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