Home

Context Navigation

timeseries.ncl @ 1985

Last change on this file since 1985 was 1559, checked in by heinze, 10 years ago
Changes to allow for using NCL version 6.2.1 and higher. Backward compatibility is also ensured.
File size: 52.2 KB

Rev	Line
[194]	1	load "$NCARG_ROOT/lib/ncarg/nclex/gsun/gsn_code.ncl"
[154]	2	load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
	3	load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"
	4	load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/shea_util.ncl"
[1559]	5
[190]	6	;***************************************************
[1559]	7	; Checking the kind of the script
[190]	8	;***************************************************
[194]	9
	10	function which_script()
	11	local script
	12	begin
[1559]	13	script = "timeseries"
[194]	14	return(script)
	15	end
[418]	16
[1559]	17	;***************************************************
	18	; Retrieving the NCL version used
	19	;***************************************************
	20
	21	ncl_version_ch = systemfunc("ncl -V")
	22	ncl_version = stringtofloat(ncl_version_ch)
	23
	24	;***************************************************
	25	; Function for checking file existence in dependence
	26	; on NCL version
	27	;***************************************************
	28
	29	function file_exists(version:string,file_name:string)
	30	begin
	31	if( version .EQ. "6.2.1" ) then
	32	existing = fileexists(file_name)
	33	else
	34	existing = isfilepresent(file_name)
	35	end if
	36	return(existing)
	37	end
	38
	39	;***************************************************
	40	; load .ncl.config or .ncl.config.default
	41	;***************************************************
	42
	43	file_name ="$PALM_BIN/../../.ncl.config"
	44	existing_file = file_exists(ncl_version_ch,file_name)
	45
	46	if (existing_file) then
[418]	47	loadscript("$PALM_BIN/../../.ncl.config")
[190]	48	else
[1559]	49	file_name = "$PALM_BIN/NCL/.ncl.config.default"
	50	existing_file = file_exists(ncl_version_ch,file_name)
	51	if (existing_file) then
	52	loadscript( "$PALM_BIN/NCL/.ncl.config.default")
	53	else
[418]	54	palm_bin_path = getenv("PALM_BIN")
[190]	55	print(" ")
[418]	56	print("Neither the personal configuration file '.ncl.config' exists in")
	57	print("~/palm/current_version")
[566]	58	print("nor the default configuration file '.ncl.config.default' "+\
	59	"exists in")
[418]	60	print(palm_bin_path + "/NCL")
[190]	61	print(" ")
	62	exit
[154]	63	end if
[1559]	64	end if
	65
[190]	66	begin
[154]	67
[190]	68	;***************************************************
[526]	69	; Retrieving the double quote character
	70	;***************************************************
	71
	72	dq=str_get_dq()
	73
	74	;***************************************************
[190]	75	; set up default parameter values and strings
	76	;***************************************************
[218]	77
[190]	78	if (file_1 .EQ. "File in") then
	79	print(" ")
[418]	80	print("Declare input file 'file_1=' in '.ncl.config' or prompt")
[190]	81	print(" ")
	82	exit
[175]	83	else
	84	file_in = file_1
[154]	85	end if
[175]	86
[566]	87	if (format_out .NE. "x11" .AND. format_out .NE. "pdf" .AND. \
	88	format_out .NE. "eps" .AND. format_out .NE. "ps" .AND. \
	89	format_out .NE. "epsi" .AND. format_out .NE. "ncgm" .AND. \
	90	format_out .NE. "png")then
[190]	91	print(" ")
	92	print("'format_out = "+format_out+"' is invalid and set to'x11'")
	93	print(" ")
	94	format_out="x11"
[154]	95	end if
[175]	96
[532]	97	if (ncl_version .LE. 5.1 .AND. format_out .EQ. "png") then
	98	print(" ")
	99	print("Output of png files not available")
	100	print("png output is avaiable with NCL version 5.2.0 and higher ")
	101	print("NCL version used: " + ncl_version_ch)
	102	print(" ")
	103	exit
	104	end if
	105
[190]	106	if (over .NE. 0 .AND. over .NE. 1) then
	107	print(" ")
	108	print("'over'= "+over+" is invalid and set to 0")
	109	print(" ")
[161]	110	over = 0
[218]	111	end if
	112
	113	if (norm_t .EQ. 0) then
	114	print(" ")
	115	print("Normalising with 0 is not allowed, 'norm_t' is set to 1.0")
	116	print(" ")
	117	norm_t = 1.0
	118	end if
[190]	119
[175]	120
[190]	121	;***************************************************
[154]	122	; open input file
[190]	123	;***************************************************
[218]	124
	125	file_in_1 = False
	126	if (isStrSubset(file_in, ".nc"))then
	127	start_f = -2
	128	end_f = -2
	129	file_in_1 = True
	130	end if
[154]	131
[218]	132	if (start_f .EQ. -1)then
	133	print(" ")
[566]	134	print("'start_f' must be one of the cyclic numbers (at least 0) of "+\
	135	"your input file(s)")
[218]	136	print(" ")
	137	exit
	138	end if
	139	if (end_f .EQ. -1)then
	140	print(" ")
[566]	141	print("'end_f' must be one of the cyclic numbers (at least 0) of "+\
	142	"your input file(s)")
[218]	143	print(" ")
	144	exit
	145	end if
[154]	146
[218]	147	files=new(end_f-start_f+1,string)
	148
	149	if (file_in_1)then
	150	if (isfilepresent(file_in))then
	151	files(0)=file_in
	152	else
	153	print(" ")
	154	print("1st input file: '"+file_in+"' does not exist")
	155	print(" ")
	156	exit
	157	end if
	158	else
	159	if (start_f .EQ. 0)then
	160	if (isfilepresent(file_in+".nc"))then
	161	files(0)=file_in+".nc"
	162	do i=1,end_f
	163	if (isfilepresent(file_in+"."+i+".nc"))then
	164	files(i)=file_in+"."+i+".nc"
	165	else
	166	print(" ")
	167	print("Input file: '"+file_in+"."+i+".nc' does not exist")
	168	print(" ")
	169	exit
	170	end if
	171	end do
	172	else
	173	print(" ")
	174	print("Input file: '"+file_in+".nc' does not exist")
	175	print(" ")
	176	exit
	177	end if
	178	else
	179	do i=start_f,end_f
	180	if (isfilepresent(file_in+"."+i+".nc"))then
	181	files(i-start_f)=file_in+"."+i+".nc"
	182	else
	183	print(" ")
	184	print("Input file: '"+file_in+"."+i+".nc' does not exist")
	185	print(" ")
	186	exit
	187	end if
	188	end do
	189	end if
	190	end if
	191
	192	f=addfiles(files,"r")
	193	f_att=addfile(files(0),"r")
	194	ListSetType(f,"cat")
	195
[585]	196	vNam =getfilevarnames(f_att)
[1248]	197	if( ncl_version .GE. 6.1 ) then
	198	vNam = vNam(::-1)
	199	end if
[585]	200	vType = getfilevartypes(f_att,vNam)
	201	vNam_static = vNam
[218]	202
[585]	203	if ((all(vType .eq. "double"))) then ;distinction if data is double or float
	204	check_vType = True
	205	else
	206	check_vType = False
	207	end if
	208
[154]	209	print(" ")
[190]	210	print("Variables in input file:")
	211	print("- "+ vNam)
[154]	212	print(" ")
[161]	213	dim = dimsizes(vNam)
	214	if (dim .EQ. 0) then
	215	print(" ")
[310]	216	print("There is no data on file")
[161]	217	print(" ")
	218	end if
	219
