| 1 | pro model_levels, plotlev |
|---|
| 2 | ;; first is one (as in gw.def) |
|---|
| 3 | |
|---|
| 4 | ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; |
|---|
| 5 | ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; |
|---|
| 6 | file = './zefolder/' |
|---|
| 7 | file = file + 'wrfout_d01_2024-03-21_00:00:00' |
|---|
| 8 | |
|---|
| 9 | file = './wrfout_d01_9999-09-09_09:00:00' |
|---|
| 10 | ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; |
|---|
| 11 | ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; |
|---|
| 12 | |
|---|
| 13 | ;; |
|---|
| 14 | ;; ini |
|---|
| 15 | ;; |
|---|
| 16 | ;if (n_elements(plotlev) eq 0) then plotlev=7 |
|---|
| 17 | |
|---|
| 18 | ;set_plot, 'ps' |
|---|
| 19 | ;device, filename='model_levels.ps' |
|---|
| 20 | |
|---|
| 21 | ; |
|---|
| 22 | ; constantes |
|---|
| 23 | ; |
|---|
| 24 | grav = 3.72 |
|---|
| 25 | |
|---|
| 26 | ; |
|---|
| 27 | ; retrieve model geopotential |
|---|
| 28 | ; |
|---|
| 29 | openr,unit,'heights.dat',/get_lun,error=err |
|---|
| 30 | ;IF (err ne 0) THEN BEGIN |
|---|
| 31 | getcdf, $ |
|---|
| 32 | file=file, $ |
|---|
| 33 | charvar='PHTOT', $ |
|---|
| 34 | invar=geop |
|---|
| 35 | getcdf, $ |
|---|
| 36 | file=file, $ |
|---|
| 37 | charvar='HGT', $ |
|---|
| 38 | invar=hgt |
|---|
| 39 | ; save, geop, filename='heights.dat' |
|---|
| 40 | ;ENDIF ELSE BEGIN |
|---|
| 41 | ; restore, filename='heights.dat' |
|---|
| 42 | ;ENDELSE |
|---|
| 43 | |
|---|
| 44 | ; |
|---|
| 45 | ; size |
|---|
| 46 | ; |
|---|
| 47 | s = size(geop) |
|---|
| 48 | nx = s[1] ;& print, nx |
|---|
| 49 | ny = s[2] ;& print, ny |
|---|
| 50 | nz = s[3] ;& print, nz |
|---|
| 51 | nt = s[4] ;& print, nt |
|---|
| 52 | |
|---|
| 53 | if (n_elements(plotlev) eq 0) then plotlev=nz-1 |
|---|
| 54 | |
|---|
| 55 | for lev = 1, nz-1 do begin |
|---|
| 56 | |
|---|
| 57 | print, '*************************************************************************' |
|---|
| 58 | print, 'model level ', lev |
|---|
| 59 | |
|---|
| 60 | ; |
|---|
| 61 | ; heights of model levels |
|---|
| 62 | ; |
|---|
| 63 | ;height = (geop(*,*,lev-1,*) - geop(*,*,0.,*)) / grav |
|---|
| 64 | height = fltarr(nx,ny,nt) |
|---|
| 65 | for i=0,nt-1 do height(*,*,i) = geop(*,*,lev-1,i) / grav - hgt(*,*) |
|---|
| 66 | |
|---|
| 67 | ; |
|---|
| 68 | ; spatial mean (temporal variations) |
|---|
| 69 | ; |
|---|
| 70 | yes = TOTAL(TOTAL(height,1),1) / float(nx) / float(ny) |
|---|
| 71 | print, 'spat. mean (temp. var.) ', mean(yes), max(yes), min(yes);, max(yes) - min(yes) |
|---|
| 72 | |
|---|
| 73 | ; |
|---|
| 74 | ; temporal mean (spatial variations) |
|---|
| 75 | ; |
|---|
| 76 | yet = TOTAL(height,3) / float(nt) |
|---|
| 77 | print, 'temp. mean (spat. var.) ', mean(yet), max(yet), min(yet);, max(yet) - min(yet) |
|---|
| 78 | |
|---|
| 79 | ; |
|---|
| 80 | ; complete mean |
|---|
| 81 | ; |
|---|
| 82 | ye = TOTAL(TOTAL(TOTAL(height,1),1),1) / float(nx) / float(ny) / float(nt) |
|---|
| 83 | print, '******* mean ', ye |
|---|
| 84 | |
|---|
| 85 | if (lev eq plotlev) then begin |
|---|
| 86 | !p.multi=[0,2,1] |
|---|
| 87 | contour, yet, /cell_fill, nlevels=30, title=string(min(yet))+string(max(yet)) |
|---|
| 88 | plot, yes, yrange=[min(yes),max(yes)] |
|---|
| 89 | stop |
|---|
| 90 | endif |
|---|
| 91 | |
|---|
| 92 | endfor |
|---|
| 93 | |
|---|
| 94 | ;device, /close |
|---|
| 95 | end |
|---|
