1 | load "$NCARG_ROOT/lib/ncarg/nclex/gsun/gsn_code.ncl" |
---|
2 | load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl" |
---|
3 | load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl" |
---|
4 | load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl" |
---|
5 | load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/shea_util.ncl" |
---|
6 | |
---|
7 | ;*************************************************** |
---|
8 | ; load ncl_preferences.ncl |
---|
9 | ;*************************************************** |
---|
10 | |
---|
11 | function which_script() |
---|
12 | local script |
---|
13 | begin |
---|
14 | script="cross_section" |
---|
15 | return(script) |
---|
16 | end |
---|
17 | |
---|
18 | if (isfilepresent("~/ncl_preferences.ncl")) then |
---|
19 | loadscript("~/ncl_preferences.ncl") |
---|
20 | else |
---|
21 | if (isfilepresent("~/palm/current_version/trunk/SCRIPTS/NCL/ncl_preferences.ncl")) then |
---|
22 | loadscript( "~/palm/current_version/trunk/SCRIPTS/NCL/ncl_preferences.ncl") |
---|
23 | else |
---|
24 | print(" ") |
---|
25 | print("'ncl_preferences.ncl' does not exist in $home or $home/palm/current_version/trunk/SCRIPTS/NCL/") |
---|
26 | print(" ") |
---|
27 | exit |
---|
28 | end if |
---|
29 | end if |
---|
30 | |
---|
31 | begin |
---|
32 | |
---|
33 | ; *************************************************** |
---|
34 | ; set default parameter values and strings if not assigned in prompt or parameter list |
---|
35 | ; *************************************************** |
---|
36 | |
---|
37 | if (file_1 .EQ. "File in") then |
---|
38 | print(" ") |
---|
39 | print("Declare input file 'file_1=' in 'ncl_preferences.ncl' or prompt") |
---|
40 | print(" ") |
---|
41 | exit |
---|
42 | else |
---|
43 | file_in = file_1 |
---|
44 | end if |
---|
45 | |
---|
46 | if (format_out .NE. "x11" .AND. format_out .NE. "pdf" .AND. format_out .NE. "eps" .AND. format_out .NE. "ps" .AND. format_out .NE. "epsi" .AND. format_out .NE. "ncgm")then |
---|
47 | print(" ") |
---|
48 | print("'format_out = "+format_out+"' is invalid and set to'x11'") |
---|
49 | print(" ") |
---|
50 | format_out="x11" |
---|
51 | end if |
---|
52 | |
---|
53 | if (sort .NE. "layer" .AND. sort .NE. "time") then |
---|
54 | print(" ") |
---|
55 | print("'sort'= "+sort+" is invalid and set to 'layer'") |
---|
56 | print(" ") |
---|
57 | sort = "layer" |
---|
58 | end if |
---|
59 | |
---|
60 | if (mode .NE. "Fill" .AND. mode .NE. "Line" .AND. mode .NE. "Both") then |
---|
61 | print(" ") |
---|
62 | print("'mode'= "+mode+" is invalid and set to 'Fill'") |
---|
63 | print(" ") |
---|
64 | mode = "Fill" |
---|
65 | end if |
---|
66 | |
---|
67 | if (fill_mode .NE. "AreaFill" .AND. fill_mode .NE. "RasterFill" .AND. fill_mode .NE. "CellFill") then |
---|
68 | print(" ") |
---|
69 | print("'fill_mode'= "+fill_mode+" is invalid and set to 'AreaFill'") |
---|
70 | print(" ") |
---|
71 | fill_mode = "AreaFill" |
---|
72 | end if |
---|
73 | |
---|
74 | if (shape .NE. 0 .AND. shape .NE. 1) then |
---|
75 | print(" ") |
---|
76 | print("'shape'= "+shape+" is invalid and set to 1") |
---|
77 | print(" ") |
---|
78 | shape = 1 |
---|
79 | end if |
---|
80 | |
---|
81 | if (xyc .NE. 0 .AND. xyc .NE. 1)then |
---|
82 | print(" ") |
---|
83 | print("'xyc'= "+xyc+" is invalid and set to 0") |
---|
84 | print(" ") |
---|
85 | xyc = 0 |
---|
86 | end if |
---|
87 | |
---|
88 | if (xzc .NE. 0 .AND. xzc .NE. 1)then |
---|
89 | print(" ") |
---|
90 | print("'xzc'= "+xzc+" is invalid and set to 0") |
---|
91 | print(" ") |
---|
92 | xyc = 0 |
---|
93 | end if |
---|
94 | |
---|
95 | if (yzc .NE. 0 .AND. yzc .NE. 1)then |
---|
96 | print(" ") |
---|
97 | print("'yzc'= "+yzc+" is invalid and set to 0") |
---|
98 | print(" ") |
---|
99 | xyc = 0 |
---|
100 | end if |
---|
101 | |
---|
102 | if (xyc .EQ. 0 .AND. xzc .EQ. 0 .AND. yzc .EQ. 0) then |
---|
103 | print(" ") |
---|
104 | print("Select one crossection (xyc=1 or xzc=1 or yzc=1)") |
---|
105 | print(" ") |
---|
106 | exit |
---|
107 | end if |
---|
108 | if (xyc .EQ. 1 ) then |
---|
109 | if (xzc .EQ. 1 .OR. yzc .EQ. 1) then |
---|
110 | print(" ") |
---|
111 | print("Select just one crossection (xyc=1 or xzc=1 or yzc=1)") |
---|
112 | print(" ") |
---|
113 | exit |
---|
114 | end if |
---|
115 | end if |
---|
116 | if (xzc .EQ. 1) then |
---|
117 | if (xyc .EQ. 1 .OR. yzc .EQ. 1) then |
---|
118 | print(" ") |
---|
119 | print("Select just one crossection (xyc=1 or xzc=1 or yzc=1)") |
---|
120 | print(" ") |
---|
121 | exit |
---|
122 | end if |
---|
123 | end if |
---|
124 | if (yzc .EQ. 1) then |
---|
125 | if (xyc .EQ. 1 .OR. xzc .EQ. 1) then |
---|
126 | print(" ") |
---|
127 | print("Select just one crossection (xyc=1 or xzc=1 or yzc=1)") |
---|
128 | print(" ") |
---|
129 | exit |
---|
130 | end if |
---|
131 | end if |
---|
132 | |
---|
133 | if (vector .NE. 0 .AND. vector .NE. 1) then |
---|
134 | print(" ") |
---|
135 | print("Set 'vector' to 0 or 1") |
---|
136 | print(" ") |
---|
137 | exit |
---|
138 | end if |
---|
139 | |
---|
140 | if (norm_x .EQ. 0) then |
---|
141 | print(" ") |
---|
142 | print("Normalising with 0 is not allowed, 'norm_x' is set to 1.0") |
---|
143 | print(" ") |
---|
144 | norm_x = 1.0 |
---|
145 | end if |
---|
146 | if (norm_y .EQ. 0) then |
---|
147 | print(" ") |
---|
148 | print("Normalising with 0 is not allowed, 'norm_y' is set to 1.0") |
---|
149 | print(" ") |
---|
150 | norm_y = 1.0 |
---|
151 | end if |
---|
152 | if (norm_z .EQ. 0) then |
---|
153 | print(" ") |
---|
154 | print("Normalising with 0 is not allowed, 'norm_z' is set to 1.0") |
---|
155 | print(" ") |
---|
156 | norm_z = 1.0 |
---|
157 | end if |
---|
158 | |
---|
159 | ; *************************************************** |
---|
160 | ; open input file |
---|
161 | ; *************************************************** |
---|
162 | |
---|
163 | file_in_1 = False |
---|
164 | if (isStrSubset(file_in, ".nc"))then |
---|
165 | start_f = -2 |
---|
166 | end_f = -2 |
---|
167 | file_in_1 = True |
---|
168 | end if |
---|
169 | |
---|
170 | if (start_f .EQ. -1)then |
---|
171 | print(" ") |
---|
172 | print("'start_f' must be one of the cyclic numbers (at least 0) of your input file(s)") |
---|
173 | print(" ") |
---|
174 | exit |
---|
175 | end if |
---|
176 | if (end_f .EQ. -1)then |
---|
177 | print(" ") |
---|
178 | print("'end_f' must be one of the cyclic numbers (at least 0) of your input file(s)") |
---|
179 | print(" ") |
---|
180 | exit |
---|
181 | end if |
---|
182 | |
---|
183 | files=new(end_f-start_f+1,string) |
---|
184 | |
---|
185 | if (file_in_1)then |
---|
186 | if (isfilepresent(file_in))then |
---|
187 | files(0)=file_in |
---|
188 | else |
---|
189 | print(" ") |
---|
190 | print("1st input file: '"+file_in+"' does not exist") |
---|
191 | print(" ") |
---|
192 | exit |
---|
193 | end if |
---|
194 | else |
---|
195 | if (start_f .EQ. 0)then |
---|
196 | if (isfilepresent(file_in+".nc"))then |
---|
197 | files(0)=file_in+".nc" |
---|
198 | do i=1,end_f |
---|
199 | if (isfilepresent(file_in+"."+i+".nc"))then |
---|
200 | files(i)=file_in+"."+i+".nc" |
---|
201 | else |
---|
202 | print(" ") |
---|
203 | print("Input file: '"+file_in+"."+i+".nc' does not exist") |
---|
204 | print(" ") |
---|
205 | exit |
---|
206 | end if |
---|
207 | end do |
---|
208 | else |
---|
209 | print(" ") |
---|
210 | print("Input file: '"+file_in+".nc' does not exist") |
---|
211 | print(" ") |
---|
212 | exit |
---|
213 | end if |
---|
214 | else |
---|
215 | do i=start_f,end_f |
---|
216 | if (isfilepresent(file_in+"."+i+".nc"))then |
---|
217 | files(i-start_f)=file_in+"."+i+".nc" |
---|
218 | else |
---|
219 | print(" ") |
---|
220 | print("Input file: '"+file_in+"."+i+".nc' does not exist") |
---|
221 | print(" ") |
---|
222 | exit |
---|
223 | end if |
---|
224 | end do |
---|
225 | end if |
---|
226 | end if |
---|
227 | |
---|
228 | f=addfiles(files,"r") |
---|
229 | f_att=addfile(files(0),"r") |
---|
230 | ListSetType(f,"cat") |
---|
231 | |
---|
232 | vNam = getfilevarnames(f_att) |
---|
233 | print(" ") |
---|
234 | print("Variables in input file:") |
---|
235 | print("- "+ vNam) |
---|
236 | print(" ") |
---|
237 | dim = dimsizes(vNam) |
---|
238 | |
---|
239 | ; *************************************************** |
---|
240 | ; check for kind of input file |
---|
241 | ; *************************************************** |
---|
242 | |
---|
243 | check_3d = True |
---|
244 | do varn=0,dim-1 |
---|
245 | checkxy = isStrSubset( vNam(varn),"xy") |
---|
246 | checkxz = isStrSubset( vNam(varn),"xz") |
---|
247 | checkyz = isStrSubset( vNam(varn),"yz") |
---|
248 | if (checkxy .OR. checkxz .OR. checkyz) then |
---|
249 | check_3d = False |
---|
250 | break |
---|
251 | end if |
---|
252 | end do |
---|
253 | if (.not. check_3d) then |
---|
254 | if (xyc .EQ. 1)then |
---|
255 | if (.not. checkxy) then |
---|
256 | print(" ") |
---|
257 | print("Your input file doesn't have values for xy cross-sections") |
---|
258 | if (checkxz)then |
---|
259 | print("Select another input file or xz cross-sections") |
---|
260 | else |
---|
261 | print("Select another input file or yz cross-sections") |
---|
262 | end if |
---|
263 | print(" ") |
---|
264 | exit |
---|
265 | else |
---|
266 | print(" ") |
---|
267 | print("Your input file contains xy data") |
---|
268 | print(" ") |
---|
269 | end if |
---|
270 | end if |
---|
271 | if (xzc .EQ. 1)then |
---|
272 | if (.not. checkxz) then |
---|
273 | print(" ") |
---|
274 | print("Your input file doesn't have values for xz cross-sections") |
---|
275 | if (checkxy)then |
---|
276 | print("Select another input file or xy cross-sections") |
---|
277 | else |
---|
278 | print("Select another input file or yz cross-sections") |
---|
279 | end if |
---|
280 | print(" ") |
---|
281 | exit |
---|
282 | else |
---|
283 | print(" ") |
---|
284 | print("Your input file contains xz data") |
---|
285 | print(" ") |
---|
286 | end if |
---|
287 | end if |
---|
288 | if (yzc .EQ. 1)then |
---|
289 | if (.not. checkyz) then |
---|
290 | print(" ") |
---|
291 | print("Your input file doesn't have values for yz cross-sections") |
---|
292 | if (checkxy)then |
---|
293 | print("Select another input file or xy cross-sections") |
---|
294 | else |
---|
295 | print("Select another input file or xz cross-sections") |
---|
296 | end if |
---|
297 | print(" ") |
---|
298 | exit |
---|
299 | else |
---|
300 | print(" ") |
---|
301 | print("Your input file contains yz data") |
---|
302 | print(" ") |
---|
303 | end if |
---|
304 | end if |
---|
305 | else |
---|
306 | print(" ") |
---|
307 | print("Your input file contains 3d or other data") |
---|
308 | print(" ") |
---|
309 | end if |
---|
310 | |
---|
311 | if (dim .EQ. 0) then |
---|
312 | print(" ") |
---|
313 | print("There is no data on file") |
---|
314 | print(" ") |
---|
315 | exit |
---|
316 | end if |
---|
317 | |
---|
318 | ; *************************************************** |
---|
319 | ; set up recourses |
---|
320 | ; *************************************************** |
---|
321 | |
---|
322 | cs_res = True |
---|
323 | vecres = True |
---|
324 | cs_res@gsnYAxisIrregular2Linear = True |
---|
325 | |
---|
326 | cs_res@gsnDraw = False |
---|
327 | cs_res@gsnFrame = False |
---|
328 | cs_res@gsnMaximize = True |
---|
329 | |
---|
330 | cs_res@tmXBLabelFontHeightF = font_size |
---|
331 | cs_res@tmYLLabelFontHeightF = font_size |
---|
332 | cs_res@tiXAxisFontHeightF = font_size |
---|
333 | cs_res@tiYAxisFontHeightF = font_size |
---|
334 | cs_res@tmXBMode ="Automatic" |
---|
335 | cs_res@tmYLMode ="Automatic" |
---|
336 | |
---|
337 | cs_res@cnLevelSelectionMode = "ManualLevels" |
---|
338 | cs_res@lbLabelFontHeightF = font_size_legend |
---|
339 | cs_res@lbLabelStride = legend_label_stride |
---|
340 | |
---|
341 | |
---|
342 | cs_resP = True |
---|
343 | cs_resP@txString = f_att@title |
---|
344 | cs_resP@txFuncCode = "~" |
---|
345 | cs_resP@txFontHeightF = 0.015 |
---|
346 | |
---|
347 | if ( mode .eq. "Fill" ) then |
---|
348 | cs_res@cnFillOn = True |
---|
349 | cs_res@gsnSpreadColors = True |
---|
350 | cs_res@cnFillMode = fill_mode |
---|
351 | cs_res@cnLinesOn = False |
---|
352 | cs_res@cnLineLabelsOn = False |
---|
353 | end if |
---|
354 | |
---|
355 | if ( mode .eq. "Both" ) then |
---|
356 | cs_res@cnFillOn = True |
---|
357 | cs_res@gsnSpreadColors = True |
---|
358 | cs_res@cnFillMode = fill_mode |
---|
359 | cs_res@cnLinesOn = True |
---|
360 | cs_res@cnLineLabelsOn = True |
---|
361 | end if |
---|
362 | |
---|
363 | ; ********************************************* |
---|
364 | ; set up of start_time_step and end_time_step |
---|
365 | ; ********************************************* |
---|
366 | |
---|
367 | t_all = f[:]->time |
---|
368 | nt = dimsizes(t_all) |
---|
369 | delta_t = t_all(nt-1)/nt |
---|
370 | |
---|
371 | ; **************************************************** |
---|
372 | ; start of time step and different types of mistakes that could be done |
---|
373 | ; **************************************************** |
---|
374 | |
---|
375 | if (start_time_step .EQ. -1.) then |
---|
376 | start_time_step=t_all(0)/3600 |
---|
377 | else |
---|
378 | if (start_time_step .GT. t_all(nt-1)/3600) |
---|
379 | print(" ") |
---|
380 | print("'start_time_step' = "+ start_time_step +"h is greater than last time step = " + t_all(nt-1)+"s = "+t_all(nt-1)/3600+"h") |
---|
381 | print(" ") |
---|
382 | print("Select another 'start_time_step'") |
---|
383 | print(" ") |
---|
384 | exit |
---|
385 | end if |
---|
386 | if (start_time_step .LT. t_all(0)/3600) |
---|
387 | print(" ") |
---|
388 | print("'start_time_step' = "+ start_time_step +"h is lower than first time step = " + t_all(0)+"s = "+t_all(0)/3600+"h") |
---|
389 | print(" ") |
---|
390 | print("Select another 'start_time_step'") |
---|
391 | print(" ") |
---|
392 | exit |
---|
393 | end if |
---|
394 | end if |
---|
395 | do i=0,nt-1 |
---|
396 | if (start_time_step .GE. (t_all(i)-delta_t/2)/3600 .AND. start_time_step .LT. (t_all(i)+delta_t/2)/3600)then |
---|
397 | st=i |
---|
398 | break |
---|
399 | end if |
---|
400 | end do |
---|
401 | |
---|
402 | if (.not. isvar("st"))then |
---|
403 | print(" ") |
---|
404 | print("'start_time_step' = "+ start_time_step +"h is invalid") |
---|
405 | print(" ") |
---|
406 | print("Select another 'start_time_step'") |
---|
407 | print(" ") |
---|
408 | exit |
---|
409 | end if |
---|
410 | |
---|
411 | ; **************************************************** |
---|
412 | ; end of time step and different types of mistakes that could be done |
---|
413 | ; **************************************************** |
---|
414 | |
---|
415 | if (end_time_step .EQ. -1.) then |
---|
416 | end_time_step = t_all(nt-1)/3600 |
---|
417 | else |
---|
418 | if (end_time_step .LT. t_all(0)/3600) |
---|
419 | print(" ") |
---|
420 | print("'end_time_step' = "+end_time_step+ "h is lower than first time step = " + t_all(0)+"s = "+t_all(0)/3600+"h") |
---|
421 | print(" ") |
---|
422 | print("Select another 'end_time_step'") |
---|
423 | print(" ") |
---|
424 | exit |
---|
425 | end if |
---|
426 | if (end_time_step .GT. t_all(nt-1)/3600) |
---|
427 | print(" ") |
---|
428 | print("'end_time_step' = "+ end_time_step +"h is greater than last time step = " + t_all(nt-1)+"s = "+t_all(nt-1)/3600+"h") |
---|
429 | print(" ") |
---|
430 | print("Select another 'end_time_step'") |
---|
431 | print(" ") |
---|
432 | exit |
---|
433 | end if |
---|
434 | if (end_time_step .LT. start_time_step/3600) |
---|
435 | print(" ") |
---|
436 | print("'end_time_step' = "+ end_time_step +"h is lower than 'start_time_step' = "+start_time_step+"h") |
---|
437 | print(" ") |
---|
438 | print("Select another 'start_time_step' or 'end_time_step'") |
---|
439 | print(" ") |
---|
440 | exit |
---|
441 | end if |
---|
442 | end if |
---|
443 | do i=0,nt-1 |
---|
444 | if (end_time_step .GE. (t_all(i)-delta_t/2)/3600 .AND. end_time_step .LT. (t_all(i)+delta_t/2)/3600)then |
---|
445 | et=i |
---|
446 | break |
---|
447 | end if |
---|
448 | end do |
---|
449 | |
---|
450 | if (.not. isvar("et"))then |
---|
451 | print(" ") |
---|
452 | print("'end_time_step' = "+ end_time_step +"h is invalid") |
---|
453 | print(" ") |
---|
454 | print("Select another 'end_time_step'") |
---|
455 | print(" ") |
---|
456 | exit |
---|
457 | end if |
---|
458 | |
---|
459 | delete(start_time_step) |
---|
460 | start_time_step=round(st,3) |
---|
461 | delete(end_time_step) |
---|
462 | end_time_step=round(et,3) |
---|
463 | |
---|
464 | print(" ") |
---|
465 | print("Output of time steps from "+t_all(start_time_step)/3600+" h = "+t_all(start_time_step)+" s => index = "+start_time_step) |
---|
466 | print(" till "+t_all(end_time_step)/3600+" h = "+t_all(end_time_step)+" s => index = "+end_time_step) |
---|
467 | print(" ") |
---|
468 | |
---|
469 | no_time=(end_time_step-start_time_step)+1 |
---|
470 | |
---|
471 | ; **************************************************** |
---|
472 | ; set up variables for vector plot if required |
---|
473 | ; **************************************************** |
---|
474 | |
---|
475 | if (vector .EQ. 1) then |
---|
476 | if (vec1 .EQ. "component1") then |
---|
477 | print(" ") |
---|
478 | print("Indicate Vector 1 ('vec1') for Vector-Plot or set 'vector' to 0") |
---|
479 | print(" ") |
---|
480 | exit |
---|
481 | end if |
---|
482 | if (vec2 .EQ. "component2") then |
---|
483 | print(" ") |
---|
484 | print("Indicate Vector 2 ('vec2') for Vector-Plot or set 'vector' to 0") |
---|
485 | print(" ") |
---|
486 | exit |
---|
487 | end if |
---|
488 | end if |
---|
489 | |
---|
490 | check_vec1 = False |
---|
491 | check_vec2 = False |
---|
492 | check_vecp = False |
---|
493 | |
---|
494 | ; **************************************************** |
---|
495 | ; get data for plots |
---|
496 | ; **************************************************** |
---|
497 | |
---|
498 | if (xyc .EQ. 1) then |
---|
499 | do varn=0,dim-1 |
---|
500 | if (vNam(varn) .eq. "xu" .OR. vNam(varn) .eq. "x")then |
---|
501 | x_d = f_att->$vNam(varn)$ |
---|
502 | xdim = dimsizes(x_d)-1 |
---|
503 | delta_x = x_d(1)-x_d(0) |
---|
504 | break |
---|
505 | end if |
---|
506 | end do |
---|
507 | do varn=0,dim-1 |
---|
508 | if (vNam(varn) .eq. "yv" .or. vNam(varn) .eq. "y")then |
---|
509 | y_d = f_att->$vNam(varn)$ |
---|
510 | ydim = dimsizes(y_d)-1 |
---|
511 | delta_y = y_d(1)-y_d(0) |
---|
512 | break |
---|
513 | end if |
---|
514 | end do |
---|
515 | do varn=0,dim-1 |
---|
516 | if (vNam(varn) .eq. "zu_3d" .OR. vNam(varn) .eq. "zw_3d")then |
---|
517 | z_d = f_att->$vNam(varn)$ |
---|
518 | zdim = dimsizes(z_d)-1 |
---|
519 | delta_z = 0 |
---|
520 | break |
---|
521 | else |
---|
522 | if (vNam(varn) .eq. "zu_xy" .OR. vNam(varn) .eq. "zw_xy") then |
---|
523 | z_d = f_att->$vNam(varn)$ |
---|
524 | zdim = dimsizes(z_d)-1 |
---|
525 | delta_z = -1.d |
---|
526 | break |
---|
527 | else |
---|
528 | zdim = 0 |
---|
529 | delta_z = -1.d |
---|
530 | end if |
---|
531 | end if |
---|
532 | end do |
---|
533 | end if |
---|
534 | if (xzc .EQ. 1) then |
---|
535 | do varn=0,dim-1 |
---|
536 | if (vNam(varn) .eq. "xu" .OR. vNam(varn) .eq. "x") then |
---|
537 | x_d = f_att->$vNam(varn)$ |
---|
538 | xdim = dimsizes(x_d)-1 |
---|
539 | delta_x = x_d(1)-x_d(0) |
---|
540 | break |
---|
541 | end if |
---|
542 | end do |
---|
543 | do varn=0,dim-1 |
---|
544 | if (vNam(varn) .eq. "zw" .or. vNam(varn) .eq. "zu" .OR. vNam(varn) .eq. "zu_3d" .or. vNam(varn) .eq. "zw_3d") then |
---|
545 | z_d = f_att->$vNam(varn)$ |
---|
546 | zdim = dimsizes(z_d)-1 |
---|
547 | delta_z = z_d(1)-z_d(0) |
---|
548 | break |
---|
549 | end if |
---|
550 | end do |
---|
551 | do varn=0,dim-1 |
---|
552 | if (vNam(varn) .eq. "yv" .or. vNam(varn) .eq. "y") then |
---|
553 | y_d = f_att->$vNam(varn)$ |
---|
554 | ydim = dimsizes(y_d)-1 |
---|
555 | delta_y = y_d(1)-y_d(0) |
---|
556 | break |
---|
557 | else |
---|
558 | if (vNam(varn) .eq. "y_xz" .or. vNam(varn) .eq. "yv_xz") then |
---|
559 | y_d = f_att->$vNam(varn)$ |
---|
560 | ydim = dimsizes(y_d)-1 |
---|
561 | delta_y = -1.d |
---|
562 | break |
---|
563 | else |
---|
564 | ydim = 0 |
---|
565 | delta_y = -1.d |
---|
566 | end if |
---|
567 | end if |
---|
568 | end do |
---|
569 | end if |
---|
570 | if (yzc .EQ. 1) then |
---|
571 | do varn=0,dim-1 |
---|
572 | if (vNam(varn) .eq. "yv" .or. vNam(varn) .eq. "y") then |
---|
573 | y_d = f_att->$vNam(varn)$ |
---|
574 | ydim = dimsizes(y_d)-1 |
---|
575 | delta_y = y_d(1)-y_d(0) |
---|
576 | break |
---|
577 | end if |
---|
578 | end do |
---|
579 | do varn=0,dim-1 |
---|
580 | if (vNam(varn) .eq. "zw" .or. vNam(varn) .eq. "zu" .OR. vNam(varn) .eq. "zu_3d" .or. vNam(varn) .eq. "zw_3d") then |
---|
581 | z_d = f_att->$vNam(varn)$ |
---|
582 | zdim = dimsizes(z_d)-1 |
---|
583 | delta_z = z_d(1)-z_d(0) |
---|
584 | break |
---|
585 | end if |
---|
586 | end do |
---|
587 | do varn=0,dim-1 |
---|
588 | if (vNam(varn) .eq. "xu" .or. vNam(varn) .eq. "x") |
---|
589 | x_d = f_att->$vNam(varn)$ |
---|
590 | xdim = dimsizes(x_d)-1 |
---|
591 | delta_x = x_d(1)-x_d(0) |
---|
592 | break |
---|
593 | else |
---|
594 | if (vNam(varn) .eq. "xu_yz" .or. vNam(varn) .eq. "x_yz" ) then |
---|
595 | x_d = f_att->$vNam(varn)$ |
---|
596 | xdim = dimsizes(x_d)-1 |
---|
597 | delta_x = -1.d |
---|
598 | break |
---|
599 | else |
---|
600 | xdim = 0 |
---|
601 | delta_x = -1.d |
---|
602 | end if |
---|
603 | end if |
---|
604 | end do |
---|
605 | end if |
---|
606 | |
---|
607 | ; **************************************************** |
---|
608 | ; set up ranges of x-, y- and z-coordinates |
---|
609 | ; **************************************************** |
---|
610 | |
---|
611 | if (xs .EQ. -1.d) then |
---|
612 | xs = x_d(0) |
---|
613 | if (delta_x .EQ. -1) then |
---|
614 | print(" ") |
---|
615 | print("You cannot choose a start value for x, there are preseted layers for x") |
---|
616 | print(" ") |
---|
617 | xstart=0 |
---|
618 | end if |
---|
619 | else |
---|
620 | if (delta_x .EQ. -1) then |
---|
621 | print(" ") |
---|
622 | print("You cannot choose a start value for x, there are preseted layers for x") |
---|
623 | print(" ") |
---|
624 | xstart=0 |
---|
625 | else |
---|
626 | if (xs .LT. 0-delta_x/2) then |
---|
627 | print(" ") |
---|
628 | print("range start for x-coordinate = "+xs+"m is lower than first value = "+(0-delta_x/2)+"m") |
---|
629 | print(" ") |
---|
630 | exit |
---|
631 | end if |
---|
632 | if (xyc .EQ. 1 .OR. xzc .EQ. 1) then |
---|
633 | if (xs .GE. xdim*delta_x) then |
---|
634 | print(" ") |
---|
635 | print("range start for x-coordinate = "+xs+"m is equal or greater than last value = "+xdim*delta_x+"m") |
---|
636 | print(" ") |
---|
637 | exit |
---|
638 | end if |
---|
639 | else |
---|
640 | if (xs .GT. xdim*delta_x) then |
---|
641 | print(" ") |
---|
642 | print("range start for x-coordinate = "+xs+"m is greater than last value = "+xdim*delta_x+"m") |
---|
643 | print(" ") |
---|
644 | exit |
---|
645 | end if |
---|
646 | end if |
---|
647 | end if |
---|
648 | end if |
---|
649 | |
---|
650 | if (ys .EQ. -1.d) then |
---|
651 | ys = y_d(0) |
---|
652 | if (delta_y .EQ. -1) then |
---|
653 | print(" ") |
---|
654 | print("You cannot choose a start value for y, there are preseted layers for y") |
---|
655 | print(" ") |
---|
656 | ystart=0 |
---|
657 | end if |
---|
658 | else |
---|
659 | if (delta_y .EQ. -1) then |
---|
660 | print(" ") |
---|
661 | print("You cannot choose a start value for y, there are preseted layers for y") |
---|
662 | print(" ") |
---|
663 | ystart=0 |
---|
664 | else |
---|
665 | if (ys .LT. 0-delta_y/2) then |
---|
666 | print(" ") |
---|
667 | print("range start for y-coordinate = "+ys+"m is lower than first value = "+0-delta_y/2+"m") |
---|
668 | print(" ") |
---|
669 | exit |
---|
670 | end if |
---|
671 | if (xyc .EQ. 1 .OR. yzc .EQ. 1) then |
---|
672 | if (ys .GE. ydim*delta_y) then |
---|
673 | print(" ") |
---|
674 | print("range start for y-coordinate = "+ys+"m is equal or greater than last value = "+ydim*delta_y+"m") |
---|
675 | print(" ") |
---|
676 | exit |
---|
677 | end if |
---|
678 | else |
---|
679 | if (ys .GT. ydim*delta_y) then |
---|
680 | print(" ") |
---|
681 | print("range start for y-coordinate = "+ys+"m is greater than last value = "+ydim*delta_y+"m") |
---|
682 | print(" ") |
---|
683 | exit |
---|
684 | end if |
---|
685 | end if |
---|
686 | end if |
---|
687 | end if |
---|
688 | |
---|
689 | if (zs .EQ. -1) then |
---|
690 | zs = 0 |
---|
691 | if (delta_z .EQ. -1) then |
---|
692 | print(" ") |
---|
693 | print("You cannot choose a start value for z, there are preseted layers for z") |
---|
694 | print(" ") |
---|
695 | end if |
---|
696 | else |
---|
697 | if (delta_z .EQ. -1) then |
---|
698 | print(" ") |
---|
699 | print("You cannot choose a start value for z, there are preseted layers for z") |
---|
700 | print(" ") |
---|
701 | else |
---|
702 | if (zs .LT. 0) then |
---|
703 | print(" ") |
---|
704 | print("range start for z-coordinate = "+zs+" is lower than first gridpoint = 0") |
---|
705 | print(" ") |
---|
706 | exit |
---|
707 | end if |
---|
708 | if (xzc .EQ. 1 .OR. yzc .EQ. 1) then |
---|
709 | if (zs .GE. zdim) then |
---|
710 | print(" ") |
---|
711 | print("range start for z-coordinate = "+zs+" is equal or greater than last gridpoint = "+zdim) |
---|
712 | print(" ") |
---|
713 | exit |
---|
714 | end if |
---|
715 | else |
---|
716 | if (zs .GT. zdim) then |
---|
717 | print(" ") |
---|
718 | print("range start for z-coordinate = "+zs+" is greater than last gridpoint = "+zdim) |
---|
719 | print(" ") |
---|
720 | exit |
---|
721 | end if |
---|
722 | end if |
---|
723 | end if |
---|
724 | end if |
---|
725 | |
---|
726 | if (xe .EQ. -1) then |
---|
727 | xe = x_d(xdim) |
---|
728 | if (delta_x .EQ. -1) then |
---|
729 | print(" ") |
---|
730 | print("You cannot choose an end value for x, there are preseted layers for x") |
---|
731 | print(" ") |
---|
732 | xend=xdim |
---|
733 | end if |
---|
734 | else |
---|
735 | if (delta_x .EQ. -1) then |
---|
736 | print(" ") |
---|
737 | print("You cannot choose an end value for x, there are preseted layers for x") |
---|
738 | print(" ") |
---|
739 | xend=xdim |
---|
740 | else |
---|
741 | if (xe .GT. xdim*delta_x) then |
---|
742 | print(" ") |
---|
743 | print("range end for x-coordinate = "+xe+"m is greater than last value = "+xdim*delta_x+"m") |
---|
744 | print(" ") |
---|
745 | exit |
---|
746 | end if |
---|
747 | if (xyc .EQ. 1 .OR. xzc .EQ. 1) then |
---|
748 | if (xe .LE. 0-delta_x/2) |
---|
749 | print(" ") |
---|
750 | print("range end for x-coordinate = "+xe+"m is equal or lower than first value = "+(0-delta_x/2)+"m") |
---|
751 | print(" ") |
---|
752 | exit |
---|
753 | end if |
---|
754 | if (xe .LE. xs) then |
---|
755 | print(" ") |
---|
756 | print("range end for x-coordinate = "+xe+"m is equal or lower than start range = "+xs+"m") |
---|
757 | print(" ") |
---|
758 | exit |
---|
759 | end if |
---|
760 | if (xe .EQ. xs+1) then |
---|
761 | print(" ") |
---|
762 | print("range end for x-coordinate = "+xe+"m must be at least two more gridpoints greater than start range = "+xs+"m") |
---|
763 | print(" ") |
---|
764 | exit |
---|
765 | end if |
---|
766 | else |
---|
767 | if (xe .LT. 0-delta_x/2) |
---|
768 | print(" ") |
---|
769 | print("range end for x-coordinate = "+xe+"m is lower than first value = "+(0-delta_x/2)+"m") |
---|
770 | print(" ") |
---|
771 | exit |
---|
772 | end if |
---|
773 | if (xe .LT. xs) then |
---|
774 | print(" ") |
---|
775 | print("range end for x-coordinate = "+xe+"m is lower than start range = "+xs+"m") |
---|
776 | print(" ") |
---|
777 | exit |
---|
778 | end if |
---|
779 | end if |
---|
780 | end if |
---|
781 | end if |
---|
782 | |
---|
783 | if (ye .EQ. -1) then |
---|
784 | ye = y_d(ydim) |
---|
785 | if (delta_y .EQ. -1) then |
---|
786 | print(" ") |
---|
787 | print("You cannot choose an end value for y, there are preseted layers for y") |
---|
788 | print(" ") |
---|
789 | yend=ydim |
---|
790 | end if |
---|
791 | else |
---|
792 | if (delta_y .EQ. -1) then |
---|
793 | print(" ") |
---|
794 | print("You cannot choose an end value for y, there are preseted layers for y") |
---|
795 | print(" ") |
---|
796 | yend=ydim |
---|
797 | else |
---|
798 | if (ye .GT. ydim*delta_y) then |
---|
799 | print(" ") |
---|
800 | print("range end for y-coordinate = "+ye+"m is greater than last value = "+ydim*delta_y+"m") |
---|
801 | print(" ") |
---|
802 | exit |
---|
803 | end if |
---|
804 | if (xyc .EQ. 1 .OR. yzc .EQ. 1) then |
---|
805 | if (ye .LE. 0-delta_y/2) |
---|
806 | print(" ") |
---|
807 | print("range end for y-coordinate = "+ye+"m is equal or lower than first value = "+(0-delta_y/2)+"m") |
---|
808 | print(" ") |
---|
809 | exit |
---|
810 | end if |
---|
811 | if (ye .LE. ys) then |
---|
812 | print(" ") |
---|
813 | print("range end for y-coordinate = "+ye+"m is equal or lower than start range = "+ys+"m") |
---|
814 | print(" ") |
---|
815 | exit |
---|
816 | end if |
---|
817 | if (ye .EQ. ys+1) then |
---|
818 | print(" ") |
---|
819 | print("range end for y-coordinate = "+ye+"m must be at least two more gridpoints greater than start range = "+ys+"m") |
---|
820 | print(" ") |
---|
821 | exit |
---|
822 | end if |
---|
823 | else |
---|
824 | if (ye .LT. 0-delta_y/2) |
---|
825 | print(" ") |
---|
826 | print("range end for y-coordinate = "+ye+"m is lower than first value = "+(0-delta_y/2)+"m") |
---|
827 | print(" ") |
---|
828 | exit |
---|
829 | end if |
---|
830 | if (ye .LT. ys) then |
---|
831 | print(" ") |
---|
832 | print("range end for y-coordinate = "+ye+"m is lower than start range = "+ys+"m") |
---|
833 | print(" ") |
---|
834 | exit |
---|
835 | end if |
---|
836 | end if |
---|
837 | end if |
---|
838 | end if |
---|
839 | |
---|
840 | if (ze .EQ. -1) then |
---|
841 | ze = zdim |
---|
842 | if (delta_z .EQ. -1) then |
---|
843 | print(" ") |
---|
844 | print("You cannot choose an end value for z, there are preseted layers for z") |
---|
845 | print(" ") |
---|
846 | ze = zdim |
---|
847 | end if |
---|
848 | else |
---|
849 | if (delta_z .EQ. -1) then |
---|
850 | print(" ") |
---|
851 | print("You cannot choose an end value for z, there are preseted layers for z") |
---|
852 | print(" ") |
---|
853 | ze = zdim |
---|
854 | else |
---|
855 | if (ze .GT. zdim) then |
---|
856 | print(" ") |
---|
857 | print("range end for z-coordinate = "+ze+" is greater than last gridpoint = "+zdim) |
---|
858 | print(" ") |
---|
859 | exit |
---|
860 | end if |
---|
861 | if (xzc .EQ. 1 .OR. yzc .EQ. 1) then |
---|
862 | if (ze .LE. 0) |
---|
863 | print(" ") |
---|
864 | print("range end for z-coordinate = "+ze+" is equal or lower than first gridpoint = 0") |
---|
865 | print(" ") |
---|
866 | exit |
---|
867 | end if |
---|
868 | if (ze .LE. zs) then |
---|
869 | print(" ") |
---|
870 | print("range end for z-coordinate = "+ze+" is equal or lower than start range = "+zs) |
---|
871 | print(" ") |
---|
872 | exit |
---|
873 | end if |
---|
874 | if (ze .EQ. zs+1) then |
---|
875 | print(" ") |
---|
876 | print("range end for z-coordinate = "+ze+" must be at least two more gridpoints greater than start range = "+zs) |
---|
877 | print(" ") |
---|
878 | exit |
---|
879 | end if |
---|
880 | else |
---|
881 | if (ze .LT. 0) |
---|
882 | print(" ") |
---|
883 | print("range end for z-coordinate = "+ze+" is lower than first gridpoint = 0") |
---|
884 | print(" ") |
---|
885 | exit |
---|
886 | end if |
---|
887 | if (ze .LT. zs) then |
---|
888 | print(" ") |
---|
889 | print("range end for z-coordinate = "+ze+" is lower than start range = "+zs) |
---|
890 | print(" ") |
---|
891 | exit |
---|
892 | end if |
---|
893 | end if |
---|
894 | end if |
---|
895 | end if |
---|
896 | |
---|
897 | if (delta_x .NE. -1) then |
---|
898 | do i=0,xdim |
---|
899 | if (xs .GT. x_d(i)-delta_x/2 .AND. xs .LE. x_d(i)+delta_x/2)then |
---|
900 | xstart=i |
---|
901 | break |
---|
902 | end if |
---|
903 | end do |
---|
904 | do i=0,xdim |
---|
905 | if (xe .GT. x_d(i)-delta_x/2 .AND. xe .LE. x_d(i)+delta_x/2)then |
---|
906 | xend=i |
---|
907 | break |
---|
908 | end if |
---|
909 | end do |
---|
910 | end if |
---|
911 | if (delta_y .NE. -1) then |
---|
912 | do i=0,ydim |
---|
913 | if (ys .GT. y_d(i)-delta_y/2 .AND. ys .LE. y_d(i)+delta_y/2)then |
---|
914 | ystart=i |
---|
915 | break |
---|
916 | end if |
---|
917 | end do |
---|
918 | do i=0,ydim |
---|
919 | if (ye .GT. y_d(i)-delta_y/2 .AND. ye .LE. y_d(i)+delta_y/2)then |
---|
920 | yend=i |
---|
921 | break |
---|
922 | end if |
---|
923 | end do |
---|
924 | end if |
---|
925 | |
---|
926 | if( shape .eq. 1 ) then |
---|
927 | if (xyc .EQ. 1)then |
---|
928 | cs_res@vpWidthF = (xe-xs)/(ye-ys) |
---|
929 | cs_res@vpHeightF = 1 |
---|
930 | vecres@vpWidthF = (xe-xs)/(ye-ys) |
---|
931 | vecres@vpHeightF = 1 |
---|
932 | if (xe-xs .GT. ye-ys)then |
---|
933 | if (no_rows .LE. no_columns)then |
---|
934 | cs_res@gsnPaperOrientation = "landscape" |
---|
935 | vecres@gsnPaperOrientation = "landscape" |
---|
936 | cs_res@lbOrientation = "Horizontal" |
---|
937 | else |
---|
938 | cs_res@gsnPaperOrientation = "portrait" |
---|
939 | vecres@gsnPaperOrientation = "portrait" |
---|
940 | cs_res@lbOrientation = "Vertical" |
---|
941 | end if |
---|
942 | else |
---|
943 | if (no_rows .GE. no_columns)then |
---|
944 | cs_res@gsnPaperOrientation = "portrait" |
---|
945 | vecres@gsnPaperOrientation = "portrait" |
---|
946 | cs_res@lbOrientation = "Vertical" |
---|
947 | else |
---|
948 | cs_res@gsnPaperOrientation = "landscape" |
---|
949 | vecres@gsnPaperOrientation = "landscape" |
---|
950 | cs_res@lbOrientation = "Horizontal" |
---|
951 | end if |
---|
952 | end if |
---|
953 | end if |
---|
954 | if (xzc .EQ. 1)then |
---|
955 | cs_res@vpWidthF = (xe-xs)/(delta_x*(ze-zs)) |
---|
956 | cs_res@vpHeightF = 1 |
---|
957 | vecres@vpWidthF = (xe-xs)/(delta_x*(ze-zs)) |
---|
958 | vecres@vpHeightF = 1 |
---|
959 | if (xe-xs .GT. (delta_x*(ze-zs)))then |
---|
960 | if (no_rows .LE. no_columns)then |
---|
961 | cs_res@gsnPaperOrientation = "landscape" |
---|
962 | vecres@gsnPaperOrientation = "landscape" |
---|
963 | cs_res@lbOrientation = "Horizontal" |
---|
964 | else |
---|
965 | cs_res@gsnPaperOrientation = "portrait" |
---|
966 | vecres@gsnPaperOrientation = "portrait" |
---|
967 | cs_res@lbOrientation = "Vertical" |
---|
968 | end if |
---|
969 | else |
---|
970 | if (no_rows .GE. no_columns)then |
---|
971 | cs_res@gsnPaperOrientation = "portrait" |
---|
972 | vecres@gsnPaperOrientation = "portrait" |
---|
973 | cs_res@lbOrientation = "Vertical" |
---|
974 | else |
---|
975 | cs_res@gsnPaperOrientation = "landscape" |
---|
976 | vecres@gsnPaperOrientation = "landscape" |
---|
977 | cs_res@lbOrientation = "Horizontal" |
---|
978 | end if |
---|
979 | end if |
---|
980 | end if |
---|
981 | if (yzc .EQ. 1)then |
---|
982 | cs_res@vpWidthF = (ye-ys)/(delta_y*(ze-zs)) |
---|
983 | cs_res@vpHeightF = 1 |
---|
984 | vecres@vpWidthF = (ye-ys)/(delta_y*(ze-zs)) |
---|
985 | vecres@vpHeightF = 1 |
---|
986 | if (ye-ys .GT. (delta_y*(ze-zs)))then |
---|
987 | if (no_rows .LE. no_columns)then |
---|
988 | cs_res@gsnPaperOrientation = "landscape" |
---|
989 | vecres@gsnPaperOrientation = "landscape" |
---|
990 | cs_res@lbOrientation = "Horizontal" |
---|
991 | else |
---|
992 | cs_res@gsnPaperOrientation = "portrait" |
---|
993 | vecres@gsnPaperOrientation = "portrait" |
---|
994 | cs_res@lbOrientation = "Vertical" |
---|
995 | end if |
---|
996 | else |
---|
997 | if (no_rows .GE. no_columns)then |
---|
998 | cs_res@gsnPaperOrientation = "portrait" |
---|
999 | vecres@gsnPaperOrientation = "portrait" |
---|
1000 | cs_res@lbOrientation = "Vertical" |
---|
1001 | else |
---|
1002 | cs_res@gsnPaperOrientation = "landscape" |
---|
1003 | vecres@gsnPaperOrientation = "landscape" |
---|
1004 | cs_res@lbOrientation = "Horizontal" |
---|
1005 | end if |
---|
1006 | end if |
---|
1007 | end if |
---|
1008 | end if |
---|
1009 | |
---|
1010 | delete(xs) |
---|
1011 | delete(xe) |
---|
1012 | delete(ys) |
---|
1013 | delete(ye) |
---|
1014 | |
---|
1015 | xs=xstart |
---|
1016 | xe=xend |
---|
1017 | ys=ystart |
---|
1018 | ye=yend |
---|
1019 | |
---|
1020 | if (xyc .EQ. 1) then |
---|
1021 | d=(ye-ys+1)/(major_ticks_y-1) |
---|
1022 | e=(xe-xs+1)/(major_ticks_x-1) |
---|
1023 | array_yl =new(major_ticks_y,integer) |
---|
1024 | array_yl_labels=new(major_ticks_y,double) |
---|
1025 | array_minor_yl =new((major_ticks_y-1)*4,double) |
---|
1026 | array_xb =new(major_ticks_x,integer) |
---|
1027 | array_xb_labels=new(major_ticks_x,double) |
---|
1028 | array_minor_xb =new((major_ticks_x-1)*4,double) |
---|
1029 | array_yl(0)=ys |
---|
1030 | array_xb(0)=xs |
---|
1031 | array_yl_labels(0)=y_d(ys) |
---|
1032 | array_xb_labels(0)=x_d(xs) |
---|
1033 | do ar=1,major_ticks_y-1 |
---|
1034 | array_yl(ar)=d*(ar-1)+d-1 |
---|
1035 | array_yl_labels(ar) = y_d(array_yl(ar)) |
---|
1036 | if (ar .GT. 0) |
---|
1037 | do min_ar=0,3 |
---|
1038 | array_minor_yl(4*(ar-1)+min_ar)= int2dble(array_yl(ar-1))+int2dble(array_yl(ar)-array_yl(ar-1))/5*(min_ar+1) |
---|
1039 | end do |
---|
1040 | end if |
---|
1041 | end do |
---|
1042 | do br=1,major_ticks_x-1 |
---|
1043 | array_xb(br)=e*(br-1)+e-1 |
---|
1044 | array_xb_labels(br) = x_d(array_xb(br)) |
---|
1045 | if (br .GT. 0) |
---|
1046 | do min_br=0,3 |
---|
1047 | array_minor_xb(4*(br-1)+min_br)= int2dble(array_xb(br-1))+int2dble(array_xb(br)-array_xb(br-1))/5*(min_br+1) |
---|
1048 | end do |
---|
1049 | end if |
---|
1050 | end do |
---|
1051 | end if |
---|
1052 | |
---|
1053 | if (xzc .EQ. 1) then |
---|
1054 | d=(ze-zs+1)/(major_ticks_y-1) |
---|
1055 | e=(xe-xs+1)/(major_ticks_x-1) |
---|
1056 | array_yl =new(major_ticks_y,integer) |
---|
1057 | array_yl_labels=new(major_ticks_y,double) |
---|
1058 | array_minor_yl =new((major_ticks_y-1)*4,double) |
---|
1059 | array_xb =new(major_ticks_x,integer) |
---|
1060 | array_xb_labels=new(major_ticks_x,double) |
---|
1061 | array_minor_xb =new((major_ticks_x-1)*4,double) |
---|
1062 | array_yl(0)=zs |
---|
1063 | array_xb(0)=xs |
---|
1064 | array_yl_labels(0)=z_d(zs) |
---|
1065 | array_xb_labels(0)=x_d(xs) |
---|
1066 | do ar=1,major_ticks_y-1 |
---|
1067 | array_yl(ar)=d*(ar-1)+d-1 |
---|
1068 | array_yl_labels(ar) = z_d(array_yl(ar)) |
---|
1069 | if (ar .GT. 0) |
---|
1070 | do min_ar=0,3 |
---|
1071 | array_minor_yl(4*(ar-1)+min_ar)= int2dble(array_yl(ar-1))+int2dble(array_yl(ar)-array_yl(ar-1))/5*(min_ar+1) |
---|
1072 | end do |
---|
1073 | end if |
---|
1074 | end do |
---|
1075 | do br=1,major_ticks_x-1 |
---|
1076 | array_xb(br)=e*(br-1)+e-1 |
---|
1077 | array_xb_labels(br) = x_d(array_xb(br)) |
---|
1078 | if (br .GT. 0) |
---|
1079 | do min_br=0,3 |
---|
1080 | array_minor_xb(4*(br-1)+min_br)= int2dble(array_xb(br-1))+int2dble(array_xb(br)-array_xb(br-1))/5*(min_br+1) |
---|
1081 | end do |
---|
1082 | end if |
---|
1083 | end do |
---|
1084 | end if |
---|
1085 | |
---|
1086 | if (yzc .EQ. 1) then |
---|
1087 | d=(ze-zs+1)/(major_ticks_y-1) |
---|
1088 | e=(ye-ys+1)/(major_ticks_x-1) |
---|
1089 | array_yl =new(major_ticks_y,integer) |
---|
1090 | array_yl_labels=new(major_ticks_y,double) |
---|
1091 | array_minor_yl =new((major_ticks_y-1)*4,double) |
---|
1092 | array_xb =new(major_ticks_x,integer) |
---|
1093 | array_xb_labels=new(major_ticks_x,double) |
---|
1094 | array_minor_xb =new((major_ticks_x-1)*4,double) |
---|
1095 | array_yl(0)=zs |
---|
1096 | array_xb(0)=ys |
---|
1097 | array_yl_labels(0)=z_d(zs) |
---|
1098 | array_xb_labels(0)=y_d(ys) |
---|
1099 | do ar=1,major_ticks_y-1 |
---|
1100 | array_yl(ar)=d*(ar-1)+d-1 |
---|
1101 | array_yl_labels(ar) = y_d(array_yl(ar)) |
---|
1102 | if (ar .GT. 0) |
---|
1103 | do min_ar=0,3 |
---|
1104 | array_minor_yl(4*(ar-1)+min_ar)= int2dble(array_yl(ar-1))+int2dble(array_yl(ar)-array_yl(ar-1))/5*(min_ar+1) |
---|
1105 | end do |
---|
1106 | end if |
---|
1107 | end do |
---|
1108 | do br=1,major_ticks_x-1 |
---|
1109 | array_xb(br)=e*(br-1)+e-1 |
---|
1110 | array_xb_labels(br) = x_d(array_xb(br)) |
---|
1111 | if (br .GT. 0) |
---|
1112 | do min_br=0,3 |
---|
1113 | array_minor_xb(4*(br-1)+min_br)= int2dble(array_xb(br-1))+int2dble(array_xb(br)-array_xb(br-1))/5*(min_br+1) |
---|
1114 | end do |
---|
1115 | end if |
---|
1116 | end do |
---|
1117 | end if |
---|
1118 | |
---|
1119 | if (axes_explicit .EQ. 1)then |
---|
1120 | cs_res@tmYLMode = "Explicit" |
---|
1121 | cs_res@tmXBMode = "Explicit" |
---|
1122 | cs_res@tmYLValues = array_yl |
---|
1123 | if (xyc .EQ. 1)then |
---|
1124 | cs_res@tmYLLabels = array_yl_labels/norm_y |
---|
1125 | cs_res@tmXBLabels = array_xb_labels/norm_x |
---|
1126 | if (norm_x .NE. 1.)then |
---|
1127 | cs_res@tiXAxisString = "x ["+unit_x+"]" |
---|
1128 | else |
---|
1129 | cs_res@tiXAxisString = "x [m]" |
---|
1130 | end if |
---|
1131 | if (norm_y .NE. 1.)then |
---|
1132 | cs_res@tiYAxisString = "y ["+unit_y+"]" |
---|
1133 | else |
---|
1134 | cs_res@tiYAxisString = "y [m]" |
---|
1135 | end if |
---|
1136 | end if |
---|
1137 | if (xzc .EQ. 1)then |
---|
1138 | cs_res@tmYLLabels = array_yl_labels/norm_z |
---|
1139 | cs_res@tmXBLabels = array_xb_labels/norm_x |
---|
1140 | if (norm_x .NE. 1.)then |
---|
1141 | cs_res@tiXAxisString = "x ["+unit_x+"]" |
---|
1142 | else |
---|
1143 | cs_res@tiXAxisString = "x [m]" |
---|
1144 | end if |
---|
1145 | if (norm_z .NE. 1.)then |
---|
1146 | cs_res@tiYAxisString = "z ["+unit_z+"]" |
---|
1147 | else |
---|
1148 | cs_res@tiYAxisString = "z [m]" |
---|
1149 | end if |
---|
1150 | end if |
---|
1151 | if (yzc .EQ. 1)then |
---|
1152 | cs_res@tmYLLabels = array_yl_labels/norm_z |
---|
1153 | cs_res@tmXBLabels = array_xb_labels/norm_y |
---|
1154 | if (norm_y .NE. 1.)then |
---|
1155 | cs_res@tiXAxisString = "y ["+unit_y+"]" |
---|
1156 | else |
---|
1157 | cs_res@tiXAxisString = "y [m]" |
---|
1158 | end if |
---|
1159 | if (norm_z .NE. 1.)then |
---|
1160 | cs_res@tiYAxisString = "z ["+unit_z+"]" |
---|
1161 | else |
---|
1162 | cs_res@tiYAxisString = "z [m]" |
---|
1163 | end if |
---|
1164 | end if |
---|
1165 | cs_res@tmXBValues = array_xb |
---|
1166 | cs_res@tmYLMinorValues = array_minor_yl |
---|
1167 | cs_res@tmXBMinorValues = array_minor_xb |
---|
1168 | else |
---|
1169 | if (axes_explicit .EQ. 0)then |
---|
1170 | cs_res@tmYLMinorPerMajor = 4 |
---|
1171 | cs_res@tmXBMinorPerMajor = 4 |
---|
1172 | else |
---|
1173 | print(" ") |
---|
1174 | print("'axes_explicit' is invalid and set to 0") |
---|
1175 | print(" ") |
---|
1176 | axes_explicit = 0 |
---|
1177 | cs_res@tmYLMinorPerMajor = 4 |
---|
1178 | cs_res@tmXBMinorPerMajor = 4 |
---|
1179 | end if |
---|
1180 | if (xyc .EQ. 1)then |
---|
1181 | cs_res@tiXAxisString = "x [m]" |
---|
1182 | cs_res@tiYAxisString = "y [m]" |
---|
1183 | end if |
---|
1184 | if (xzc .EQ. 1)then |
---|
1185 | cs_res@tiXAxisString = "x [m]" |
---|
1186 | cs_res@tiYAxisString = "z [m]" |
---|
1187 | end if |
---|
1188 | if (yzc .EQ. 1)then |
---|
1189 | cs_res@tiXAxisString = "y [m]" |
---|
1190 | cs_res@tiYAxisString = "z [m]" |
---|
1191 | end if |
---|
1192 | end if |
---|
1193 | |
---|
1194 | if (xyc .EQ. 1 .OR. xzc .EQ.1)then |
---|
1195 | if (xe .EQ. xs+1) then |
---|
1196 | print(" ") |
---|
1197 | print("range end for x-coordinate="+xe*delta_x+"m) must be at least two") |
---|
1198 | print("more gridpoints(="+2*delta_x+"m) greater than start range="+xs*delta_x+"m)") |
---|
1199 | print(" ") |
---|
1200 | exit |
---|
1201 | end if |
---|
1202 | end if |
---|
1203 | if (xyc .EQ. 1 .OR. yzc .EQ.1)then |
---|
1204 | if (ye .EQ. ys+1) then |
---|
1205 | print(" ") |
---|
1206 | print("range end for y-coordinate="+ye*delta_y+"m) must be at least two") |
---|
1207 | print("more gridpoints(="+2*delta_y+"m) greater than start range="+ys*delta_y+"m)") |
---|
1208 | print(" ") |
---|
1209 | exit |
---|
1210 | end if |
---|
1211 | end if |
---|
1212 | if (xzc .EQ. 1 .OR. yzc .EQ.1)then |
---|
1213 | if (ze .EQ. zs+1) then |
---|
1214 | print(" ") |
---|
1215 | print("range end for x-coordinate="+ze+") must be at least two") |
---|
1216 | print("more gridpoints greater than start range="+zs+")") |
---|
1217 | print(" ") |
---|
1218 | exit |
---|
1219 | end if |
---|
1220 | end if |
---|
1221 | |
---|
1222 | if (xyc .EQ. 1) then |
---|
1223 | no_layer = (ze-zs)+1 |
---|
1224 | data = new((/dim,nt,(ze-zs)+1,(ye-ys)+1,(xe-xs)+1/),float) |
---|
1225 | end if |
---|
1226 | if (xzc .EQ. 1) then |
---|
1227 | no_layer = (ye-ys)+1 |
---|
1228 | data = new((/dim,nt,(ze-zs)+1,(ye-ys)+1,(xe-xs)+1/),float) |
---|
1229 | end if |
---|
1230 | if (yzc .EQ. 1) then |
---|
1231 | no_layer = (xe-xs)+1 |
---|
1232 | data = new((/dim,nt,(ze-zs)+1,(ye-ys)+1,(xe-xs)+1/),float) |
---|
1233 | end if |
---|
1234 | |
---|
1235 | MinVal = new(dim,float) |
---|
1236 | MaxVal = new(dim,float) |
---|
1237 | unit = new(dim,string) |
---|
1238 | |
---|
1239 | ; **************************************************** |
---|
1240 | ; define inner and outer loops depending on "sort" |
---|
1241 | ; **************************************************** |
---|
1242 | |
---|
1243 | if ( xyc .eq. 1 ) then |
---|
1244 | lays = zs |
---|
1245 | laye = ze |
---|
1246 | end if |
---|
1247 | if ( xzc .eq. 1 ) then |
---|
1248 | lays = ys |
---|
1249 | laye = ye |
---|
1250 | end if |
---|
1251 | if ( yzc .eq. 1) then |
---|
1252 | lays = xs |
---|
1253 | laye = xe |
---|
1254 | end if |
---|
1255 | |
---|
1256 | if ( sort .eq. "time" ) then |
---|
1257 | lis = start_time_step |
---|
1258 | lie = end_time_step |
---|
1259 | los = 0 |
---|
1260 | loe = laye-lays |
---|
1261 | else |
---|
1262 | lis = 0 |
---|
1263 | lie = laye-lays |
---|
1264 | los = start_time_step |
---|
1265 | loe = end_time_step |
---|
1266 | end if |
---|
1267 | |
---|
1268 | check = True |
---|
1269 | v1=0 |
---|
1270 | v2=0 |
---|
1271 | no_var=0 |
---|
1272 | n=0 |
---|
1273 | |
---|
1274 | do varn=dim-1,0,1 |
---|
1275 | |
---|
1276 | if (vector .EQ. 1) then |
---|
1277 | check_vec1 = isStrSubset( vec1,","+vNam(varn)+",") |
---|
1278 | check_vec2 = isStrSubset( vec2,","+vNam(varn)+",") |
---|
1279 | end if |
---|
1280 | |
---|
1281 | if ( vNam(varn) .eq. "zw" .or. vNam(varn) .eq. "yv" .or. vNam(varn) .eq. "xu" .or. vNam(varn) .eq. "zu" .or. vNam(varn) .eq. "zwwi" .or. vNam(varn) .eq. "zusi" .or. vNam(varn) .eq. "time" .or. vNam(varn) .eq. "zu_3d" .or. vNam(varn) .eq. "zw_3d" .or. vNam(varn) .eq. "x".or. vNam(varn) .eq. "y" .or. vNam(varn) .eq. "zu_xy" .or. vNam(varn) .eq. "zw_xy" .or. vNam(varn) .eq. "zu1_xy" .or. vNam(varn) .eq. "ind_z_xy" .or. vNam(varn) .eq. "y_xz" .or. vNam(varn) .eq. "yv_xz" .or. vNam(varn) .eq. "ind_y_xz" .or. vNam(varn) .eq. "x_yz" .or. vNam(varn) .eq. "xu_yz" .or. vNam(varn) .eq. "ind_x_yz") then |
---|
1282 | check = False |
---|
1283 | else |
---|
1284 | check = True |
---|
1285 | end if |
---|
1286 | |
---|
1287 | if(check) then |
---|
1288 | |
---|
1289 | no_var=no_var+1 |
---|
1290 | |
---|
1291 | if (vector .EQ. 1) then |
---|
1292 | if (","+vNam(varn)+"," .EQ. vec1) then |
---|
1293 | v1=v1+1 |
---|
1294 | end if |
---|
1295 | if (","+vNam(varn)+"," .EQ. vec2) then |
---|
1296 | v2=v2+1 |
---|
1297 | end if |
---|
1298 | end if |
---|
1299 | |
---|
1300 | if (xyc .EQ. 1) then |
---|
1301 | temp = f[:]->$vNam(varn)$ |
---|
1302 | data_att = f_att->$vNam(varn)$ |
---|
1303 | data(varn,:,:,:,:)=temp(:,zs:ze,ys:ye,xs:xe) |
---|
1304 | end if |
---|
1305 | if ( xzc .eq. 1 ) then |
---|
1306 | temp = f[:]->$vNam(varn)$ |
---|
1307 | data_att = f_att->$vNam(varn)$ |
---|
1308 | data(varn,:,:,:,:)=temp(:,zs:ze,ys:ye,xs:xe) |
---|
1309 | end if |
---|
1310 | if ( yzc .eq. 1) then |
---|
1311 | temp = f[:]->$vNam(varn)$ |
---|
1312 | data_att = f_att->$vNam(varn)$ |
---|
1313 | data(varn,:,:,:,:)=temp(:,zs:ze,ys:ye,xs:xe) |
---|
1314 | end if |
---|
1315 | if (check_vec1) then |
---|
1316 | vect1=data(varn,:,:,:,:) |
---|
1317 | end if |
---|
1318 | if (check_vec2) then |
---|
1319 | vect2=data(varn,:,:,:,:) |
---|
1320 | end if |
---|
1321 | |
---|
1322 | data!0 = "var" |
---|
1323 | data!1 = "t" |
---|
1324 | data!2 = "z" |
---|
1325 | data!3 = "y" |
---|
1326 | data!4 = "x" |
---|
1327 | |
---|
1328 | MinVal(varn) = min(data(varn,:,:,:,:)) |
---|
1329 | MaxVal(varn) = max(data(varn,:,:,:,:)) |
---|
1330 | |
---|
1331 | unit(varn) = data_att@units |
---|
1332 | |
---|
1333 | end if |
---|
1334 | |
---|
1335 | end do |
---|
1336 | |
---|
1337 | var_input=new(no_var,string) |
---|
1338 | numb=0 |
---|
1339 | no_var=0 |
---|
1340 | |
---|
1341 | do varn=dim-1,0,1 |
---|
1342 | |
---|
1343 | if (vector .EQ. 1) then |
---|
1344 | check_vec1 = isStrSubset( vec1,","+vNam(varn)+",") |
---|
1345 | check_vec2 = isStrSubset( vec2,","+vNam(varn)+",") |
---|
1346 | end if |
---|
1347 | |
---|
1348 | if ( vNam(varn) .eq. "zw" .or. vNam(varn) .eq. "yv" .or. vNam(varn) .eq. "xu" .or. vNam(varn) .eq. "zu" .or. vNam(varn) .eq. "zwwi" .or. vNam(varn) .eq. "zusi" .or. vNam(varn) .eq. "time" .or. vNam(varn) .eq. "zu_3d" .or. vNam(varn) .eq. "zw_3d" .or. vNam(varn) .eq. "x".or. vNam(varn) .eq. "y" .or. vNam(varn) .eq. "zu_xy" .or. vNam(varn) .eq. "zw_xy" .or. vNam(varn) .eq. "zu1_xy" .or. vNam(varn) .eq. "ind_z_xy" .or. vNam(varn) .eq. "y_xz" .or. vNam(varn) .eq. "yv_xz" .or. vNam(varn) .eq. "ind_y_xz" .or. vNam(varn) .eq. "x_yz" .or. vNam(varn) .eq. "xu_yz" .or. vNam(varn) .eq. "ind_x_yz") then |
---|
1349 | check = False |
---|
1350 | else |
---|
1351 | check = True |
---|
1352 | end if |
---|
1353 | |
---|
1354 | if(check) then |
---|
1355 | var_input(numb)=vNam(varn) |
---|
1356 | numb=numb+1 |
---|
1357 | end if |
---|
1358 | |
---|
1359 | if (var .NE. "all") then |
---|
1360 | check = isStrSubset( var,","+vNam(varn)+"," ) |
---|
1361 | end if |
---|
1362 | |
---|
1363 | if (check) then |
---|
1364 | no_var = no_var+1 |
---|
1365 | end if |
---|
1366 | |
---|
1367 | end do |
---|
1368 | |
---|
1369 | if (no_var .EQ. 0) then |
---|
1370 | print(" ") |
---|
1371 | print("The variables var='"+var+"' do not exist on your input file;") |
---|
1372 | print("be sure to have one comma berfore and after each variable") |
---|
1373 | print(" ") |
---|
1374 | exit |
---|
1375 | end if |
---|
1376 | |
---|
1377 | if (vector .EQ. 1) then |
---|
1378 | if (v1 .EQ. 0)then |
---|
1379 | print(" ") |
---|
1380 | print("Component 1 for the vector-plot ('vec1') must be one of the varibles on the input file:") |
---|
1381 | print("- "+var_input) |
---|
1382 | print("be sure to have one comma berfore and after the variable") |
---|
1383 | print(" ") |
---|
1384 | exit |
---|
1385 | end if |
---|
1386 | |
---|
1387 | if (v2 .EQ. 0)then |
---|
1388 | print(" ") |
---|
1389 | print("Component 2 for the vector-plot ('vec2') must be one of the varibles on the input file:") |
---|
1390 | print("- "+var_input) |
---|
1391 | print("be sure to have one comma berfore and after the variable") |
---|
1392 | print(" ") |
---|
1393 | exit |
---|
1394 | end if |
---|
1395 | end if |
---|
1396 | |
---|
1397 | ; *************************************************** |
---|
1398 | ; open workstation(s) |
---|
1399 | ; *************************************************** |
---|
1400 | |
---|
1401 | wks_ps = gsn_open_wks(format_out,file_out) |
---|
1402 | gsn_define_colormap(wks_ps,"rainbow+white") |
---|
1403 | |
---|
1404 | if (vector .EQ. 1 .AND. plotvec .EQ. "plotvec") then |
---|
1405 | plot=new((/no_time*no_layer*no_var+no_time*no_layer/),graphic) |
---|
1406 | else |
---|
1407 | plot=new((/no_time*no_layer*no_var/),graphic) |
---|
1408 | end if |
---|
1409 | |
---|
1410 | page = 0 |
---|
1411 | n=0 |
---|
1412 | print(" ") |
---|
1413 | print("Plots sorted by '"+sort+"'") |
---|
1414 | print(" ") |
---|
1415 | |
---|
1416 | ; *************************************************** |
---|
1417 | ; create plots |
---|
1418 | ; *************************************************** |
---|
1419 | |
---|
1420 | if (vector .EQ. 1 .AND. plotvec .EQ. "plotvec") then |
---|
1421 | do lo = los, loe |
---|
1422 | do li = lis, lie |
---|
1423 | if (xyc .EQ. 1)then |
---|
1424 | if (sort .EQ. "time")then |
---|
1425 | level = "z=" + z_d(lo) + "m" |
---|
1426 | else |
---|
1427 | level = "z=" + z_d(li) + "m" |
---|
1428 | end if |
---|
1429 | end if |
---|
1430 | if (xzc .EQ. 1)then |
---|
1431 | if (sort .EQ. "time")then |
---|
1432 | level = "y=" + y_d(lo) + "m" |
---|
1433 | else |
---|
1434 | level = "y=" + y_d(li) + "m" |
---|
1435 | end if |
---|
1436 | end if |
---|
1437 | if (yzc .EQ. 1)then |
---|
1438 | if (sort .EQ. "time")then |
---|
1439 | level = "x=" + x_d(lo) + "m" |
---|
1440 | else |
---|
1441 | level = "x=" + x_d(li) + "m" |
---|
1442 | end if |
---|
1443 | end if |
---|
1444 | vecres = True ; vector only resources |
---|
1445 | vecres@gsnDraw = False ; don't draw |
---|
1446 | vecres@gsnFrame = False ; don't advance frame |
---|
1447 | vecres@vcGlyphStyle = "CurlyVector" ; curly vectors |
---|
1448 | vecres@vcRefMagnitudeF = ref_mag ; define vector ref mag |
---|
1449 | vecres@vcRefLengthF = 0.05 ; define length of vec ref |
---|
1450 | vecres@gsnRightString = "Vector plot of "+vec1+" and "+vec2 |
---|
1451 | vecres@tmXBLabelFontHeightF = font_size |
---|
1452 | vecres@tmYLLabelFontHeightF = font_size |
---|
1453 | vecres@tiXAxisFontHeightF = font_size |
---|
1454 | vecres@tiYAxisFontHeightF = font_size |
---|
1455 | if (sort .EQ. "time")then |
---|
1456 | vecres@gsnLeftString = "t=" + decimalPlaces(t_all(li)/3600,2,True) +"h "+level |
---|
1457 | else |
---|
1458 | vecres@gsnLeftString = "t=" + decimalPlaces(t_all(lo)/3600,2,True) +"h "+level |
---|
1459 | end if |
---|
1460 | vecres@tiXAxisString = " " |
---|
1461 | if (xyc .EQ. 1)then |
---|
1462 | if (sort .EQ. "time")then |
---|
1463 | plot(n) = gsn_csm_vector(wks_ps,vect1(li,lo,:,:),vect2(li,lo,:,:),vecres) |
---|
1464 | else |
---|
1465 | plot(n) = gsn_csm_vector(wks_ps,vect1(lo,li,:,:),vect2(lo,li,:,:),vecres) |
---|
1466 | end if |
---|
1467 | end if |
---|
1468 | if (xzc .EQ. 1) then |
---|
1469 | if (sort .EQ. "time")then |
---|
1470 | plot(n) = gsn_csm_vector(wks_ps,vect1(li,:,lo,:),vect2(li,:,lo,:),vecres) |
---|
1471 | else |
---|
1472 | plot(n) = gsn_csm_vector(wks_ps,vect1(lo,:,li,:),vect2(lo,:,li,:),vecres) |
---|
1473 | end if |
---|
1474 | end if |
---|
1475 | if (yzc .EQ. 1) then |
---|
1476 | if (sort .EQ. "time")then |
---|
1477 | plot(n) = gsn_csm_vector(wks_ps,vect1(li,:,:,lo),vect2(li,:,:,lo),vecres) |
---|
1478 | else |
---|
1479 | plot(n) = gsn_csm_vector(wks_ps,vect1(lo,:,:,li),vect2(lo,:,:,li),vecres) |
---|
1480 | end if |
---|
1481 | end if |
---|
1482 | n=n+1 |
---|
1483 | end do |
---|
1484 | end do |
---|
1485 | end if |
---|
1486 | |
---|
1487 | do varn=dim-1,0,1 |
---|
1488 | |
---|
1489 | if (vector .EQ. 1) then |
---|
1490 | check_vecp = isStrSubset( plotvec,","+vNam(varn)+",") |
---|
1491 | end if |
---|
1492 | |
---|
1493 | if ( vNam(varn) .eq. "zw" .or. vNam(varn) .eq. "yv" .or. vNam(varn) .eq. "xu" .or. vNam(varn) .eq. "zu" .or. vNam(varn) .eq. "zwwi" .or. vNam(varn) .eq. "zusi" .or. vNam(varn) .eq. "time" .or. vNam(varn) .eq. "zu_3d" .or. vNam(varn) .eq. "zw_3d" .or. vNam(varn) .eq. "x".or. vNam(varn) .eq. "y" .or. vNam(varn) .eq. "zu_xy" .or. vNam(varn) .eq. "zw_xy" .or. vNam(varn) .eq. "zu1_xy" .or. vNam(varn) .eq. "ind_z_xy" .or. vNam(varn) .eq. "y_xz" .or. vNam(varn) .eq. "yv_xz" .or. vNam(varn) .eq. "ind_y_xz" .or. vNam(varn) .eq. "x_yz" .or. vNam(varn) .eq. "xu_yz" .or. vNam(varn) .eq. "ind_x_yz") then |
---|
1494 | check = False |
---|
1495 | else |
---|
1496 | check = True |
---|
1497 | end if |
---|
1498 | |
---|
1499 | if (var .NE. "all") then |
---|
1500 | check = isStrSubset( var,","+vNam(varn)+"," ) |
---|
1501 | end if |
---|
1502 | |
---|
1503 | if(check) then |
---|
1504 | |
---|
1505 | space=(MaxVal(varn)-MinVal(varn))/24 |
---|
1506 | |
---|
1507 | cs_res@cnMinLevelValF = MinVal(varn) |
---|
1508 | cs_res@cnMaxLevelValF = MaxVal(varn) |
---|
1509 | |
---|
1510 | cs_res@cnLevelSpacingF = space |
---|
1511 | |
---|
1512 | ; **************************************************** |
---|
1513 | ; loops over time and layer |
---|
1514 | ; **************************************************** |
---|
1515 | |
---|
1516 | |
---|
1517 | do lo = los, loe |
---|
1518 | do li = lis, lie |
---|
1519 | |
---|
1520 | ; **************************************************** |
---|
1521 | ; xy cross section |
---|
1522 | ; **************************************************** |
---|
1523 | |
---|
1524 | if ( xyc .eq. 1 ) then |
---|
1525 | |
---|
1526 | cs_res@gsnLeftString = "Plot of "+vNam(varn) |
---|
1527 | cs_res@gsnRightString = unit(varn) |
---|
1528 | |
---|
1529 | if ( sort .eq. "time" ) then |
---|
1530 | if ( z_d(zs+lo) .eq. -1)then |
---|
1531 | if (delta_z .EQ. -1) then |
---|
1532 | level = "-average" |
---|
1533 | else |
---|
1534 | level = "=" + z_d(zs+lo) + "m" |
---|
1535 | end if |
---|
1536 | else |
---|
1537 | level = "=" + z_d(zs+lo) + "m" |
---|
1538 | end if |
---|
1539 | cs_res@gsnCenterString = "t=" + decimalPlaces(t_all(li)/3600,2,True) +"h z"+level |
---|
1540 | plot(n) = gsn_csm_contour(wks_ps,data(varn,li,lo,:,:),cs_res) |
---|
1541 | if (vector .EQ. 1 .AND. check_vecp) then |
---|
1542 | vecres = True ; vector only resources |
---|
1543 | vecres@gsnDraw = False ; don't draw |
---|
1544 | vecres@gsnFrame = False ; don't advance frame |
---|
1545 | vecres@vcGlyphStyle = "CurlyVector" ; curly vectors |
---|
1546 | vecres@vcRefMagnitudeF = ref_mag ; define vector ref mag |
---|
1547 | vecres@vcRefLengthF = 0.05 ; define length of vec ref |
---|
1548 | vecres@gsnRightString = " " |
---|
1549 | vecres@gsnLeftString = " " |
---|
1550 | vecres@tiXAxisString = " " |
---|
1551 | plot_vec=gsn_csm_vector(wks_ps,vect1(li,lo,:,:),vect2(li,lo,:,:),vecres) |
---|
1552 | overlay(plot(n), plot_vec) |
---|
1553 | end if |
---|
1554 | end if |
---|
1555 | |
---|
1556 | if ( sort .eq. "layer" ) then |
---|
1557 | if (z_d(zs+li) .eq. -1) then |
---|
1558 | if (delta_z .EQ. -1) then |
---|
1559 | level = "-average" |
---|
1560 | else |
---|
1561 | level = "=" + z_d(zs+li) + "m" |
---|
1562 | end if |
---|
1563 | else |
---|
1564 | level = "=" + z_d(zs+li) + "m" |
---|
1565 | end if |
---|
1566 | cs_res@gsnCenterString = "t=" + decimalPlaces(t_all(lo)/3600,2,True) + "h z"+ level |
---|
1567 | plot(n) = gsn_csm_contour(wks_ps,data(varn,lo,li,:,:),cs_res) |
---|
1568 | if (vector .EQ. 1 .AND. check_vecp) then |
---|
1569 | vecres = True ; vector only resources |
---|
1570 | vecres@gsnDraw = False ; don't draw |
---|
1571 | vecres@gsnFrame = False ; don't advance frame |
---|
1572 | vecres@vcGlyphStyle = "CurlyVector" ; curly vectors |
---|
1573 | vecres@vcRefMagnitudeF = ref_mag ; define vector ref mag |
---|
1574 | vecres@vcRefLengthF = 0.05 ; define length of vec ref |
---|
1575 | vecres@gsnRightString = " " ; turn off right string |
---|
1576 | vecres@gsnLeftString = " " ; turn off left string |
---|
1577 | vecres@tiXAxisString = " " |
---|
1578 | plot_vec=gsn_csm_vector(wks_ps,vect1(lo,li,:,:),vect2(lo,li,:,:),vecres) |
---|
1579 | overlay(plot(n), plot_vec) |
---|
1580 | end if |
---|
1581 | end if |
---|
1582 | end if |
---|
1583 | |
---|
1584 | ; **************************************************** |
---|
1585 | ; xz cross section |
---|
1586 | ; **************************************************** |
---|
1587 | |
---|
1588 | if ( xzc .eq. 1 ) then |
---|
1589 | |
---|
1590 | cs_res@gsnLeftString = "Plot of "+vNam(varn) |
---|
1591 | |
---|
1592 | if ( sort .eq. "time" ) then |
---|
1593 | if ( y_d(ys+lo) .eq. -1 ) then |
---|
1594 | level = "-average" |
---|
1595 | else |
---|
1596 | level = "=" + y_d(ys+lo) + "m" |
---|
1597 | end if |
---|
1598 | cs_res@gsnCenterString = "t=" + decimalPlaces(t_all(li)/3600,2,True) + "h y"+ level |
---|
1599 | plot(n) = gsn_csm_contour(wks_ps,data(varn,li,:,lo,:),cs_res) |
---|
1600 | if (vector .EQ. 1 .AND. check_vecp) then |
---|
1601 | vecres = True ; vector only resources |
---|
1602 | vecres@gsnDraw = False ; don't draw |
---|
1603 | vecres@gsnFrame = False ; don't advance frame |
---|
1604 | vecres@vcGlyphStyle = "CurlyVector" ; curly vectors |
---|
1605 | vecres@vcRefMagnitudeF = ref_mag ; define vector ref mag |
---|
1606 | vecres@vcRefLengthF = 0.05 ; define length of vec ref |
---|
1607 | vecres@gsnRightString = " " ; turn off right string |
---|
1608 | vecres@gsnLeftString = " " ; turn off left string |
---|
1609 | vecres@tiXAxisString = " " |
---|
1610 | plot_vec=gsn_csm_vector(wks_ps,vect1(li,:,lo,:),vect2(li,:,lo,:),vecres) |
---|
1611 | overlay(plot(n), plot_vec) |
---|
1612 | end if |
---|
1613 | end if |
---|
1614 | |
---|
1615 | if ( sort .eq. "layer" ) then |
---|
1616 | if ( y_d(ys+li) .eq. -1 ) then |
---|
1617 | level = "-average" |
---|
1618 | else |
---|
1619 | level = "=" + y_d(ys+li) + "m" |
---|
1620 | end if |
---|
1621 | cs_res@gsnCenterString = "t=" + decimalPlaces(t_all(lo)/3600,2,True) + "h y"+ level |
---|
1622 | plot(n) = gsn_csm_contour(wks_ps,data(varn,lo,:,li,:),cs_res) |
---|
1623 | if (vector .EQ. 1 .AND. check_vecp) then |
---|
1624 | vecres = True ; vector only resources |
---|
1625 | vecres@gsnDraw = False ; don't draw |
---|
1626 | vecres@gsnFrame = False ; don't advance frame |
---|
1627 | vecres@vcGlyphStyle = "CurlyVector" ; curly vectors |
---|
1628 | vecres@vcRefMagnitudeF = ref_mag ; define vector ref mag |
---|
1629 | vecres@vcRefLengthF = 0.05 ; define length of vec ref |
---|
1630 | vecres@gsnRightString = " " ; turn off right string |
---|
1631 | vecres@gsnLeftString = " " ; turn off left string |
---|
1632 | vecres@tiXAxisString = " " |
---|
1633 | plot_vec=gsn_csm_vector(wks_ps,vect1(lo,:,li,:),vect2(lo,:,li,:),vecres) |
---|
1634 | overlay(plot(n), plot_vec) |
---|
1635 | end if |
---|
1636 | end if |
---|
1637 | end if |
---|
1638 | |
---|
1639 | ; **************************************************** |
---|
1640 | ; yz cross section |
---|
1641 | ; **************************************************** |
---|
1642 | |
---|
1643 | if ( yzc .eq. 1 ) then |
---|
1644 | |
---|
1645 | cs_res@gsnLeftString = "Plot of "+vNam(varn) |
---|
1646 | |
---|
1647 | if ( sort .eq. "time" ) then |
---|
1648 | if ( x_d(xs+lo) .eq. -1 ) then |
---|
1649 | level = "-average" |
---|
1650 | else |
---|
1651 | level = "=" + x_d(xs+lo) + "m" |
---|
1652 | end if |
---|
1653 | cs_res@gsnCenterString = "t=" + decimalPlaces(t_all(li)/3600,2,True) + "h x"+ level |
---|
1654 | plot(n) = gsn_csm_contour(wks_ps,data(varn,li,:,:,lo),cs_res) |
---|
1655 | if (vector .EQ. 1 .AND. check_vecp) then |
---|
1656 | vecres = True ; vector only resources |
---|
1657 | vecres@gsnDraw = False ; don't draw |
---|
1658 | vecres@gsnFrame = False ; don't advance frame |
---|
1659 | vecres@vcGlyphStyle = "CurlyVector" ; curly vectors |
---|
1660 | vecres@vcRefMagnitudeF = ref_mag ; define vector ref mag |
---|
1661 | vecres@vcRefLengthF = 0.05 ; define length of vec ref |
---|
1662 | vecres@gsnRightString = " " ; turn off right string |
---|
1663 | vecres@gsnLeftString = " " ; turn off left string |
---|
1664 | vecres@tiXAxisString = " " |
---|
1665 | plot_vec=gsn_csm_vector(wks_ps,vect1(li,:,:,lo),vect2(li,:,:,lo),vecres) |
---|
1666 | overlay(plot(n), plot_vec) |
---|
1667 | end if |
---|
1668 | end if |
---|
1669 | |
---|
1670 | if ( sort .eq. "layer" ) then |
---|
1671 | if ( x_d(xs+li) .eq. -1 ) then |
---|
1672 | level = "-average" |
---|
1673 | else |
---|
1674 | level = "=" + x_d(xs+li) + "m" |
---|
1675 | end if |
---|
1676 | cs_res@gsnCenterString = "t=" + decimalPlaces(t_all(lo)/3600,2,True) + "h x"+ level |
---|
1677 | plot(n) = gsn_csm_contour(wks_ps,data(varn,lo,:,:,li),cs_res) |
---|
1678 | if (vector .EQ. 1 .AND. check_vecp)then |
---|
1679 | vecres = True ; vector only resources |
---|
1680 | vecres@gsnDraw = False ; don't draw |
---|
1681 | vecres@gsnFrame = False ; don't advance frame |
---|
1682 | vecres@vcGlyphStyle = "CurlyVector" ; curly vectors |
---|
1683 | vecres@vcRefMagnitudeF = ref_mag ; define vector ref mag |
---|
1684 | vecres@vcRefLengthF = 0.05 ; define length of vec ref |
---|
1685 | vecres@gsnRightString = " " ; turn off right string |
---|
1686 | vecres@gsnLeftString = " " ; turn off left string |
---|
1687 | vecres@tiXAxisString = " " |
---|
1688 | plot_vec=gsn_csm_vector(wks_ps,vect1(lo,:,:,li),vect2(lo,:,:,li),vecres) |
---|
1689 | overlay(plot(n), plot_vec) |
---|
1690 | end if |
---|
1691 | end if |
---|
1692 | end if |
---|
1693 | n=n+1 |
---|
1694 | end do |
---|
1695 | end do |
---|
1696 | end if |
---|
1697 | end do |
---|
1698 | |
---|
1699 | ; *************************************************** |
---|
1700 | ; merge plots onto one page |
---|
1701 | ; *************************************************** |
---|
1702 | |
---|
1703 | if (vector .EQ. 1 .AND. plotvec .EQ. "plotvec")then |
---|
1704 | if (format_out .EQ. "eps" .OR. format_out .EQ. "epsi") then |
---|
1705 | gsn_panel(wks_ps,plot(0:(no_time*no_layer*(no_var+1))-1),(/no_var+1,no_layer*no_time/),cs_resP) |
---|
1706 | print(" ") |
---|
1707 | print("Outputs to .eps or .epsi have only one frame") |
---|
1708 | print(" ") |
---|
1709 | else |
---|
1710 | do np = 0,no_layer*no_time*(no_var+1)-1,no_rows*no_columns |
---|
1711 | if ( np + no_rows*no_columns .gt. (no_layer*no_time*(no_var+1))-1) then |
---|
1712 | gsn_panel(wks_ps, plot(np:(no_layer*no_time*(no_var+1))-1),(/no_rows,no_columns/),cs_resP) |
---|
1713 | else |
---|
1714 | gsn_panel(wks_ps, plot(np:np+no_rows*no_columns-1),(/no_rows,no_columns/),cs_resP) |
---|
1715 | end if |
---|
1716 | end do |
---|
1717 | end if |
---|
1718 | else |
---|
1719 | if (format_out .EQ. "eps" .OR. format_out .EQ. "epsi") then |
---|
1720 | gsn_panel(wks_ps,plot(0:(no_time*no_layer*no_var)-1),(/no_var,no_layer*no_time/),cs_resP) |
---|
1721 | print(" ") |
---|
1722 | print("Outputs to .eps or .epsi have only one frame") |
---|
1723 | print(" ") |
---|
1724 | else |
---|
1725 | do np = 0,no_layer*no_time*no_var-1,no_rows*no_columns |
---|
1726 | if ( np + no_rows*no_columns .gt. (no_layer*no_time*no_var)-1) then |
---|
1727 | gsn_panel(wks_ps, plot(np:(no_layer*no_time*no_var)-1),(/no_rows,no_columns/),cs_resP) |
---|
1728 | else |
---|
1729 | gsn_panel(wks_ps, plot(np:np+no_rows*no_columns-1),(/no_rows,no_columns/),cs_resP) |
---|
1730 | end if |
---|
1731 | end do |
---|
1732 | end if |
---|
1733 | end if |
---|
1734 | |
---|
1735 | print(" ") |
---|
1736 | print("Output to: " + file_out +"."+ format_out) |
---|
1737 | print(" ") |
---|
1738 | |
---|
1739 | end |
---|