[310]	220	t_all = f[:]->time
	221	nt = dimsizes(t_all)
[218]	222
[827]	223	if (nt .EQ. 1)then
	224	delta_t=t_all(nt-1)/nt
	225	else
	226	delta_t_array = new(nt-1,double)
	227
	228	do i=0,nt-2
	229	delta_t_array(i) = t_all(i+1)-t_all(i)
	230	end do
	231
	232	delta_t = min(delta_t_array)
	233	delete(delta_t_array)
	234	end if
	235
[769]	236	if (nt .LE. 1) then
	237	print(" ")
	238	print("Input file contains only one time step -> " +\
	239	"plot of timeseries is not possible ")
	240	print(" ")
	241	exit
	242	end if
	243
[190]	244	;****************************************************
[154]	245	; start of time step and different types of mistakes that could be done
[190]	246	;****************************************************
[154]	247
[190]	248	if (start_time_step .EQ. -1.) then
	249	start_time_step=t_all(0)/3600
[154]	250	else
[162]	251	if (start_time_step .GE. t_all(nt-1)/3600)
[154]	252	print(" ")
[566]	253	print("'start_time_step' = "+ start_time_step +\
	254	"h is equal or greater than last time step = " + \
	255	t_all(nt-1)+"s = "+t_all(nt-1)/3600+"h")
[154]	256	print(" ")
[218]	257	print("Select another 'start_time_step'")
[162]	258	print(" ")
[154]	259	exit
	260	end if
[162]	261	if (start_time_step .LT. t_all(0)/3600)
[154]	262	print(" ")
[566]	263	print("'start_time_step' = "+ start_time_step +\
	264	"h is lower than first time step = " + t_all(0)+"s = "+\
	265	t_all(0)/3600+"h")
[154]	266	print(" ")
[218]	267	print("Select another 'start_time_step'")
[162]	268	print(" ")
[154]	269	exit
	270	end if
	271	end if
[162]	272	do i=0,nt-2
[566]	273	if (start_time_step .GE. (t_all(i)-delta_t/2)/3600 .AND. \
	274	start_time_step .LT. (t_all(i)+delta_t/2)/3600)then
[162]	275	st=i
	276	break
[769]	277	else
	278	st=0
[162]	279	end if
	280	end do
[566]	281	if (start_time_step .GE. t_all(nt-1)-delta_t/2 .AND. \
	282	start_time_step .LT. t_all(nt-1)) then
[162]	283	st=nt-2
[769]	284	else
	285	st=0
[162]	286	end if
	287
[190]	288	;****************************************************
[154]	289	; end of time step and different types of mistakes that could be done
[190]	290	;****************************************************
[154]	291
[190]	292	if (end_time_step .EQ. -1.) then
	293	end_time_step = t_all(nt-1)/3600
[154]	294	else
[162]	295	if (end_time_step .LE. t_all(0)/3600)
[154]	296	print(" ")
[566]	297	print("'end_time_step' = "+end_time_step+ \
	298	"h is lower or equal than first time step = " + \
	299	t_all(0)+"s = "+t_all(0)/3600+"h")
[154]	300	print(" ")
[218]	301	print("Select another 'end_time_step'")
[162]	302	print(" ")
[154]	303	exit
	304	end if
[162]	305	if (end_time_step .GT. t_all(nt-1)/3600)
[154]	306	print(" ")
[566]	307	print("'end_time_step' = "+ end_time_step +\
	308	"h is greater than last time step = " + t_all(nt-1)+"s = "+\
	309	t_all(nt-1)/3600+"h")
[154]	310	print(" ")
[218]	311	print("Select another 'end_time_step'")
[162]	312	print(" ")
[154]	313	exit
	314	end if
[162]	315	if (end_time_step .LE. start_time_step/3600)
[161]	316	print(" ")
[566]	317	print("'end_time_step' = "+ end_time_step +\
	318	"h is equal or lower than 'start_time_step' = "+\
	319	start_time_step+"h")
[161]	320	print(" ")
[218]	321	print("Select another 'start_time_step' or 'end_time_step'")
[162]	322	print(" ")
[161]	323	exit
[154]	324	end if
[162]	325	end if
	326	do i=0,nt-1
[566]	327	if (end_time_step .GE. (t_all(i)-delta_t/2)/3600 .AND. \
	328	end_time_step .LT. (t_all(i)+delta_t/2)/3600)then
[162]	329	et=i
	330	break
[769]	331	else
	332	et=0
[162]	333	end if
	334	end do
	335
	336	delete(start_time_step)
	337	start_time_step=round(st,3)
	338	delete(end_time_step)
	339	end_time_step=round(et,3)
[154]	340
[175]	341	print(" ")
[566]	342	print("Output of time steps from "+t_all(start_time_step)/3600+" h = "+\
	343	t_all(start_time_step)+" s => index = "+start_time_step)
	344	print(" till "+t_all(end_time_step)/3600+" h = "+\
	345	t_all(end_time_step)+" s => index = "+end_time_step)
[175]	346	print(" ")
	347
[218]	348	t = t_all(start_time_step:end_time_step)/norm_t
[154]	349
	350	; ***************************************************
	351	; set up recourses
	352	; ***************************************************
	353
	354	res = True
	355	res@gsnDraw = False
	356	res@gsnFrame = False
	357	res@tmXBMode = True
	358	res@tmYLMode = True
	359	res@txFont = "helvetica"
	360	res@tiMainFont = "helvetica"
	361	res@tiXAxisFont = "helvetica"
	362	res@tiYAxisFont = "helvetica"
	363	res@tmXBLabelFont = "helvetica"
[161]	364	res@tmYLLabelFont = "helvetica"
	365	res@xyLineColors = (/237/)
[758]	366	res@trXMaxF = t(start_time_step)
	367	res@trXMinF = t(end_time_step)
[161]	368
[218]	369	res@lgLabelFontHeightF = 0.02
[161]	370
[983]	371	resP = True
	372	resP@gsnMaximize = True
	373	resP@gsnPanelXWhiteSpacePercent = 4.0
	374	resP@gsnPanelYWhiteSpacePercent = 4.0
	375	resP@txFont = "helvetica"
	376	resP@txString = f_att@title+" time series "
	377	resP@txFuncCode = "~"
	378	if (format_out .NE. "x11") then
	379	if (no_rows .GE. 5) then
	380	resP@txFontHeightF = 0.0105
	381	else
	382	resP@txFontHeightF = 0.015
	383	end if
	384	else
	385	resP@txFontHeightF = 0.015
	386	end if
[218]	387	res@tmXBMinorPerMajor = 4
	388	res@tmYLMinorPerMajor = 4
	389
[154]	390	res@vpWidthF=4
	391
[161]	392	txres = True
[154]	393
	394	; ***************************************************
	395	; read data and create plots
	396	; ***************************************************
	397
[532]	398	if ( format_out .EQ. "pdf" .OR. format_out .EQ. "ps" ) then
	399	format_out@wkPaperSize = "A4"
	400	end if
	401	if ( format_out .EQ. "png" ) then
	402	format_out@wkWidth = 1000
	403	format_out@wkHeight = 1000
	404	end if
	405
[154]	406	wks_ps = gsn_open_wks(format_out,file_out)
	407	gsn_define_colormap(wks_ps,"rainbow+white")
[161]	408	plot_ps=new(dim,graphic)
	409
