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)= y_d(array_yl(ar-1))+(y_d(array_yl(ar))-y_d(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)= x_d(array_xb(br-1))+(x_d(array_xb(br))-x_d(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)= z_d(array_yl(ar-1))+(z_d(array_yl(ar))-z_d(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)= x_d(array_xb(br-1))+(x_d(array_xb(br))-x_d(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) = z_d(array_yl(ar)) |
---|
1102 | if (ar .GT. 0) |
---|
1103 | do min_ar=0,3 |
---|
1104 | array_minor_yl(4*(ar-1)+min_ar)= z_d(array_yl(ar-1))+(z_d(array_yl(ar))-z_d(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) = y_d(array_xb(br)) |
---|
1111 | if (br .GT. 0) |
---|
1112 | do min_br=0,3 |
---|
1113 | array_minor_xb(4*(br-1)+min_br)= y_d(array_xb(br-1))+(y_d(array_xb(br))-y_d(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 | if (xyc .EQ. 1)then |
---|
1123 | cs_res@tmYLValues = y_d(array_yl) |
---|
1124 | cs_res@tmXBValues = x_d(array_xb) |
---|
1125 | cs_res@tmYLLabels = array_yl_labels/norm_y |
---|
1126 | cs_res@tmXBLabels = array_xb_labels/norm_x |
---|
1127 | if (norm_x .NE. 1.)then |
---|
1128 | cs_res@tiXAxisString = "x ["+unit_x+"]" |
---|
1129 | else |
---|
1130 | cs_res@tiXAxisString = "x [m]" |
---|
1131 | end if |
---|
1132 | if (norm_y .NE. 1.)then |
---|
1133 | cs_res@tiYAxisString = "y ["+unit_y+"]" |
---|
1134 | else |
---|
1135 | cs_res@tiYAxisString = "y [m]" |
---|
1136 | end if |
---|
1137 | end if |
---|
1138 | if (xzc .EQ. 1)then |
---|
1139 | cs_res@tmYLValues = z_d(array_yl) |
---|
1140 | cs_res@tmXBValues = x_d(array_xb) |
---|
1141 | cs_res@tmYLLabels = array_yl_labels/norm_z |
---|
1142 | cs_res@tmXBLabels = array_xb_labels/norm_x |
---|
1143 | if (norm_x .NE. 1.)then |
---|
1144 | cs_res@tiXAxisString = "x ["+unit_x+"]" |
---|
1145 | else |
---|
1146 | cs_res@tiXAxisString = "x [m]" |
---|
1147 | end if |
---|
1148 | if (norm_z .NE. 1.)then |
---|
1149 | cs_res@tiYAxisString = "z ["+unit_z+"]" |
---|
1150 | else |
---|
1151 | cs_res@tiYAxisString = "z [m]" |
---|
1152 | end if |
---|
1153 | end if |
---|
1154 | if (yzc .EQ. 1)then |
---|
1155 | cs_res@tmYLValues = z_d(array_yl) |
---|
1156 | cs_res@tmXBValues = y_d(array_xb) |
---|
1157 | cs_res@tmYLLabels = array_yl_labels/norm_z |
---|
1158 | cs_res@tmXBLabels = array_xb_labels/norm_y |
---|
1159 | if (norm_y .NE. 1.)then |
---|
1160 | cs_res@tiXAxisString = "y ["+unit_y+"]" |
---|
1161 | else |
---|
1162 | cs_res@tiXAxisString = "y [m]" |
---|
1163 | end if |
---|
1164 | if (norm_z .NE. 1.)then |
---|
1165 | cs_res@tiYAxisString = "z ["+unit_z+"]" |
---|
1166 | else |
---|
1167 | cs_res@tiYAxisString = "z [m]" |
---|
1168 | end if |
---|
1169 | end if |
---|
1170 | cs_res@tmYLMinorValues = array_minor_yl |
---|
1171 | cs_res@tmXBMinorValues = array_minor_xb |
---|
1172 | else |
---|
1173 | if (axes_explicit .EQ. 0)then |
---|
1174 | cs_res@tmYLMinorPerMajor = 4 |
---|
1175 | cs_res@tmXBMinorPerMajor = 4 |
---|
1176 | else |
---|
1177 | print(" ") |
---|
1178 | print("'axes_explicit' is invalid and set to 0") |
---|
1179 | print(" ") |
---|
1180 | axes_explicit = 0 |
---|
1181 | cs_res@tmYLMinorPerMajor = 4 |
---|
1182 | cs_res@tmXBMinorPerMajor = 4 |
---|
1183 | end if |
---|
1184 | if (xyc .EQ. 1)then |
---|
1185 | cs_res@tiXAxisString = "x [m]" |
---|
1186 | cs_res@tiYAxisString = "y [m]" |
---|
1187 | end if |
---|
1188 | if (xzc .EQ. 1)then |
---|
1189 | cs_res@tiXAxisString = "x [m]" |
---|
1190 | cs_res@tiYAxisString = "z [m]" |
---|
1191 | end if |
---|
1192 | if (yzc .EQ. 1)then |
---|
1193 | cs_res@tiXAxisString = "y [m]" |
---|
1194 | cs_res@tiYAxisString = "z [m]" |
---|
1195 | end if |
---|
1196 | end if |
---|
1197 | |
---|
1198 | if (xyc .EQ. 1 .OR. xzc .EQ.1)then |
---|
1199 | if (xe .EQ. xs+1) then |
---|
1200 | print(" ") |
---|
1201 | print("range end for x-coordinate="+xe*delta_x+"m) must be at least two") |
---|
1202 | print("more gridpoints(="+2*delta_x+"m) greater than start range="+xs*delta_x+"m)") |
---|
1203 | print(" ") |
---|
1204 | exit |
---|
1205 | end if |
---|
1206 | end if |
---|
1207 | if (xyc .EQ. 1 .OR. yzc .EQ.1)then |
---|
1208 | if (ye .EQ. ys+1) then |
---|
1209 | print(" ") |
---|
1210 | print("range end for y-coordinate="+ye*delta_y+"m) must be at least two") |
---|
1211 | print("more gridpoints(="+2*delta_y+"m) greater than start range="+ys*delta_y+"m)") |
---|
1212 | print(" ") |
---|
1213 | exit |
---|
1214 | end if |
---|
1215 | end if |
---|
1216 | if (xzc .EQ. 1 .OR. yzc .EQ.1)then |
---|
1217 | if (ze .EQ. zs+1) then |
---|
1218 | print(" ") |
---|
1219 | print("range end for x-coordinate="+ze+") must be at least two") |
---|
1220 | print("more gridpoints greater than start range="+zs+")") |
---|
1221 | print(" ") |
---|
1222 | exit |
---|
1223 | end if |
---|
1224 | end if |
---|
1225 | |
---|
1226 | if (xyc .EQ. 1) then |
---|
1227 | no_layer = (ze-zs)+1 |
---|
1228 | data = new((/dim,nt,(ze-zs)+1,(ye-ys)+1,(xe-xs)+1/),float) |
---|
1229 | end if |
---|
1230 | if (xzc .EQ. 1) then |
---|
1231 | no_layer = (ye-ys)+1 |
---|
1232 | data = new((/dim,nt,(ze-zs)+1,(ye-ys)+1,(xe-xs)+1/),float) |
---|
1233 | end if |
---|
1234 | if (yzc .EQ. 1) then |
---|
1235 | no_layer = (xe-xs)+1 |
---|
1236 | data = new((/dim,nt,(ze-zs)+1,(ye-ys)+1,(xe-xs)+1/),float) |
---|
1237 | end if |
---|
1238 | |
---|
1239 | MinVal = new(dim,float) |
---|
1240 | MaxVal = new(dim,float) |
---|
1241 | unit = new(dim,string) |
---|
1242 | |
---|
1243 | ; **************************************************** |
---|
1244 | ; define inner and outer loops depending on "sort" |
---|
1245 | ; **************************************************** |
---|
1246 | |
---|
1247 | if ( xyc .eq. 1 ) then |
---|
1248 | lays = zs |
---|
1249 | laye = ze |
---|
1250 | end if |
---|
1251 | if ( xzc .eq. 1 ) then |
---|
1252 | lays = ys |
---|
1253 | laye = ye |
---|
1254 | end if |
---|
1255 | if ( yzc .eq. 1) then |
---|
1256 | lays = xs |
---|
1257 | laye = xe |
---|
1258 | end if |
---|
1259 | |
---|
1260 | if ( sort .eq. "time" ) then |
---|
1261 | lis = start_time_step |
---|
1262 | lie = end_time_step |
---|
1263 | los = 0 |
---|
1264 | loe = laye-lays |
---|
1265 | else |
---|
1266 | lis = 0 |
---|
1267 | lie = laye-lays |
---|
1268 | los = start_time_step |
---|
1269 | loe = end_time_step |
---|
1270 | end if |
---|
1271 | |
---|
1272 | check = True |
---|
1273 | v1=0 |
---|
1274 | v2=0 |
---|
1275 | no_var=0 |
---|
1276 | n=0 |
---|
1277 | |
---|
1278 | do varn=dim-1,0,1 |
---|
1279 | |
---|
1280 | if (vector .EQ. 1) then |
---|
1281 | check_vec1 = isStrSubset( vec1,","+vNam(varn)+",") |
---|
1282 | check_vec2 = isStrSubset( vec2,","+vNam(varn)+",") |
---|
1283 | end if |
---|
1284 | |
---|
1285 | 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 |
---|
1286 | check = False |
---|
1287 | else |
---|
1288 | check = True |
---|
1289 | end if |
---|
1290 | |
---|
1291 | if(check) then |
---|
1292 | |
---|
1293 | no_var=no_var+1 |
---|
1294 | |
---|
1295 | if (vector .EQ. 1) then |
---|
1296 | if (","+vNam(varn)+"," .EQ. vec1) then |
---|
1297 | v1=v1+1 |
---|
1298 | end if |
---|
1299 | if (","+vNam(varn)+"," .EQ. vec2) then |
---|
1300 | v2=v2+1 |
---|
1301 | end if |
---|
1302 | end if |
---|
1303 | |
---|
1304 | if (xyc .EQ. 1) then |
---|
1305 | temp = f[:]->$vNam(varn)$ |
---|
1306 | data_att = f_att->$vNam(varn)$ |
---|
1307 | if(vNam(varn) .eq. "ts_xy" .or. vNam(varn) .eq. "us_xy" .or. vNam(varn) .eq. "z0s_xy") |
---|
1308 | ;these variables depend von zu1_xy and that's why they have only one z-layer |
---|
1309 | data(varn,:,0,:,:)=temp(:,0,ys:ye,xs:xe) |
---|
1310 | else |
---|
1311 | data(varn,:,:,:,:)=temp(:,zs:ze,ys:ye,xs:xe) |
---|
1312 | end if |
---|
1313 | delete(temp) |
---|
1314 | end if |
---|
1315 | if ( xzc .eq. 1 ) then |
---|
1316 | temp = f[:]->$vNam(varn)$ |
---|
1317 | data_att = f_att->$vNam(varn)$ |
---|
1318 | data(varn,:,:,:,:)=temp(:,zs:ze,ys:ye,xs:xe) |
---|
1319 | end if |
---|
1320 | if ( yzc .eq. 1) then |
---|
1321 | temp = f[:]->$vNam(varn)$ |
---|
1322 | data_att = f_att->$vNam(varn)$ |
---|
1323 | data(varn,:,:,:,:)=temp(:,zs:ze,ys:ye,xs:xe) |
---|
1324 | end if |
---|
1325 | if (check_vec1) then |
---|
1326 | vect1=data(varn,:,:,:,:) |
---|
1327 | end if |
---|
1328 | if (check_vec2) then |
---|
1329 | vect2=data(varn,:,:,:,:) |
---|
1330 | end if |
---|
1331 | |
---|
1332 | data!0 = "var" |
---|
1333 | data!1 = "t" |
---|
1334 | data!2 = "z" |
---|
1335 | data!3 = "y" |
---|
1336 | data!4 = "x" |
---|
1337 | |
---|
1338 | MinVal(varn) = min(data(varn,:,:,:,:)) |
---|
1339 | MaxVal(varn) = max(data(varn,:,:,:,:)) |
---|
1340 | |
---|
1341 | unit(varn) = data_att@units |
---|
1342 | delete(data_att) |
---|
1343 | |
---|
1344 | end if |
---|
1345 | |
---|
1346 | end do |
---|
1347 | |
---|
1348 | var_input=new(no_var,string) |
---|
1349 | numb=0 |
---|
1350 | no_var=0 |
---|
1351 | |
---|
1352 | do varn=dim-1,0,1 |
---|
1353 | |
---|
1354 | if (vector .EQ. 1) then |
---|
1355 | check_vec1 = isStrSubset( vec1,","+vNam(varn)+",") |
---|
1356 | check_vec2 = isStrSubset( vec2,","+vNam(varn)+",") |
---|
1357 | end if |
---|
1358 | |
---|
1359 | 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 |
---|
1360 | check = False |
---|
1361 | else |
---|
1362 | check = True |
---|
1363 | end if |
---|
1364 | |
---|
1365 | if(check) then |
---|
1366 | var_input(numb)=vNam(varn) |
---|
1367 | numb=numb+1 |
---|
1368 | end if |
---|
1369 | |
---|
1370 | if (var .NE. "all") then |
---|
1371 | check = isStrSubset( var,","+vNam(varn)+"," ) |
---|
1372 | end if |
---|
1373 | |
---|
1374 | if (check) then |
---|
1375 | no_var = no_var+1 |
---|
1376 | end if |
---|
1377 | |
---|
1378 | end do |
---|
1379 | |
---|
1380 | if (no_var .EQ. 0) then |
---|
1381 | print(" ") |
---|
1382 | print("The variables var='"+var+"' do not exist on your input file;") |
---|
1383 | print("be sure to have one comma berfore and after each variable") |
---|
1384 | print(" ") |
---|
1385 | exit |
---|
1386 | end if |
---|
1387 | |
---|
1388 | if (vector .EQ. 1) then |
---|
1389 | if (v1 .EQ. 0)then |
---|
1390 | print(" ") |
---|
1391 | print("Component 1 for the vector-plot ('vec1') must be one of the varibles on the input file:") |
---|
1392 | print("- "+var_input) |
---|
1393 | print("be sure to have one comma berfore and after the variable") |
---|
1394 | print(" ") |
---|
1395 | exit |
---|
1396 | end if |
---|
1397 | |
---|
1398 | if (v2 .EQ. 0)then |
---|
1399 | print(" ") |
---|
1400 | print("Component 2 for the vector-plot ('vec2') must be one of the varibles on the input file:") |
---|
1401 | print("- "+var_input) |
---|
1402 | print("be sure to have one comma berfore and after the variable") |
---|
1403 | print(" ") |
---|
1404 | exit |
---|
1405 | end if |
---|
1406 | end if |
---|
1407 | |
---|
1408 | ; *************************************************** |
---|
1409 | ; open workstation(s) |
---|
1410 | ; *************************************************** |
---|
1411 | |
---|
1412 | wks_ps = gsn_open_wks(format_out,file_out) |
---|
1413 | gsn_define_colormap(wks_ps,"rainbow+white") |
---|
1414 | |
---|
1415 | if (vector .EQ. 1 .AND. plotvec .EQ. "plotvec") then |
---|
1416 | plot=new((/no_time*no_layer*no_var+no_time*no_layer/),graphic) |
---|
1417 | else |
---|
1418 | plot=new((/no_time*no_layer*no_var/),graphic) |
---|
1419 | end if |
---|
1420 | |
---|
1421 | page = 0 |
---|
1422 | n=0 |
---|
1423 | print(" ") |
---|
1424 | print("Plots sorted by '"+sort+"'") |
---|
1425 | print(" ") |
---|
1426 | |
---|
1427 | ; *************************************************** |
---|
1428 | ; create plots |
---|
1429 | ; *************************************************** |
---|
1430 | |
---|
1431 | if (vector .EQ. 1 .AND. plotvec .EQ. "plotvec") then |
---|
1432 | do lo = los, loe |
---|
1433 | do li = lis, lie |
---|
1434 | if (xyc .EQ. 1)then |
---|
1435 | if (sort .EQ. "time")then |
---|
1436 | level = "z=" + z_d(lo) + "m" |
---|
1437 | else |
---|
1438 | level = "z=" + z_d(li) + "m" |
---|
1439 | end if |
---|
1440 | end if |
---|
1441 | if (xzc .EQ. 1)then |
---|
1442 | if (sort .EQ. "time")then |
---|
1443 | level = "y=" + y_d(lo) + "m" |
---|
1444 | else |
---|
1445 | level = "y=" + y_d(li) + "m" |
---|
1446 | end if |
---|
1447 | end if |
---|
1448 | if (yzc .EQ. 1)then |
---|
1449 | if (sort .EQ. "time")then |
---|
1450 | level = "x=" + x_d(lo) + "m" |
---|
1451 | else |
---|
1452 | level = "x=" + x_d(li) + "m" |
---|
1453 | end if |
---|
1454 | end if |
---|
1455 | vecres = True ; vector only resources |
---|
1456 | vecres@gsnDraw = False ; don't draw |
---|
1457 | vecres@gsnFrame = False ; don't advance frame |
---|
1458 | vecres@vcGlyphStyle = "CurlyVector" ; curly vectors |
---|
1459 | vecres@vcRefMagnitudeF = ref_mag ; define vector ref mag |
---|
1460 | vecres@vcRefLengthF = 0.05 ; define length of vec ref |
---|
1461 | vecres@gsnRightString = "Vector plot of "+vec1+" and "+vec2 |
---|
1462 | vecres@tmXBLabelFontHeightF = font_size |
---|
1463 | vecres@tmYLLabelFontHeightF = font_size |
---|
1464 | vecres@tiXAxisFontHeightF = font_size |
---|
1465 | vecres@tiYAxisFontHeightF = font_size |
---|
1466 | if (sort .EQ. "time")then |
---|
1467 | vecres@gsnLeftString = "t=" + decimalPlaces(t_all(li)/3600,2,True) +"h "+level |
---|
1468 | else |
---|
1469 | vecres@gsnLeftString = "t=" + decimalPlaces(t_all(lo)/3600,2,True) +"h "+level |
---|
1470 | end if |
---|
1471 | vecres@tiXAxisString = " " |
---|
1472 | if (xyc .EQ. 1)then |
---|
1473 | if (sort .EQ. "time")then |
---|
1474 | plot(n) = gsn_csm_vector(wks_ps,vect1(li,lo,:,:),vect2(li,lo,:,:),vecres) |
---|
1475 | else |
---|
1476 | plot(n) = gsn_csm_vector(wks_ps,vect1(lo,li,:,:),vect2(lo,li,:,:),vecres) |
---|
1477 | end if |
---|
1478 | end if |
---|
1479 | if (xzc .EQ. 1) then |
---|
1480 | if (sort .EQ. "time")then |
---|
1481 | plot(n) = gsn_csm_vector(wks_ps,vect1(li,:,lo,:),vect2(li,:,lo,:),vecres) |
---|
1482 | else |
---|
1483 | plot(n) = gsn_csm_vector(wks_ps,vect1(lo,:,li,:),vect2(lo,:,li,:),vecres) |
---|
1484 | end if |
---|
1485 | end if |
---|
1486 | if (yzc .EQ. 1) then |
---|
1487 | if (sort .EQ. "time")then |
---|
1488 | plot(n) = gsn_csm_vector(wks_ps,vect1(li,:,:,lo),vect2(li,:,:,lo),vecres) |
---|
1489 | else |
---|
1490 | plot(n) = gsn_csm_vector(wks_ps,vect1(lo,:,:,li),vect2(lo,:,:,li),vecres) |
---|
1491 | end if |
---|
1492 | end if |
---|
1493 | n=n+1 |
---|
1494 | end do |
---|
1495 | end do |
---|
1496 | end if |
---|
1497 | |
---|
1498 | do varn=dim-1,0,1 |
---|
1499 | |
---|
1500 | if (vector .EQ. 1) then |
---|
1501 | check_vecp = isStrSubset( plotvec,","+vNam(varn)+",") |
---|
1502 | end if |
---|
1503 | |
---|
1504 | 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 |
---|
1505 | check = False |
---|
1506 | else |
---|
1507 | check = True |
---|
1508 | end if |
---|
1509 | |
---|
1510 | if (var .NE. "all") then |
---|
1511 | check = isStrSubset( var,","+vNam(varn)+"," ) |
---|
1512 | end if |
---|
1513 | |
---|
1514 | if(check) then |
---|
1515 | |
---|
1516 | space=(MaxVal(varn)-MinVal(varn))/24 |
---|
1517 | |
---|
1518 | cs_res@cnMinLevelValF = MinVal(varn) |
---|
1519 | cs_res@cnMaxLevelValF = MaxVal(varn) |
---|
1520 | |
---|
1521 | cs_res@cnLevelSpacingF = space |
---|
1522 | |
---|
1523 | ; **************************************************** |
---|
1524 | ; loops over time and layer |
---|
1525 | ; **************************************************** |
---|
1526 | |
---|
1527 | |
---|
1528 | do lo = los, loe |
---|
1529 | do li = lis, lie |
---|
1530 | |
---|
1531 | ; **************************************************** |
---|
1532 | ; xy cross section |
---|
1533 | ; **************************************************** |
---|
1534 | |
---|
1535 | if ( xyc .eq. 1 ) then |
---|
1536 | |
---|
1537 | cs_res@gsnLeftString = "Plot of "+vNam(varn) |
---|
1538 | cs_res@gsnRightString = unit(varn) |
---|
1539 | |
---|
1540 | if ( sort .eq. "time" ) then |
---|
1541 | if ( z_d(zs+lo) .eq. -1)then |
---|
1542 | if (delta_z .EQ. -1) then |
---|
1543 | level = "-average" |
---|
1544 | else |
---|
1545 | level = "=" + z_d(zs+lo) + "m" |
---|
1546 | end if |
---|
1547 | else |
---|
1548 | level = "=" + z_d(zs+lo) + "m" |
---|
1549 | end if |
---|
1550 | cs_res@gsnCenterString = "t=" + decimalPlaces(t_all(li)/3600,2,True) +"h z"+level |
---|
1551 | plot(n) = gsn_csm_contour(wks_ps,data(varn,li,lo,:,:),cs_res) |
---|
1552 | if (vector .EQ. 1 .AND. check_vecp) then |
---|
1553 | vecres = True ; vector only resources |
---|
1554 | vecres@gsnDraw = False ; don't draw |
---|
1555 | vecres@gsnFrame = False ; don't advance frame |
---|
1556 | vecres@vcGlyphStyle = "CurlyVector" ; curly vectors |
---|
1557 | vecres@vcRefMagnitudeF = ref_mag ; define vector ref mag |
---|
1558 | vecres@vcRefLengthF = 0.05 ; define length of vec ref |
---|
1559 | vecres@gsnRightString = " " |
---|
1560 | vecres@gsnLeftString = " " |
---|
1561 | vecres@tiXAxisString = " " |
---|
1562 | plot_vec=gsn_csm_vector(wks_ps,vect1(li,lo,:,:),vect2(li,lo,:,:),vecres) |
---|
1563 | overlay(plot(n), plot_vec) |
---|
1564 | end if |
---|
1565 | end if |
---|
1566 | |
---|
1567 | if ( sort .eq. "layer" ) then |
---|
1568 | if (z_d(zs+li) .eq. -1) then |
---|
1569 | if (delta_z .EQ. -1) then |
---|
1570 | level = "-average" |
---|
1571 | else |
---|
1572 | level = "=" + z_d(zs+li) + "m" |
---|
1573 | end if |
---|
1574 | else |
---|
1575 | level = "=" + z_d(zs+li) + "m" |
---|
1576 | end if |
---|
1577 | cs_res@gsnCenterString = "t=" + decimalPlaces(t_all(lo)/3600,2,True) + "h z"+ level |
---|
1578 | plot(n) = gsn_csm_contour(wks_ps,data(varn,lo,li,:,:),cs_res) |
---|
1579 | if (vector .EQ. 1 .AND. check_vecp) then |
---|
1580 | vecres = True ; vector only resources |
---|
1581 | vecres@gsnDraw = False ; don't draw |
---|
1582 | vecres@gsnFrame = False ; don't advance frame |
---|
1583 | vecres@vcGlyphStyle = "CurlyVector" ; curly vectors |
---|
1584 | vecres@vcRefMagnitudeF = ref_mag ; define vector ref mag |
---|
1585 | vecres@vcRefLengthF = 0.05 ; define length of vec ref |
---|
1586 | vecres@gsnRightString = " " ; turn off right string |
---|
1587 | vecres@gsnLeftString = " " ; turn off left string |
---|
1588 | vecres@tiXAxisString = " " |
---|
1589 | plot_vec=gsn_csm_vector(wks_ps,vect1(lo,li,:,:),vect2(lo,li,:,:),vecres) |
---|
1590 | overlay(plot(n), plot_vec) |
---|
1591 | end if |
---|
1592 | end if |
---|
1593 | end if |
---|
1594 | |
---|
1595 | ; **************************************************** |
---|
1596 | ; xz cross section |
---|
1597 | ; **************************************************** |
---|
1598 | |
---|
1599 | if ( xzc .eq. 1 ) then |
---|
1600 | |
---|
1601 | cs_res@gsnLeftString = "Plot of "+vNam(varn) |
---|
1602 | |
---|
1603 | if ( sort .eq. "time" ) then |
---|
1604 | if ( y_d(ys+lo) .eq. -1 ) then |
---|
1605 | level = "-average" |
---|
1606 | else |
---|
1607 | level = "=" + y_d(ys+lo) + "m" |
---|
1608 | end if |
---|
1609 | cs_res@gsnCenterString = "t=" + decimalPlaces(t_all(li)/3600,2,True) + "h y"+ level |
---|
1610 | plot(n) = gsn_csm_contour(wks_ps,data(varn,li,:,lo,:),cs_res) |
---|
1611 | if (vector .EQ. 1 .AND. check_vecp) then |
---|
1612 | vecres = True ; vector only resources |
---|
1613 | vecres@gsnDraw = False ; don't draw |
---|
1614 | vecres@gsnFrame = False ; don't advance frame |
---|
1615 | vecres@vcGlyphStyle = "CurlyVector" ; curly vectors |
---|
1616 | vecres@vcRefMagnitudeF = ref_mag ; define vector ref mag |
---|
1617 | vecres@vcRefLengthF = 0.05 ; define length of vec ref |
---|
1618 | vecres@gsnRightString = " " ; turn off right string |
---|
1619 | vecres@gsnLeftString = " " ; turn off left string |
---|
1620 | vecres@tiXAxisString = " " |
---|
1621 | plot_vec=gsn_csm_vector(wks_ps,vect1(li,:,lo,:),vect2(li,:,lo,:),vecres) |
---|
1622 | overlay(plot(n), plot_vec) |
---|
1623 | end if |
---|
1624 | end if |
---|
1625 | |
---|
1626 | if ( sort .eq. "layer" ) then |
---|
1627 | if ( y_d(ys+li) .eq. -1 ) then |
---|
1628 | level = "-average" |
---|
1629 | else |
---|
1630 | level = "=" + y_d(ys+li) + "m" |
---|
1631 | end if |
---|
1632 | cs_res@gsnCenterString = "t=" + decimalPlaces(t_all(lo)/3600,2,True) + "h y"+ level |
---|
1633 | plot(n) = gsn_csm_contour(wks_ps,data(varn,lo,:,li,:),cs_res) |
---|
1634 | if (vector .EQ. 1 .AND. check_vecp) then |
---|
1635 | vecres = True ; vector only resources |
---|
1636 | vecres@gsnDraw = False ; don't draw |
---|
1637 | vecres@gsnFrame = False ; don't advance frame |
---|
1638 | vecres@vcGlyphStyle = "CurlyVector" ; curly vectors |
---|
1639 | vecres@vcRefMagnitudeF = ref_mag ; define vector ref mag |
---|
1640 | vecres@vcRefLengthF = 0.05 ; define length of vec ref |
---|
1641 | vecres@gsnRightString = " " ; turn off right string |
---|
1642 | vecres@gsnLeftString = " " ; turn off left string |
---|
1643 | vecres@tiXAxisString = " " |
---|
1644 | plot_vec=gsn_csm_vector(wks_ps,vect1(lo,:,li,:),vect2(lo,:,li,:),vecres) |
---|
1645 | overlay(plot(n), plot_vec) |
---|
1646 | end if |
---|
1647 | end if |
---|
1648 | end if |
---|
1649 | |
---|
1650 | ; **************************************************** |
---|
1651 | ; yz cross section |
---|
1652 | ; **************************************************** |
---|
1653 | |
---|
1654 | if ( yzc .eq. 1 ) then |
---|
1655 | |
---|
1656 | cs_res@gsnLeftString = "Plot of "+vNam(varn) |
---|
1657 | |
---|
1658 | if ( sort .eq. "time" ) then |
---|
1659 | if ( x_d(xs+lo) .eq. -1 ) then |
---|
1660 | level = "-average" |
---|
1661 | else |
---|
1662 | level = "=" + x_d(xs+lo) + "m" |
---|
1663 | end if |
---|
1664 | cs_res@gsnCenterString = "t=" + decimalPlaces(t_all(li)/3600,2,True) + "h x"+ level |
---|
1665 | plot(n) = gsn_csm_contour(wks_ps,data(varn,li,:,:,lo),cs_res) |
---|
1666 | if (vector .EQ. 1 .AND. check_vecp) then |
---|
1667 | vecres = True ; vector only resources |
---|
1668 | vecres@gsnDraw = False ; don't draw |
---|
1669 | vecres@gsnFrame = False ; don't advance frame |
---|
1670 | vecres@vcGlyphStyle = "CurlyVector" ; curly vectors |
---|
1671 | vecres@vcRefMagnitudeF = ref_mag ; define vector ref mag |
---|
1672 | vecres@vcRefLengthF = 0.05 ; define length of vec ref |
---|
1673 | vecres@gsnRightString = " " ; turn off right string |
---|
1674 | vecres@gsnLeftString = " " ; turn off left string |
---|
1675 | vecres@tiXAxisString = " " |
---|
1676 | plot_vec=gsn_csm_vector(wks_ps,vect1(li,:,:,lo),vect2(li,:,:,lo),vecres) |
---|
1677 | overlay(plot(n), plot_vec) |
---|
1678 | end if |
---|
1679 | end if |
---|
1680 | |
---|
1681 | if ( sort .eq. "layer" ) then |
---|
1682 | if ( x_d(xs+li) .eq. -1 ) then |
---|
1683 | level = "-average" |
---|
1684 | else |
---|
1685 | level = "=" + x_d(xs+li) + "m" |
---|
1686 | end if |
---|
1687 | cs_res@gsnCenterString = "t=" + decimalPlaces(t_all(lo)/3600,2,True) + "h x"+ level |
---|
1688 | plot(n) = gsn_csm_contour(wks_ps,data(varn,lo,:,:,li),cs_res) |
---|
1689 | if (vector .EQ. 1 .AND. check_vecp)then |
---|
1690 | vecres = True ; vector only resources |
---|
1691 | vecres@gsnDraw = False ; don't draw |
---|
1692 | vecres@gsnFrame = False ; don't advance frame |
---|
1693 | vecres@vcGlyphStyle = "CurlyVector" ; curly vectors |
---|
1694 | vecres@vcRefMagnitudeF = ref_mag ; define vector ref mag |
---|
1695 | vecres@vcRefLengthF = 0.05 ; define length of vec ref |
---|
1696 | vecres@gsnRightString = " " ; turn off right string |
---|
1697 | vecres@gsnLeftString = " " ; turn off left string |
---|
1698 | vecres@tiXAxisString = " " |
---|
1699 | plot_vec=gsn_csm_vector(wks_ps,vect1(lo,:,:,li),vect2(lo,:,:,li),vecres) |
---|
1700 | overlay(plot(n), plot_vec) |
---|
1701 | end if |
---|
1702 | end if |
---|
1703 | end if |
---|
1704 | n=n+1 |
---|
1705 | end do |
---|
1706 | end do |
---|
1707 | end if |
---|
1708 | end do |
---|
1709 | |
---|
1710 | ; *************************************************** |
---|
1711 | ; merge plots onto one page |
---|
1712 | ; *************************************************** |
---|
1713 | |
---|
1714 | if (vector .EQ. 1 .AND. plotvec .EQ. "plotvec")then |
---|
1715 | if (format_out .EQ. "eps" .OR. format_out .EQ. "epsi") then |
---|
1716 | gsn_panel(wks_ps,plot(0:(no_time*no_layer*(no_var+1))-1),(/no_var+1,no_layer*no_time/),cs_resP) |
---|
1717 | print(" ") |
---|
1718 | print("Outputs to .eps or .epsi have only one frame") |
---|
1719 | print(" ") |
---|
1720 | else |
---|
1721 | do np = 0,no_layer*no_time*(no_var+1)-1,no_rows*no_columns |
---|
1722 | if ( np + no_rows*no_columns .gt. (no_layer*no_time*(no_var+1))-1) then |
---|
1723 | gsn_panel(wks_ps, plot(np:(no_layer*no_time*(no_var+1))-1),(/no_rows,no_columns/),cs_resP) |
---|
1724 | else |
---|
1725 | gsn_panel(wks_ps, plot(np:np+no_rows*no_columns-1),(/no_rows,no_columns/),cs_resP) |
---|
1726 | end if |
---|
1727 | end do |
---|
1728 | end if |
---|
1729 | else |
---|
1730 | if (format_out .EQ. "eps" .OR. format_out .EQ. "epsi") then |
---|
1731 | gsn_panel(wks_ps,plot(0:(no_time*no_layer*no_var)-1),(/no_var,no_layer*no_time/),cs_resP) |
---|
1732 | print(" ") |
---|
1733 | print("Outputs to .eps or .epsi have only one frame") |
---|
1734 | print(" ") |
---|
1735 | else |
---|
1736 | do np = 0,no_layer*no_time*no_var-1,no_rows*no_columns |
---|
1737 | if ( np + no_rows*no_columns .gt. (no_layer*no_time*no_var)-1) then |
---|
1738 | gsn_panel(wks_ps, plot(np:(no_layer*no_time*no_var)-1),(/no_rows,no_columns/),cs_resP) |
---|
1739 | else |
---|
1740 | gsn_panel(wks_ps, plot(np:np+no_rows*no_columns-1),(/no_rows,no_columns/),cs_resP) |
---|
1741 | end if |
---|
1742 | end do |
---|
1743 | end if |
---|
1744 | end if |
---|
1745 | |
---|
1746 | print(" ") |
---|
1747 | print("Output to: " + file_out +"."+ format_out) |
---|
1748 | print(" ") |
---|
1749 | |
---|
1750 | end |
---|