1 | pro map_uvt |
---|
2 | ; |
---|
3 | ; |
---|
4 | ; |
---|
5 | what_I_plot=0. |
---|
6 | overcontour=0 |
---|
7 | no3d='true' |
---|
8 | @map_uvt_inc.pro |
---|
9 | SPAWN, '\rm param_plot.idl ; cp map_uvt_inc.pro param_plot.idl ; cp -f map_uvt_inc.pro '+save_ps+'.map_uvt_inc.pro' |
---|
10 | if (n_elements(coord2d) eq 0) then coord2d='false' |
---|
11 | if (n_elements(overvector_x) eq 0) then overvector_x = 1 ;; initialisation (!! map_uvt_inc.pro utilise dans map_uvt et map_latlon) |
---|
12 | ; |
---|
13 | ; |
---|
14 | ; |
---|
15 | if (n_elements(field1) ne 0) then getcdf, file=filename, charvar=field1, invar=cfield1 |
---|
16 | if (n_elements(field2) ne 0) then getcdf, file=filename, charvar=field2, invar=cfield2 |
---|
17 | if ( ( overvector_x ne 0 ) or ( n_elements(field1) eq 0 ) ) then begin |
---|
18 | u = getget(filename, 'Um', count=[0,0,1,0], offset=[0,0,nlevel,0]) |
---|
19 | v = getget(filename, 'Vm', count=[0,0,1,0], offset=[0,0,nlevel,0]) |
---|
20 | help, u |
---|
21 | endif |
---|
22 | getcdf, file=filename, charvar='XLONG', invar=longi |
---|
23 | getcdf, file=filename, charvar='XLAT', invar=lati |
---|
24 | getcdf, file=filename, charvar='HGT', invar=hgt |
---|
25 | ; |
---|
26 | ; |
---|
27 | ; |
---|
28 | sp = 5 ;; relaxation width |
---|
29 | nx = n_elements(longi(0,*)) |
---|
30 | ny = n_elements(longi(*,0)) |
---|
31 | if (n_elements(field1) ne 0) then begin |
---|
32 | if (no3d ne 'true') then cfield1 = cfield1 (sp:nx-sp-1,sp:ny-sp-1,*,*) else cfield1 = cfield1 (sp:nx-sp-1,sp:ny-sp-1,*) |
---|
33 | endif |
---|
34 | if (n_elements(field2) ne 0) then begin |
---|
35 | if (no3d ne 'true') then cfield2 = cfield2 (sp:nx-sp-1,sp:ny-sp-1,*,*) else cfield2 = cfield2 (sp:nx-sp-1,sp:ny-sp-1,*) |
---|
36 | endif |
---|
37 | if ( ( overvector_x ne 0 ) or ( n_elements(field1) eq 0 ) ) then begin |
---|
38 | u = reform (u (sp:nx-sp-1,sp:ny-sp-1,0,*) ) |
---|
39 | v = reform (v (sp:nx-sp-1,sp:ny-sp-1,0,*) ) |
---|
40 | endif |
---|
41 | longi = longi (sp:nx-sp-1,sp:ny-sp-1,*) |
---|
42 | lati = lati (sp:nx-sp-1,sp:ny-sp-1,*) |
---|
43 | hgt = hgt (sp:nx-sp-1,sp:ny-sp-1,*) |
---|
44 | nx = n_elements(longi(*,0,0)) |
---|
45 | ny = n_elements(longi(0,*,0)) |
---|
46 | nt = n_elements(longi(0,0,*)) |
---|
47 | if ( overvector_x ne 0 ) then begin |
---|
48 | overvector_x = u |
---|
49 | overvector_y = v |
---|
50 | endif |
---|
51 | ; |
---|
52 | ; |
---|
53 | ; |
---|
54 | if (ntime eq 99) then begin |
---|
55 | PRINT, '-- ALL TIME STEPS', nt-1 |
---|
56 | ntstart = 0 & ntend = nt-1 |
---|
57 | endif else begin |
---|
58 | PRINT, '-- ONLY TIME STEP ', string(ntime,'(I0)') |
---|
59 | ntstart = ntime & ntend = ntime |
---|
60 | endelse |
---|
61 | for ntime = ntstart,ntend do begin |
---|
62 | help, u |
---|
63 | ; |
---|
64 | ; |
---|
65 | ; |
---|
66 | zefile=save_ps+string(100+ntime,'(I0)') |
---|
67 | PS_Start, filename=zefile+'.ps' |
---|
68 | print, zefile+'.ps' |
---|
69 | !P.Charsize = 1.2 |
---|
70 | !p.charthick = 2.0 |
---|
71 | !p.thick = 2.0 |
---|
72 | !x.thick = 2.0 |
---|
73 | !y.thick = 2.0 |
---|
74 | ; |
---|
75 | ; |
---|
76 | ; |
---|
77 | lon = reform(longi(*,*,ntime)) |
---|
78 | lat = reform(lati(*,*,ntime)) |
---|
79 | if (n_elements(field1) eq 0) then begin |
---|
80 | print, 'field1: horizontal velocity' |
---|
81 | zevel = overvector_x^2 + overvector_y^2 ;; attention il faut que les tableaux soient de la meme taille, OK si uvmet |
---|
82 | what_I_plot = sqrt(zevel) |
---|
83 | endif else begin |
---|
84 | if (no3d ne 'true') then what_I_plot = reform(cfield1(*,*,nlevel,ntime)) else what_I_plot = reform(cfield1(*,*,ntime)) |
---|
85 | if (n_elements(u) ne 0) then overvector_x = reform(u(*,*,ntime)) ;; ne pas utiliser test overvector_x a cause de la boucle temps |
---|
86 | if (n_elements(v) ne 0) then overvector_y = reform(v(*,*,ntime)) ;; ne pas utiliser test overvector_y a cause de la boucle temps |
---|
87 | endelse |
---|
88 | if (n_elements(field2) eq 0) then begin |
---|
89 | print, 'field2: topography' |
---|
90 | overcontour = reform(hgt(*,*,ntime)) |
---|
91 | endif else begin |
---|
92 | if (no3d ne 'true') then overcontour = reform(cfield2(*,*,nlevel,ntime)) else overcontour = reform(cfield2(*,*,ntime)) |
---|
93 | endelse |
---|
94 | ; |
---|
95 | ; |
---|
96 | ; |
---|
97 | if ( coord2d eq 'polar' ) then begin |
---|
98 | print, 'OK YOU USE MAP_SET with POLAR PROJECTION. VECTORS ARE NOT SUPPORTED. USE polar_uv OR ADAPT THIS SCRIPT.' |
---|
99 | overvector_x = 0 |
---|
100 | overvector_y = 0 |
---|
101 | if (n_elements(windowx) ne 0) then begin |
---|
102 | latmin = windowy(0) & latmax = windowy(1) & lonmin = windowx(0) & lonmax = windowx(1) |
---|
103 | endif else begin |
---|
104 | latmin = 65. & latmax = 90. & lonmin = -180. & lonmax = 180. |
---|
105 | endelse |
---|
106 | print, 'latmin,lonmin,latmax,lonmax', latmin,lonmin,latmax,lonmax |
---|
107 | map_set, latmax, 0., /isotropic, /azimuthal, /noborder, limit=[latmin,lonmin,latmax,lonmax], title=title_user, /advance |
---|
108 | endif else begin |
---|
109 | minlat=min(lat) & maxlat=max(lat) & minlon=min(lon) & maxlon=max(lon) |
---|
110 | if (coord2d eq 'true') then begin |
---|
111 | ;;;; CECI EST DESORMAIS FAIT DANS MAP_LATLON |
---|
112 | ;;npoints=n_elements(lon(*,0)) + n_elements(lon(0,*)) ;; trop de points, mais au moins on ne perd rien |
---|
113 | ;;TRIANGULATE, lon, lat, tr |
---|
114 | ;;what_I_plot = GRIDDATA( lon, lat, what_I_plot, /LINEAR, triangles=tr, dimension=npoints, MISSING=!VALUES.F_NAN ) |
---|
115 | ;;if ( overvector_x ne 0 ) then begin |
---|
116 | ;; overvector_x = GRIDDATA( lon, lat, overvector_x, /LINEAR, triangles=tr, dimension=npoints, MISSING=!VALUES.F_NAN ) |
---|
117 | ;; overvector_y = GRIDDATA( lon, lat, overvector_y, /LINEAR, triangles=tr, dimension=npoints, MISSING=!VALUES.F_NAN ) |
---|
118 | ;;endif |
---|
119 | ;;overcontour = GRIDDATA( lon, lat, overcontour, /LINEAR, triangles=tr, dimension=npoints, MISSING=!VALUES.F_NAN ) |
---|
120 | ;; if ( overvector_x ne 0 ) then begin |
---|
121 | ;; ; sale sale sale |
---|
122 | ;; if (minlat lt min(lat(*,0))) then overvector_y=-overvector_y |
---|
123 | ;; if (minlon lt min(lon(0,*))) then overvector_x=-overvector_x |
---|
124 | ;; endif |
---|
125 | ;;lon = minlon + (maxlon - minlon)*findgen(npoints)/float(npoints-1) |
---|
126 | ;;lat = minlat + (maxlat - minlat)*findgen(npoints)/float(npoints-1) |
---|
127 | endif else begin |
---|
128 | lon = reform(lon(*,0)) |
---|
129 | lat = reform(lat(0,*)) |
---|
130 | ;;npoints=n_elements(lon(*,0)) + n_elements(lon(0,*)) |
---|
131 | ;;what_I_plot = REBIN( what_I_plot, npoints, npoints ) |
---|
132 | ;;overvector_x = REBIN( overvector_x, npoints, npoints ) |
---|
133 | ;;overvector_y = REBIN( overvector_y, npoints, npoints ) |
---|
134 | ;;overcontour = REBIN( overcontour, npoints, npoints ) |
---|
135 | ;;lon = minlon + (maxlon - minlon)*findgen(npoints)/float(npoints-1) |
---|
136 | ;;lat = minlat + (maxlat - minlat)*findgen(npoints)/float(npoints-1) |
---|
137 | ;lon = findgen( n_elements(what_I_plot(*,0)) ) |
---|
138 | ;lat = findgen( n_elements(what_I_plot(0,*)) ) |
---|
139 | ;print, 'ATTENTION!!!! NOUS SOMMES BIEN D ACCORD QUE VOUS NE CHARGEZ PAS UMET et VMET ??' |
---|
140 | endelse |
---|
141 | ;;;;; trouve dans polar_uv.pro |
---|
142 | ;;if (coord2d eq 'regular') then begin ;; carte avec mercator |
---|
143 | ;; lon = reform(lon(*,0)) |
---|
144 | ;; lat = reform(lat(0,*)) |
---|
145 | ;;endif |
---|
146 | endelse |
---|
147 | help, what_I_plot, lon, lat |
---|
148 | print, min(what_I_plot) |
---|
149 | print, max(what_I_plot) |
---|
150 | print, min(overcontour) |
---|
151 | print, max(overcontour) |
---|
152 | print, min(overvector_x) |
---|
153 | print, min(overvector_y) |
---|
154 | ; |
---|
155 | ; |
---|
156 | ; |
---|
157 | map_latlon, $ |
---|
158 | what_I_plot, $ ; 2D field |
---|
159 | lon, $ ; 1D latitude |
---|
160 | lat, $ ; 1D longitude |
---|
161 | minfield=minfield_init, $ ; minimum value of plotted field (=0: calculate) |
---|
162 | maxfield=maxfield_init, $ ; maximum value of plotted field (=0: calculate) |
---|
163 | overcontour=overcontour, $ ; another 2D field to overplot with contour lines (=0: no) |
---|
164 | overvector_x=overvector_x, $ ; wind vector - x component (=0: no) |
---|
165 | overvector_y=overvector_y, $ ; wind vector - y component (=0: no) |
---|
166 | ct=pal, $ ; color table (33-rainbow is default) |
---|
167 | colors=colors, $ ; number of colors/levels (32 is default) |
---|
168 | title=title_user, $ ; title of the plot ('' is default) |
---|
169 | format=format ; format of colorbar annotations ('(F6.2)' is default) |
---|
170 | ; |
---|
171 | ; |
---|
172 | ; |
---|
173 | if ( coord2d eq 'polar' ) then begin ;;; pourrait aller dans map_latlon |
---|
174 | loadct, 0 |
---|
175 | MAP_GRID, $ |
---|
176 | CHARSIZE = 1., $ |
---|
177 | COLOR = 0, $ |
---|
178 | ;lats=-60, $ |
---|
179 | LABEL = 1, $ ;; /LABEL or LABEL=2 (one label any 2 grid lines) |
---|
180 | LATDEL = intervaly, $ ;;5 10 |
---|
181 | LONDEL = intervalx, $ ;;15 |
---|
182 | LONLAB = latmin + intervaly/2., $ ;5. + (latmin+latmax)/2., $ ;0. |
---|
183 | LATLAB = -0.001, $ |
---|
184 | GLINESTYLE = 2, $ |
---|
185 | GLINETHICK = 0.3 |
---|
186 | ;LONALIGN = 0., $ |
---|
187 | ;LATALIGN = 1. |
---|
188 | endif |
---|
189 | ; |
---|
190 | ; |
---|
191 | ; |
---|
192 | PS_End, /PNG |
---|
193 | endfor |
---|
194 | end |
---|