[154]	410	n=0
[585]	411	if (check_vType) then
	412	minE = 100000.d
	413	maxE = -100000.d
	414	minus = 100000.d
	415	maxus = -100000.d
	416	minu = 100000.d
	417	maxu = -100000.d
	418	minz = 100000.d
	419	maxz = -100000.d
	420	minw = 100000.d
	421	maxw = -100000.d
	422	minp = 100000.d
	423	maxp = -100000.d
	424	mins = 100000.d
	425	maxs = -100000.d
	426	else
	427	minE = 100000.
	428	maxE = -100000.
	429	minus = 100000.
	430	maxus = -100000.
	431	minu = 100000.
	432	maxu = -100000.
	433	minz = 100000.
	434	maxz = -100000.
	435	minw = 100000.
	436	maxw = -100000.
	437	minp = 100000.
	438	maxp = -100000.
	439	mins = 100000.
	440	maxs = -100000.
	441	end if
[161]	442
[585]	443	if (check_vType) then
	444	data = new((/dim,(end_time_step-start_time_step)+1/),double)
	445	data_0 = new((end_time_step-start_time_step)+1,double)
	446	mini = new(dim,double)
	447	maxi = new(dim,double)
	448	else
	449	data = new((/dim,(end_time_step-start_time_step)+1/),float)
	450	data_0 = new((end_time_step-start_time_step)+1,float)
	451	mini = new(dim,float)
	452	maxi = new(dim,float)
	453	end if
	454	unit = new(dim,string)
[161]	455	data_0 = 0.0
[585]	456
[161]	457
	458	if (over .EQ. 1) then
	459	plot_E = gsn_csm_xy(wks_ps,t,data_0(:),res)
	460	plot_Es = gsn_csm_xy(wks_ps,t,data_0(:),res)
	461	plot_us = gsn_csm_xy(wks_ps,t,data_0(:),res)
	462	plot_ws = gsn_csm_xy(wks_ps,t,data_0(:),res)
	463	plot_umax = gsn_csm_xy(wks_ps,t,data_0(:),res)
	464	plot_vmax = gsn_csm_xy(wks_ps,t,data_0(:),res)
	465	plot_wmax = gsn_csm_xy(wks_ps,t,data_0(:),res)
	466	plot_z_i_wpt = gsn_csm_xy(wks_ps,t,data_0(:),res)
	467	plot_z_i_pt = gsn_csm_xy(wks_ps,t,data_0(:),res)
	468	plot_wpptp0 = gsn_csm_xy(wks_ps,t,data_0(:),res)
	469	plot_wpptp = gsn_csm_xy(wks_ps,t,data_0(:),res)
	470	plot_wpt = gsn_csm_xy(wks_ps,t,data_0(:),res)
	471	plot_pt_0_ = gsn_csm_xy(wks_ps,t,data_0(:),res)
	472	plot_pt_zp_ = gsn_csm_xy(wks_ps,t,data_0(:),res)
	473	end if
[162]	474
[161]	475
[310]	476	if (var .NE. "all") then
	477	vNam_temp = new((/dim/),string)
	478	comma = 0
	479	var_char = stringtocharacter(var)
	480	no_char = dimsizes(var_char)
	481
	482	do j = 0, no_char-1
	483	if(var_char(j) .eq. ",")
	484	comma = comma + 1
	485	end if
	486	end do
	487
[317]	488	if(comma .eq. 0 .or. comma .eq. 1)
	489	print(" ")
	490	print("The variables 'var="+var+"'" )
	491	print("do not exist on your input file;")
	492	print("be sure to have one comma before and after each variable")
	493	print(" ")
	494	exit
	495	end if
	496
[310]	497	indices = new((/comma/), integer)
	498	comma = 0
	499
	500	do j = 0, no_char-1
	501	if(var_char(j) .eq. ",")
	502	indices(comma) = j
	503	comma = comma + 1
	504	end if
	505	end do
	506
	507	do j = 0, comma -2
[566]	508	vNam_temp(j) = charactertostring(\
	509	var_char(indices(j)+1:indices(j+1)-1))
[310]	510	end do
	511
	512	delete(vNam)
	513	vNam = new((/comma-1/),string)
	514	do varn = 0, comma - 2
	515	vNam(varn) = vNam_temp(comma - 2 - varn)
	516	end do
	517
	518	do j = 0, comma-2
	519	count_check = 0
	520	do varn = 0, dim-1
	521	if(vNam(j) .ne. vNam_static(varn))
	522	count_check=count_check+1
	523	end if
[317]	524	end do
	525
[310]	526	if(count_check .eq. dim)
	527	print(" ")
	528	print("The variables 'var="+var+"'" )
	529	print("do not exist on your input file;")
	530	print("be sure to have one comma before and after each variable")
	531	print(" ")
	532	exit
	533	end if
	534	end do
	535
	536	dim = comma-1
	537	end if
	538
	539
	540	do varn = dim-1, 0, 1
	541
[161]	542	if( isStrSubset (vNam(varn), "time") )
	543	check = False
	544	else
	545	check = True
[310]	546	end if
[161]	547
	548	if(check) then
[162]	549
[218]	550	data_all = f[:]->$vNam(varn)$
	551	data_att = f_att->$vNam(varn)$
	552	unit(varn) = data_att@units
[310]	553
	554
[161]	555	data(varn,:)=data_all(start_time_step:end_time_step)
	556
	557	if (over .EQ. 1) then
	558
	559	mini(varn) = min(data(varn,:))
	560	maxi(varn) = max(data(varn,:))
	561
[566]	562	if (vNam(varn) .EQ. "E" .OR. vNam(varn) .EQ. "Es" .OR. \
	563	vNam(varn) .EQ. "E*") then
[161]	564	if (mini(varn) .EQ. maxi(varn)) then
	565	if (min(data(varn,:)) .EQ. 0)then
[218]	566	mini(varn)= mini(varn)-0.1
	567	maxi(varn)= maxi(varn)+0.1
[161]	568	end if
	569	if (min(data(varn,:)) .LT. 0)then
	570	mini(varn)= mini(varn)-1.+(mini(varn))/2
	571	maxi(varn)= maxi(varn)+1.-(maxi(varn))/2
	572	end if
	573	if (min(data(varn,:)) .GT. 0)then
	574	mini(varn)= mini(varn)-1.-(mini(varn))/2
	575	maxi(varn)= maxi(varn)+1.+(maxi(varn))/2
	576	end if
	577	end if
	578	minE=min((/minE,mini(varn)/))
	579	maxE=max((/maxE,maxi(varn)/))
	580	end if
	581
[526]	582	if (vNam(varn) .EQ. "us" .OR. vNam(varn) .EQ. "ws"\
	583	.OR. vNam(varn) .EQ. "u" .OR. vNam(varn) .EQ. "w") then
[161]	584	if (mini(varn) .EQ. maxi(varn)) then
	585	if (min(data(varn,:)) .EQ. 0)then
[218]	586	mini(varn)= mini(varn)-0.1
	587	maxi(varn)= maxi(varn)+0.1
[161]	588	end if
	589	if (min(data(varn,:)) .LT. 0)then
	590	mini(varn)= mini(varn)-1.+(mini(varn))/2
	591	maxi(varn)= maxi(varn)+1.-(maxi(varn))/2
	592	end if
	593	if (min(data(varn,:)) .GT. 0)then
	594	mini(varn)= mini(varn)-1.-(mini(varn))/2
	595	maxi(varn)= maxi(varn)+1.+(maxi(varn))/2
	596	end if
	597	end if
	598	minus=min((/minus,mini(varn)/))
	599	maxus=max((/maxus,maxi(varn)/))
	600	end if
	601
[566]	602	if (vNam(varn) .EQ. "umax" .OR. vNam(varn) .EQ. "vmax" .OR. \
	603	vNam(varn) .EQ. "wmax") then
[161]	604	if (mini(varn) .EQ. maxi(varn)) then
	605	if (mini(varn) .EQ. 0)then
[218]	606	mini(varn)= mini(varn)-0.1
	607	maxi(varn)= maxi(varn)+0.1
[161]	608	end if
	609	if (mini(varn) .LT. 0)then
	610	mini(varn)= mini(varn)-1.+(mini(varn))/2
	611	maxi(varn)= maxi(varn)+1.-(maxi(varn))/2
	612	end if
	613	if (mini(varn) .GT. 0)then
	614	mini(varn)= mini(varn)-1.-(mini(varn))/2
	615	maxi(varn)= maxi(varn)+1.+(maxi(varn))/2
	616	end if
	617	end if
	618	minu=min((/minu,mini(varn)/))
	619	maxu=max((/maxu,maxi(varn)/))
	620	end if
	621
	622	if (vNam(varn) .EQ. "z_i_wpt" .OR. vNam(varn) .EQ. "z_i_pt") then
	623	if (mini(varn) .EQ. maxi(varn)) then
	624	if (min(data(varn,:)) .EQ. 0)then
[218]	625	mini(varn)= mini(varn)-0.1
	626	maxi(varn)= maxi(varn)+0.1
[161]	627	end if
	628	if (min(data(varn,:)) .LT. 0)then
	629	mini(varn)= mini(varn)-1.+(mini(varn))/2
	630	maxi(varn)= maxi(varn)+1.-(maxi(varn))/2
	631	end if
	632	if (min(data(varn,:)) .GT. 0)then
	633	mini(varn)= mini(varn)-1.-(mini(varn))/2
	634	maxi(varn)= maxi(varn)+1.+(maxi(varn))/2
	635	end if
	636	end if
	637	minz=min((/minz,mini(varn)/))
	638	maxz=max((/maxz,maxi(varn)/))
	639	end if
	640
[566]	641	if (vNam(varn) .EQ. "wpptp0" .OR. vNam(varn) .EQ. "wpptp" .OR. \
	642	vNam(varn) .EQ. "wpt" .OR. vNam(varn) .EQ. \
	643	"w"+dq+"pt"+dq+"0" .OR. vNam(varn) .EQ. "w"+dq+"pt"+dq ) then
[161]	644	if (mini(varn) .EQ. maxi(varn)) then
	645	if (min(data(varn,:)) .EQ. 0)then
[218]	646	mini(varn)= mini(varn)-0.1
	647	maxi(varn)= maxi(varn)+0.1
[161]	648	end if
	649	if (min(data(varn,:)) .LT. 0)then
	650	mini(varn)= mini(varn)-1.+(mini(varn))/2
	651	maxi(varn)= maxi(varn)+1.-(maxi(varn))/2
	652	end if
	653	if (min(data(varn,:)) .GT. 0)then
	654	mini(varn)= mini(varn)-1.-(mini(varn))/2
	655	maxi(varn)= maxi(varn)+1.+(maxi(varn))/2
	656	end if
	657	end if
	658	minw=min((/minw,mini(varn)/))
	659	maxw=max((/maxw,maxi(varn)/))
	660	end if
	661
[526]	662	if (vNam(varn) .EQ. "pt_0_" .OR. vNam(varn) .EQ. "pt_zp_"\
	663	.OR. vNam(varn) .EQ. "pt(0)" .OR. vNam(varn) .EQ. "pt(zp)") then
[161]	664	if (mini(varn) .EQ. maxi(varn)) then
	665	if (min(data(varn,:)) .EQ. 0)then
[218]	666	mini(varn)= mini(varn)-0.1
	667	maxi(varn)= maxi(varn)+0.1
[161]	668	end if
	669	if (min(data(varn,:)) .LT. 0)then
	670	mini(varn)= mini(varn)-1.+(mini(varn))/2
	671	maxi(varn)= maxi(varn)+1.-(maxi(varn))/2
	672	end if
	673	if (min(data(varn,:)) .GT. 0)then
	674	mini(varn)= mini(varn)-1.-(mini(varn))/2
	675	maxi(varn)= maxi(varn)+1.+(maxi(varn))/2
	676	end if
	677	end if
	678	minp=min((/minp,mini(varn)/))
	679	maxp=max((/maxp,maxi(varn)/))
	680	end if
	681
	682	end if
	683	end if
	684	end do
[310]	685
	686	if (isStrSubset(data@long_name," SR " ))then
[162]	687	print(" ")
[566]	688	print("If you have outputs of statistic regions you cannot overlay "+\
	689	"variables;")
[310]	690	print("'over' is set to 0" )
[162]	691	print(" ")
[310]	692	over = 0
[175]	693	end if
[162]	694
[310]	695	do varn = dim-1, 0, 1
[161]	696
[154]	697	if( isStrSubset (vNam(varn), "time") )
	698	check = False
[161]	699	else
	700	check = True
	701	end if
[190]	702
	703	if (var .NE. "all") then
[161]	704	check = isStrSubset( var,","+vNam(varn)+"," )
	705	end if
[154]	706
	707	if(check) then
[310]	708
[161]	709	if (over .EQ. 1) then
[154]	710
[161]	711	res@gsnLeftString = "overlayed plot"
	712	res@gsnRightString = unit(varn)
[218]	713	if (norm_t .NE. 1.)then
[534]	714	res@tiXAxisString = "t ("+unit_t+")"
[218]	715	else
[534]	716	res@tiXAxisString = "t (s)"
[218]	717	end if
	718	res@tiYAxisString = " "
	719	res@tiXAxisFontHeightF = font_size
	720	res@txFontHeightF = font_size
	721	res@tiYAxisFontHeightF = font_size
[154]	722
[161]	723	if (vNam(varn) .EQ. "E")
	724	E=0
	725	res@xyLineColors = (/237/)
	726	res@xyLineLabelFontHeightF = 0.05
	727	res@xyLineLabelFontColor = 237
[310]	728	res@trYMaxF = maxE
	729	res@trYMinF = minE
[161]	730	plot_E = gsn_csm_xy(wks_ps,t,data(varn,:),res)
	731	end if
[526]	732	if (vNam(varn) .EQ. "Es" .OR. vNam(varn) .EQ. "E*")
[161]	733	Es=0
	734	res@xyLineColors = (/144/)
	735	res@xyLineLabelFontHeightF = 0.05
	736	res@xyLineLabelFontColor = 144
	737	plot_Es = gsn_csm_xy(wks_ps,t,data(varn,:),res)
	738	end if
	739
[526]	740	if (vNam(varn) .EQ. "us" .OR. vNam(varn) .EQ. "u*")
[161]	741	us=0
	742	res@xyLineColors = (/237/)
	743	res@xyLineLabelFontHeightF = 0.05
	744	res@xyLineLabelFontColor = 237
[310]	745	res@trYMaxF = maxus
	746	res@trYMinF = minus
[161]	747	plot_us = gsn_csm_xy(wks_ps,t,data(varn,:),res)
	748	end if
[526]	749	if (vNam(varn) .EQ. "ws" .OR. vNam(varn) .EQ. "w*")
[161]	750	ws=0
	751	res@xyLineColors = (/144/)
	752	res@xyLineLabelFontHeightF = 0.05
	753	res@xyLineLabelFontColor = 144
	754	plot_ws = gsn_csm_xy(wks_ps,t,data(varn,:),res)
	755	end if
	756
	757	if (vNam(varn) .EQ. "umax")
	758	u=0
	759	res@xyLineColors = (/237/)
	760	res@xyLineLabelFontHeightF = 0.05
	761	res@xyLineLabelFontColor = 237
[310]	762	res@trYMaxF = maxu
	763	res@trYMinF = minu
[161]	764	plot_umax = gsn_csm_xy(wks_ps,t,data(varn,:),res)
	765	end if
	766	if (vNam(varn) .EQ. "vmax")
	767	v=0
	768	res@xyLineColors = (/144/)
	769	res@xyLineLabelFontHeightF = 0.05
	770	res@xyLineLabelFontColor = 144
	771	plot_vmax = gsn_csm_xy(wks_ps,t,data(varn,:),res)
	772	end if
	773	if (vNam(varn) .EQ. "wmax")
	774	w=0
	775	res@xyLineColors = (/80/)
	776	res@xyLineLabelFontHeightF = 0.05
	777	res@xyLineLabelFontColor = 80
	778	plot_wmax = gsn_csm_xy(wks_ps,t,data(varn,:),res)
	779	end if
[154]	780
[161]	781	if (vNam(varn) .EQ. "z_i_wpt")
	782	zw=0
	783	res@xyLineColors = (/237/)
	784	res@xyLineLabelFontHeightF = 0.05
	785	res@xyLineLabelFontColor = 237
[310]	786	res@trYMaxF = maxz
	787	res@trYMinF = minz
[161]	788	plot_z_i_wpt = gsn_csm_xy(wks_ps,t,data(varn,:),res)
	789	end if
	790	if (vNam(varn) .EQ. "z_i_pt")
	791	z=0
	792	res@xyLineColors = (/144/)
	793	res@xyLineLabelFontHeightF = 0.05
	794	res@xyLineLabelFontColor = 144
	795	plot_z_i_pt = gsn_csm_xy(wks_ps,t,data(varn,:),res)
	796	end if
	797
[566]	798	if (vNam(varn) .EQ. "wpptp0" .OR. \
	799	vNam(varn) .EQ. "w"+dq+"pt"+dq+"0" )
[161]	800	w0=0
	801	res@xyLineColors = (/237/)
	802	res@xyLineLabelFontHeightF = 0.05
	803	res@xyLineLabelFontColor = 237
[310]	804	res@trYMaxF = maxw
	805	res@trYMinF = minw
[161]	806	plot_wpptp0 = gsn_csm_xy(wks_ps,t,data(varn,:),res)
	807	end if
[526]	808	if (vNam(varn) .EQ. "wpptp" .OR. vNam(varn) .EQ. "w"+dq+"pt"+dq)
[161]	809	wp=0
	810	res@xyLineColors = (/144/)
	811	res@xyLineLabelFontHeightF = 0.05
	812	res@xyLineLabelFontColor = 144
	813	plot_wpptp = gsn_csm_xy(wks_ps,t,data(varn,:),res)
	814	end if
	815	if (vNam(varn) .EQ. "wpt")
	816	wt=0
	817	res@xyLineColors = (/80/)
	818	res@xyLineLabelFontHeightF = 0.05
	819	res@xyLineLabelFontColor = 80
	820	plot_wpt = gsn_csm_xy(wks_ps,t,data(varn,:),res)
	821	end if
	822
[526]	823	if (vNam(varn) .EQ. "pt_0_" .OR. vNam(varn) .EQ. "pt(0)")
[161]	824	p=0
	825	res@xyLineColors = (/237/)
	826	res@xyLineLabelFontHeightF = 0.05
	827	res@xyLineLabelFontColor = 237
[310]	828	res@trYMaxF = maxp
	829	res@trYMinF = minp
[161]	830	plot_pt_0_ = gsn_csm_xy(wks_ps,t,data(varn,:),res)
	831	end if
[526]	832	if (vNam(varn) .EQ. "pt_zp_" .OR. vNam(varn) .EQ. "pt(zp)")
[161]	833	pz=0
	834	res@xyLineColors = (/144/)
	835	res@xyLineLabelFontHeightF = 0.05
	836	res@xyLineLabelFontColor = 144
	837	plot_pt_zp_ = gsn_csm_xy(wks_ps,t,data(varn,:),res)
	838	end if
	839
	840	end if
[154]	841	end if
[161]	842	end do
[310]	843
	844
	845	do varn = dim-1, 0, 1
[161]	846
	847	if( isStrSubset (vNam(varn), "time") )
	848	check = False
	849	else
	850	check = True
	851	end if
[190]	852
	853	if (var.NE. "all") then
[161]	854	check = isStrSubset( var,","+vNam(varn)+"," )
	855	end if
	856
	857	if(check) then
	858
	859	if (over .EQ. 1) then
	860
	861	if (vNam(varn) .EQ. "E" .AND. Es .NE. 1) then
	862	E=1
	863	overlay(plot_E,plot_Es)
	864	n=n+1
[162]	865	plot_ps(n) = plot_E
	866
	867	; ***************************************************
[566]	868	; legend for overlaid plot
[162]	869	; ***************************************************
	870
	871	lgres = True
	872	lgMonoDashIndex = False
	873	lgres@lgLabelFont = "helvetica"
	874	lgres@lgLabelFontHeightF = .1
	875	lgres@vpWidthF = 0.4
	876	lgres@vpHeightF = 0.4
	877	lgres@lgDashIndexes = (/0,0,0/)
	878	lgres@lgLineColors = (/237,144,80/)
[526]	879	lbid = gsn_create_legend(wks_ps,2,(/"E","E*"/),lgres)
[162]	880
	881	amres = True
	882	amres@amParallelPosF = 0.6
	883	amres@amOrthogonalPosF = -0.2
	884	annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)
[161]	885	end if
[566]	886	if ((vNam(varn) .EQ. "Es" .OR. vNam(varn) .EQ. "E*") .AND. \
	887	E .NE. 1) then
[161]	888	Es=1
	889	overlay(plot_E,plot_Es)
	890	n=n+1
[162]	891	plot_ps(n) = plot_E
	892
	893	; ***************************************************
[566]	894	; legend for overlaid plot
[162]	895	; ***************************************************
	896
	897	lgres = True
	898	lgMonoDashIndex = False
	899	lgres@lgLabelFont = "helvetica"
	900	lgres@lgLabelFontHeightF = .1
	901	lgres@vpWidthF = 0.4
	902	lgres@vpHeightF = 0.4
	903	lgres@lgDashIndexes = (/0,0,0/)
	904	lgres@lgLineColors = (/237,144,80/)
[526]	905	lbid = gsn_create_legend(wks_ps,2,(/"E","E*"/),lgres)
[162]	906
	907	amres = True
	908	amres@amParallelPosF = 0.6
	909	amres@amOrthogonalPosF = -0.2
	910	annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)
[161]	911	end if
	912
[566]	913	if ((vNam(varn) .EQ. "us" .OR. vNam(varn) .EQ. "u*") .AND. \
	914	ws .NE. 1) then
[161]	915	us=1
	916	overlay(plot_us,plot_ws)
	917	n=n+1
	918	plot_ps(n) = plot_us
[162]	919
	920	; ***************************************************
[566]	921	; legend for overlaid plot
[162]	922	; ***************************************************
	923
	924	lgres = True
	925	lgMonoDashIndex = False
	926	lgres@lgLabelFont = "helvetica"
	927	lgres@lgLabelFontHeightF = .1
	928	lgres@vpWidthF = 0.4
	929	lgres@vpHeightF = 0.4
	930	lgres@lgDashIndexes = (/0,0,0/)
	931	lgres@lgLineColors = (/237,144,80/)
[526]	932	lbid = gsn_create_legend(wks_ps,2,(/"u","w"/),lgres)
[162]	933
	934	amres = True
	935	amres@amParallelPosF = 0.6
	936	amres@amOrthogonalPosF = -0.2
	937	annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)
[161]	938	end if
[566]	939	if ((vNam(varn) .EQ. "ws" .OR. vNam(varn) .EQ. "w*") .AND. \
	940	us .NE. 1) then
[161]	941	ws=1
	942	overlay(plot_us,plot_ws)
	943	n=n+1
	944	plot_ps(n) = plot_us
[162]	945
	946	; ***************************************************
[566]	947	; legend for overlaid plot
[162]	948	; ***************************************************
	949
	950	lgres = True
	951	lgMonoDashIndex = False
	952	lgres@lgLabelFont = "helvetica"
	953	lgres@lgLabelFontHeightF = .1
	954	lgres@vpWidthF = 0.4
	955	lgres@vpHeightF = 0.4
	956	lgres@lgDashIndexes = (/0,0,0/)
	957	lgres@lgLineColors = (/237,144,80/)
[526]	958	lbid = gsn_create_legend(wks_ps,2,(/"u","w"/),lgres)
[162]	959
	960	amres = True
	961	amres@amParallelPosF = 0.6
	962	amres@amOrthogonalPosF = -0.2
	963	annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)
[161]	964	end if
	965
	966	if (vNam(varn) .EQ. "umax" .AND. v .NE. 1)
	967	if (w .NE. 1) then
	968	u=1
	969	overlay(plot_umax,plot_vmax)
	970	overlay(plot_umax,plot_wmax)
[310]	971
[161]	972	n=n+1
	973	plot_ps(n) = plot_umax
[162]	974
	975	; ***************************************************
[566]	976	; legend for overlaid plot
[162]	977	; ***************************************************
	978
	979	lgres = True
	980	lgMonoDashIndex = False
	981	lgres@lgLabelFont = "helvetica"
	982	lgres@lgLabelFontHeightF = .1
	983	lgres@vpWidthF = 0.4
	984	lgres@vpHeightF = 0.4
	985	lgres@lgDashIndexes = (/0,0,0/)
	986	lgres@lgLineColors = (/237,144,80/)
[566]	987	lbid = gsn_create_legend(\
	988	wks_ps,3,(/"umax","vmax","wmax"/),lgres)
[162]	989
	990	amres = True
	991	amres@amParallelPosF = 0.6
	992	amres@amOrthogonalPosF = -0.2
	993	annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)
[161]	994	end if
	995	end if
	996	if (vNam(varn) .EQ. "vmax" .AND. u .NE. 1)
	997	if (w .NE. 1) then
	998	v=1
	999	overlay(plot_umax,plot_vmax)
	1000	overlay(plot_umax,plot_wmax)
	1001	n=n+1
	1002	plot_ps(n) = plot_umax
[162]	1003
	1004	; ***************************************************
[566]	1005	; legend for overlaid plot
[162]	1006	; ***************************************************
	1007
	1008	lgres = True
	1009	lgMonoDashIndex = False
	1010	lgres@lgLabelFont = "helvetica"
	1011	lgres@lgLabelFontHeightF = .1
	1012	lgres@vpWidthF = 0.4
	1013	lgres@vpHeightF = 0.4
	1014	lgres@lgDashIndexes = (/0,0,0/)
	1015	lgres@lgLineColors = (/237,144,80/)
[566]	1016	lbid = gsn_create_legend(\
	1017	wks_ps,3,(/"umax","vmax","wmax"/),lgres)
[162]	1018	amres = True
	1019	amres@amParallelPosF = 0.6
	1020	amres@amOrthogonalPosF = -0.2
	1021	annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)
[161]	1022	end if
	1023	end if
	1024	if (vNam(varn) .EQ. "wmax" .AND. v .NE. 1)
	1025	if(u .NE. 1) then
	1026	w=1
	1027	overlay(plot_umax,plot_vmax)
	1028	overlay(plot_umax,plot_wmax)
	1029	n=n+1
	1030	plot_ps(n) = plot_umax
[162]	1031
	1032	; ***************************************************
[566]	1033	; legend for overlaid plot
[162]	1034	; ***************************************************
	1035
	1036	lgres = True
	1037	lgMonoDashIndex = False
	1038	lgres@lgLabelFont = "helvetica"
	1039	lgres@lgLabelFontHeightF = .1
	1040	lgres@vpWidthF = 0.4
	1041	lgres@vpHeightF = 0.4
	1042	lgres@lgDashIndexes = (/0,0,0/)
	1043	lgres@lgLineColors = (/237,144,80/)
[566]	1044	lbid = gsn_create_legend(\
	1045	wks_ps,3,(/"umax","vmax","wmax"/),lgres)
[162]	1046	amres = True
	1047	amres@amParallelPosF = 0.6
	1048	amres@amOrthogonalPosF = -0.2
	1049	annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)
[161]	1050	end if
	1051	end if
	1052
	1053	if (vNam(varn) .EQ. "z_i_wpt" .AND. z .NE. 1) then
	1054	zw=1
	1055	overlay(plot_z_i_wpt,plot_z_i_pt)
	1056	n=n+1
[162]	1057	plot_ps(n) = plot_z_i_wpt
	1058
	1059	; ***************************************************
[566]	1060	; legend for overlaid plot
[162]	1061	; ***************************************************
	1062
	1063	lgres = True
	1064	lgMonoDashIndex = False
	1065	lgres@lgLabelFont = "helvetica"
	1066	lgres@lgLabelFontHeightF = .1
	1067	lgres@vpWidthF = 0.4
	1068	lgres@vpHeightF = 0.4
	1069	lgres@lgDashIndexes = (/0,0,0/)
	1070	lgres@lgLineColors = (/237,144,80/)
[566]	1071	lbid = gsn_create_legend(wks_ps,2,(/"z_i_wpt","z_i_pt"/),lgres)
[162]	1072
	1073	amres = True
	1074	amres@amParallelPosF = 0.6
	1075	amres@amOrthogonalPosF = -0.2
	1076	annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)
[161]	1077	end if
	1078	if (vNam(varn) .EQ. "z_i_pt" .AND. zw .NE. 1) then
	1079	z=1
	1080	overlay(plot_z_i_wpt,plot_z_i_pt)
	1081	n=n+1
[162]	1082	plot_ps(n) = plot_z_i_wpt
	1083
	1084	; ***************************************************
[566]	1085	; legend for overlaid plot
[162]	1086	; ***************************************************
	1087
	1088	lgres = True
	1089	lgMonoDashIndex = False
	1090	lgres@lgLabelFont = "helvetica"
	1091	lgres@lgLabelFontHeightF = .1
	1092	lgres@vpWidthF = 0.4
	1093	lgres@vpHeightF = 0.4
	1094	lgres@lgDashIndexes = (/0,0,0/)
	1095	lgres@lgLineColors = (/237,144,80/)
[566]	1096	lbid = gsn_create_legend(wks_ps,2,(/"z_i_wpt","z_i_pt"/),lgres)
[162]	1097
	1098	amres = True
	1099	amres@amParallelPosF = 0.6
	1100	amres@amOrthogonalPosF = -0.2
	1101	annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)
[161]	1102	end if
	1103
[526]	1104	if ((vNam(varn) .EQ. "wpptp0" .OR. vNam(varn) .EQ. "w"+dq+"pt"+dq+"0") .AND. wp .NE. 1)
[161]	1105	if (wt .NE. 1) then
	1106	w0=1
	1107	overlay(plot_wpptp0,plot_wpptp)
	1108	overlay(plot_wpptp0,plot_wpt)
	1109	n=n+1
	1110	plot_ps(n) = plot_wpptp0
[162]	1111
	1112	; ***************************************************
[566]	1113	; legend for overlaid plot
[162]	1114	; ***************************************************
	1115
	1116	lgres = True
	1117	lgMonoDashIndex = False
	1118	lgres@lgLabelFont = "helvetica"
	1119	lgres@lgLabelFontHeightF = .1
	1120	lgres@vpWidthF = 0.4
	1121	lgres@vpHeightF = 0.4
	1122	lgres@lgDashIndexes = (/0,0,0/)
	1123	lgres@lgLineColors = (/237,144,80/)
[566]	1124	lbid = gsn_create_legend(\
	1125	wks_ps,3,(/"wpptp0","wpptp","wpt"/),lgres)
[162]	1126	amres = True
	1127	amres@amParallelPosF = 0.6
	1128	amres@amOrthogonalPosF = -0.2
	1129	annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)
[161]	1130	end if
	1131	end if
[566]	1132	if ((vNam(varn) .EQ. "wpptp" .OR. \
	1133	vNam(varn) .EQ. "w"+dq+"pt"+dq) .AND. w0 .NE. 1)
[161]	1134	if (wt .NE. 1) then
	1135	wp=1
	1136	overlay(plot_wpptp0,plot_wpptp)
	1137	overlay(plot_wpptp0,plot_wpt)
	1138	n=n+1
[162]	1139	plot_ps(n) = plot_wpptp0
	1140
	1141	; ***************************************************
[566]	1142	; legend for overlaid plot
[162]	1143	; ***************************************************
	1144
	1145	lgres = True
	1146	lgMonoDashIndex = False
	1147	lgres@lgLabelFont = "helvetica"
	1148	lgres@lgLabelFontHeightF = .1
	1149	lgres@vpWidthF = 0.4
	1150	lgres@vpHeightF = 0.4
	1151	lgres@lgDashIndexes = (/0,0,0/)
	1152	lgres@lgLineColors = (/237,144,80/)
[566]	1153	lbid = gsn_create_legend(\
	1154	wks_ps,3,(/"wpptp0","wpptp","wpt"/),lgres)
[162]	1155
	1156	amres = True
	1157	amres@amParallelPosF = 0.6
	1158	amres@amOrthogonalPosF = -0.2
	1159	annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)
[161]	1160	end if
	1161	end if
	1162	if (vNam(varn) .EQ. "wpt" .AND. wp .NE. 1)
	1163	if (w0 .NE. 1) then
	1164	wt=1
	1165	overlay(plot_wpptp0,plot_wpptp)
	1166	overlay(plot_wpptp0,plot_wpt)
	1167	n=n+1
[162]	1168	plot_ps(n) = plot_wpptp0
	1169
	1170	; ***************************************************
[566]	1171	; legend for overlaid plot
[162]	1172	; ***************************************************
	1173
	1174	lgres = True
	1175	lgMonoDashIndex = False
	1176	lgres@lgLabelFont = "helvetica"
	1177	lgres@lgLabelFontHeightF = .1
	1178	lgres@vpWidthF = 0.4
	1179	lgres@vpHeightF = 0.4
	1180	lgres@lgDashIndexes = (/0,0,0/)
	1181	lgres@lgLineColors = (/237,144,80/)
[566]	1182	lbid = gsn_create_legend(\
	1183	wks_ps,3,(/"wpptp0","wpptp","wpt"/),lgres)
[162]	1184	amres = True
	1185	amres@amParallelPosF = 0.6
	1186	amres@amOrthogonalPosF = -0.2
	1187	annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)
[161]	1188	end if
	1189	end if
	1190
[566]	1191	if ((vNam(varn) .EQ. "pt_0_" .OR. vNam(varn) .EQ. "pt(0)") .AND. \
	1192	pz .NE. 1) then
[161]	1193	p=1
	1194	overlay(plot_pt_0_,plot_pt_zp_)
	1195	n=n+1
[162]	1196	plot_ps(n) = plot_pt_0_
	1197
	1198	; ***************************************************
[566]	1199	; legend for overlaid plot
[162]	1200	; ***************************************************
	1201
	1202	lgres = True
	1203	lgMonoDashIndex = False
	1204	lgres@lgLabelFont = "helvetica"
	1205	lgres@lgLabelFontHeightF = .1
	1206	lgres@vpWidthF = 0.4
	1207	lgres@vpHeightF = 0.4
	1208	lgres@lgDashIndexes = (/0,0,0/)
	1209	lgres@lgLineColors = (/237,144,80/)
[566]	1210	lbid = gsn_create_legend(wks_ps,2,(/"pt_0_","pt_zp_"/),lgres)
[162]	1211
	1212	amres = True
	1213	amres@amParallelPosF = 0.6
	1214	amres@amOrthogonalPosF = -0.2
	1215	annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)
[161]	1216	end if
[566]	1217	if ((vNam(varn) .EQ. "pt_zp_" .OR. vNam(varn) .EQ. "pt(zp)") .AND.\
	1218	p .NE. 1) then
[161]	1219	pz=1
	1220	overlay(plot_pt_0_,plot_pt_zp_)
	1221	n=n+1
[162]	1222	plot_ps(n) = plot_pt_0_
	1223
	1224	; ***************************************************
[566]	1225	; legend for overlaid plot
[162]	1226	; ***************************************************
	1227
	1228	lgres = True
	1229	lgMonoDashIndex = False
	1230	lgres@lgLabelFont = "helvetica"
	1231	lgres@lgLabelFontHeightF = .1
	1232	lgres@vpWidthF = 0.4
	1233	lgres@vpHeightF = 0.4
	1234	lgres@lgDashIndexes = (/0,0,0/)
	1235	lgres@lgLineColors = (/237,144,80/)
[566]	1236	lbid = gsn_create_legend(wks_ps,2,(/"pt_0_","pt_zp_"/),lgres)
[162]	1237
	1238	amres = True
	1239	amres@amParallelPosF = 0.6
	1240	amres@amOrthogonalPosF = -0.2
	1241	annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)
[161]	1242	end if
[162]	1243
[526]	1244	if(vNam(varn) .NE. "pt_zp_" .AND. vNam(varn) .NE. "pt(zp)" .AND. \
[566]	1245	vNam(varn) .NE. "pt_0_" .AND. vNam(varn) .NE. "pt(0)" .AND. \
	1246	vNam(varn) .NE. "wpt" .AND. vNam(varn) .NE. "wpptp" .AND. \
	1247	vNam(varn) .NE. "w"+dq+"pt"+dq .AND. \
	1248	vNam(varn) .NE. "wpptp0" .AND. \
	1249	vNam(varn) .NE. "w"+dq+"pt"+dq+"0" .AND. \
	1250	vNam(varn) .NE. "z_i_pt" .AND. \
	1251	vNam(varn) .NE. "z_i_wpt" .AND. vNam(varn) .NE. "wmax" .AND. \
	1252	vNam(varn) .NE. "vmax" .AND. vNam(varn) .NE. "umax" .AND. \
	1253	vNam(varn) .NE. "ws" .AND. vNam(varn) .NE. "w*" .AND. \
	1254	vNam(varn) .NE. "us" .AND. vNam(varn) .NE. "u*" .AND. \
	1255	vNam(varn) .NE. "Es" .AND. vNam(varn) .NE. "E*" .AND. \
[526]	1256	vNam(varn) .NE. "E") then
[162]	1257
[161]	1258	n=n+1
	1259	res@xyLineColors = (/237/)
	1260	res@xyLabelMode = False
	1261	res@gsnLeftString = vNam(varn)
	1262	res@gsnRightString = unit(varn)
[310]	1263	res@trYMaxF = max(data(varn,:))
	1264	res@trYMinF = min(data(varn,:))
[161]	1265	if (min(data(varn,:)) .EQ. max(data(varn,:))) then
	1266	if (min(data(varn,:)) .EQ. 0)then
[310]	1267	res@trYMaxF = max(data(varn,:))-0.1
	1268	res@trYMinF = min(data(varn,:))+0.1
[161]	1269	end if
	1270	if (min(data(varn,:)) .LT. 0)then
[310]	1271	res@trYMaxF = max(data(varn,:))+(min(data(varn,:)))/2
	1272	res@trYMinF = min(data(varn,:))-(max(data(varn,:)))/2
[161]	1273	end if
	1274	if (min(data(varn,:)) .GT. 0)then
[310]	1275	res@trYMaxF = max(data(varn,:))-(min(data(varn,:)))/2
	1276	res@trYMinF = min(data(varn,:))+(max(data(varn,:)))/2
[161]	1277	end if
	1278	end if
	1279	plot_ps(n) = gsn_csm_xy(wks_ps,t,data(varn,:),res)
	1280	end if
[310]	1281
	1282	; ***************************************************
	1283	; transparent legend as placeholder for single plots
	1284	; ***************************************************
	1285	lgresmono = True
	1286	lgMonoDashIndex = False
	1287	lgresmono@lgMonoLineColor = True
	1288	lgresmono@lgLabelFont = "helvetica"
	1289	lgresmono@lgLabelFontHeightF = .1
	1290	lgresmono@vpWidthF = 0.4
	1291	lgresmono@vpHeightF = 0.4
	1292	lgresmono@lgLineColor = 0
	1293	lgresmono@lgPerimOn = False
	1294	lbid = gsn_create_legend(wks_ps,1,(/""/),lgresmono)
	1295
	1296	amres = True
	1297	amres@amParallelPosF = 0.6
	1298	amres@amOrthogonalPosF = -0.2
	1299	annoid1 = gsn_add_annotation(plot_ps(n),lbid,amres)
[161]	1300
	1301	else
	1302
	1303	print("plot of " + vNam(varn))
	1304
	1305	n=n+1
	1306	res@xyLineColors = (/237/)
	1307	res@gsnLeftString = vNam(varn)
	1308	res@gsnRightString = unit(varn)
[310]	1309
[218]	1310	if (norm_t .NE. 1.)then
[534]	1311	res@tiXAxisString = "t ("+unit_t+")"
[218]	1312	else
[534]	1313	res@tiXAxisString = "t (s)"
[218]	1314	end if
	1315	res@tiYAxisString = " "
	1316	res@tiXAxisFontHeightF = font_size
	1317	res@txFontHeightF = font_size
	1318	res@tiYAxisFontHeightF = font_size
	1319	res@tmXBLabelFontHeightF = font_size
	1320	res@tmYLLabelFontHeightF = font_size
[310]	1321	res@trYMaxF = max(data(varn,:))
	1322	res@trYMinF = min(data(varn,:))
[161]	1323	if (min(data(varn,:)) .EQ. max(data(varn,:))) then
	1324	if (min(data(varn,:)) .EQ. 0)then
[310]	1325	res@trYMaxF = max(data(varn,:))-0.1
	1326	res@trYMinF = min(data(varn,:))+0.1
[161]	1327	end if
	1328	if (min(data(varn,:)) .LT. 0)then
[310]	1329	res@trYMaxF = max(data(varn,:))+(min(data(varn,:)))/2
	1330	res@trYMinF = min(data(varn,:))-(max(data(varn,:)))/2
[161]	1331	end if
	1332	if (min(data(varn,:)) .GT. 0)then
[310]	1333	res@trYMaxF = max(data(varn,:))-(min(data(varn,:)))/2
	1334	res@trYMinF = min(data(varn,:))+(max(data(varn,:)))/2
[161]	1335	end if
	1336	end if
	1337	plot_ps(n) = gsn_csm_xy(wks_ps,t,data(varn,:),res)
	1338
	1339	end if
	1340	end if
[154]	1341	end do
[161]	1342
[154]	1343	; ***************************************************
	1344	; merge plots onto one page
	1345	; ***************************************************
[532]	1346
	1347	no_frames = 0
[162]	1348
[566]	1349	if ((format_out .EQ. "eps" .OR. format_out .EQ. "epsi") .AND. \
	1350	n .gt. no_rows*no_columns) then
[161]	1351	gsn_panel(wks_ps,plot_ps(1:n),(/n,1/),resP)
[250]	1352	print(" ")
	1353	print("Outputs to .eps or .epsi have only one frame")
	1354	print(" ")
[161]	1355	else
[218]	1356	do np = 1,n,no_rows*no_columns
	1357	if ( np + no_rows*no_columns .gt. n) then
	1358	gsn_panel(wks_ps, plot_ps(np:n),(/no_rows,no_columns/),resP)
[532]	1359	no_frames = no_frames + 1
[161]	1360	else
[566]	1361	gsn_panel(wks_ps, plot_ps(np:np+no_rows*no_columns-1),\
	1362	(/no_rows,no_columns/),resP)
[532]	1363	no_frames = no_frames + 1
[161]	1364	end if
	1365	end do
	1366	end if
[154]	1367
[532]	1368	if (format_out .EQ. "png" ) then
	1369	png_output = new((/no_frames/), string)
	1370	j = 0
	1371
[969]	1372	if (no_frames .eq. 1) then
	1373	if (ncl_version .GE. 6 ) then
	1374	png_output(0) = file_out+".png"
	1375	else
	1376	png_output(0) = file_out+".00000"+1+".png"
	1377	end if
	1378	;using imagemagick's convert for reducing the white
	1379	;space around the plot
	1380	cmd = "convert -geometry 1000x1000 -density 300 -trim " + \
	1381	png_output(0) + " " + png_output(0)
[532]	1382	system(cmd)
[969]	1383	else
[532]	1384
[969]	1385	do i=0, no_frames-1
	1386	j = i + 1
	1387	if (j .LE. 9) then
	1388	png_output(i) = file_out+".00000"+j+".png"
	1389	end if
	1390	if (j .GT. 9 .AND. j .LE. 99) then
	1391	png_output(i) = file_out+".0000"+j+".png"
	1392	end if
	1393	if (j .GT. 99 .AND. j .LE. 999) then
	1394	png_output(i) = file_out+".000"+j+".png"
	1395	end if
	1396	if (j .GT. 999) then
	1397	png_output(i) = file_out+".00"+j+".png"
	1398	end if
	1399
	1400	;using imagemagick's convert for reducing the white
	1401	;space around the plot
	1402	cmd = "convert -geometry 1000x1000 -density 300 -trim " + \
	1403	png_output(i) + " " + png_output(i)
	1404	system(cmd)
	1405	end do
	1406	end if
	1407
[532]	1408	print(" ")
	1409	print("Output to: "+ png_output)
	1410	print(" ")
	1411	else
	1412	print(" ")
	1413	print("Output to: " + file_out +"."+ format_out)
	1414	print(" ")
	1415	end if
[969]	1416
[154]	1417	end

Note: See TracBrowser for help on using the repository browser.

Download in other formats:

Original Format