Context Navigation

map_latlon.pro @ 113

Last change on this file since 113 was 91, checked in by aslmd, 14 years ago
mars: test outliers [dans initracer.F, commente] LMD_MM_MARS: modifications mineures [retrocompatibilite ancienne physique, callphys.def test pour nouvelle physique] PLOT: generalisation de la routine map_uvt pour pouvoir tracer des figures en projection polaire complete
File size: 17.1 KB

Line
1	pro map_latlon, $
2	what_I_plot, $ ; 2D field
3	lon, $ ; 1D latitude OR 2D
4	lat, $ ; 1D longitude OR 2D
5	minfield=minfield_init, $ ; minimum value of plotted field (=0: calculate)
6	maxfield=maxfield_init, $ ; maximum value of plotted field (=0: calculate)
7	overcontour=overcontour, $ ; another 2D field to overplot with contour lines (=0: no)
8	overvector_x=overvector_x, $ ; wind vector - x component (=0: no)
9	overvector_y=overvector_y, $ ; wind vector - y component (=0: no)
10	ct=pal, $ ; color table (33-rainbow is default)
11	colors=colors, $ ; number of colors/levels (32 is default)
12	title=title_user, $ ; title of the plot ('' is default)
13	format=format ; format of colorbar annotations ('(F6.2)' is default)
14
15	;--------------------------------------------
16	;
17	; A general routine to plot lat/lon maps
18	;
19	; USE: - map_latlon, what_I_plot
20	; - map_latlon, what_I_plot, lon, lat
21	; - ... and see options above
22	;
23	;--------------------------------------------
24	; A. Spiga, August 2007
25	;--------------------------------------------
26
27
28
29
30	;---------------------------
31	; a few permanent settings
32	;---------------------------
33	;
34	missing_value=1.e30
35	flag_cb='true'
36	;-- optional:
37	;-- user defined parameters
38	; do not forget to type '.compile map_latlon' to update changes
39	SPAWN, 'touch param_plot.idl'
40	@param_plot.idl
41	;
42	;------------------
43
44
45	;------------------
46	; user settings
47	;------------------
48	if (n_elements(lon) eq 0) then begin
49	lon=findgen(n_elements(what_I_plot(*,0)))
50	endif
51	if (n_elements(lat) eq 0) then lat=findgen(n_elements(what_I_plot(0,*)))
52	if (n_elements(minfield_init) eq 0) then minfield_init=0.
53	if (n_elements(maxfield_init) eq 0) then maxfield_init=0.
54	if (n_elements(overcontour) eq 0) then overcontour=0.
55	if ((n_elements(overvector_x) eq 0) or (n_elements(overvector_y) eq 0)) then begin
56	overvector_x=0.
57	overvector_y=0.
58	endif
59	if (n_elements(colors) eq 0) then colors=32
60	if (n_elements(title_user) eq 0) then title_user=''
61	if (n_elements(subtitle_user) eq 0) then subtitle_user=''
62	if (n_elements(pal) eq 0) then pal=33
63	if (n_elements(format) eq 0) then format='(F6.2)'
64	if (format eq '') then format='(F6.2)'
65	;if (n_elements(xtitleset) eq 0) then xtitleset='longitude'
66	;if (n_elements(ytitleset) eq 0) then ytitleset='latitude'
67	if (n_elements(title_axis) ne 2) then begin
68	xtitleset='longitude'
69	ytitleset='latitude'
70	endif else begin
71	xtitleset=title_axis(0)
72	ytitleset=title_axis(1)
73	endelse
74	if (n_elements(isotropic) eq 0) then isotropic='true'
75	if (n_elements(coord2d) eq 0) then coord2d='false'
76
77	;------------------
78	; limits
79	;------------------
80	;if (minfield_init*maxfield_init eq 0) then begin
81	if (minfield_init eq maxfield_init) then begin
82	;---different min/max for each layer
83	w=where(abs(what_I_plot) lt missing_value)
84	if (w(0) ne -1) then begin
85	minfield=min(what_I_plot[w])
86	maxfield=max(what_I_plot[w])
87	;goodscale, $
88	; what_I_plot, $
89	; minfield, $
90	; maxfield
91	endif else begin
92	what_I_plot=0.
93	print, 'no valid values in this field'
94	endelse
95	endif else begin
96	;---user-defined limits
97	minfield=minfield_init
98	maxfield=maxfield_init
99	;;w=where(what_I_plot lt minfield) & if (w(0) ne -1) then what_I_plot[w]=minfield
100	;;w=where(what_I_plot gt maxfield) & if (w(0) ne -1) then what_I_plot[w]=maxfield
101	endelse
102
103	if ((minfield+maxfield eq 0) and (abs(minfield) ne abs(maxfield))) then begin
104	print, 'nothing to plot ... skipping this plot ...'
105	endif else begin
106
107
108	;-----------------------------------------
109	; fix the possible longitude ugliness :)
110	; ...get rid of the -180/180 limit
111	;-----------------------------------------
112
113	;;;;; PB WITH 2D COORDINATES
114	;w=where(lon eq min(lon))
115	;if (w(0) ne 0) then begin
116	; lon[0:w(0)-1]=lon[0:w(0)-1]-360
117	;endif
118
119
120	;---------------
121	; adapt ticks
122	;---------------
123	if (n_elements(windowx) eq 0) then windowx=[min(lon),max(lon)]
124	if (n_elements(windowy) eq 0) then windowy=[min(lat),max(lat)]
125	if (n_elements(intervalx) eq 0) then intervalx=round((windowx(1)-windowx(0))/8.)
126	if (n_elements(intervaly) eq 0) then intervaly=round((windowy(1)-windowy(0))/5.)
127
128
129	;------------------
130	; plot window
131	;------------------
132
133	;lons=lon & lats=lat
134	;;--------------------------------------------
135	;;
136	;; ANOMALY
137	;;
138	;w=where( (lon ge windowx(0)) and (lon le windowx(1)) ) & w2=where( (lat ge windowy(0)) and (lat le windowy(1)) )
139	;yeah=fltarr(n_elements(w),n_elements(w2)) & for i=0,n_elements(w)-1 do for j=0,n_elements(w2)-1 do yeah(i,j)=what_I_plot(w(i),w2(j))
140	;what_I_plot = yeah - mean(yeah)
141	;lon=lon[w] & lat=lat[w2]
142	;;--------------------------------------------
143
144
145	; to ensure the right limits
146	dumb_what_I_plot=what_I_plot
147	dumb_what_I_plot(,)=minfield & print, minfield
148	dumb_what_I_plot(0,0)=maxfield & print, maxfield
149
150	;;
151	;;
152	;subtitle_hunhun = subtitle_user
153	;subtitle_user = ''
154	;;
155	;;
156
157
158	pospos = [0.10, 0.12, 0.90, 0.92] ;; iso
159	;pospos = [0.15, 0.35, 0.95, 0.95] ;; non-iso
160
161	if (isotropic eq 'true') then begin
162	loadct, 0, /silent
163	if ( coord2d ne 'polar' ) then begin
164	;;trace non polaire
165	contour, /nodata, $
166	/isotropic, $
167	position=pospos, $
168	dumb_what_I_plot, lon, lat, $
169	title=title_user, xtitle=xtitleset, ytitle=ytitleset, $
170	xrange=windowx, yrange=windowy, $
171	xtickinterval=intervalx, ytickinterval=intervaly, $
172	subtitle=subtitle_user, color=0
173	endif else begin
174	;;trace avec MAP_SET -- MAP_SET -- MAP_SET
175	print,'------------------- Polar plot'
176	contour, /nodata, $
177	/isotropic, $
178	/overplot, $
179	dumb_what_I_plot, lon, lat, $
180	title=title_user, xtitle=xtitleset, ytitle=ytitleset, $
181	xrange=windowx, yrange=windowy, $
182	xtickinterval=intervalx, ytickinterval=intervaly, $
183	subtitle=subtitle_user, color=0
184	endelse
185	subtitle_user = ''
186	endif else begin
187	pospos = [0.15, 0.35, 0.95, 0.95] ;; non-iso
188	loadct, 0, /silent
189	if ( coord2d ne 'polar' ) then begin
190	;;trace non polaire
191	contour, /nodata, $
192	position=pospos, $
193	dumb_what_I_plot, lon, lat, $
194	title=title_user, xtitle=xtitleset, ytitle=ytitleset, $
195	xrange=windowx, yrange=windowy, $
196	xtickinterval=intervalx, ytickinterval=intervaly, $
197	color=0,subtitle=subtitle_user
198	endif else begin
199	;;trace avec MAP_SET -- MAP_SET -- MAP_SET
200	print,'------------------- Polar plot'
201	contour, /nodata, $
202	/overplot, $
203	dumb_what_I_plot, lon, lat, $
204	title=title_user, xtitle=xtitleset, ytitle=ytitleset, $
205	xrange=windowx, yrange=windowy, $
206	xtickinterval=intervalx, ytickinterval=intervaly, $
207	color=0,subtitle=subtitle_user
208	endelse
209	endelse
210
211	;------------------
212	; plot field
213	;------------------
214
215	; to avoid spurious blank zones
216	nlevels=colors
217	levp=minfield+(findgen(nlevels)/float(nlevels-1))*(maxfield-minfield)
218
219	loadct,pal,/silent
220	;reverse_ct
221	;if (pal eq 0) then reverse_ct
222	contour, what_I_plot, $
223	lon,lat, $
224	; nlevels=colors, $
225	levels=levp, $
226	/cell_fill, $
227	max_value=maxfield, $
228	min_value=minfield, $
229	/overplot
230
231
232	;nlat=128.
233	;oplot, lat0.+nlat100./1000., linestyle=1
234
235	;;; ANOMALY
236	;lon=lons
237	;lat=lats
238
239	;------------------
240	; colorbar
241	;------------------
242
243	;if (flag_cb eq 'true') then begin
244	; format_colorbar=format
245	; colorbar, $
246	; position=[0.15, 0.85, 0.95, 0.90], $
247	; divisions=8, $
248	; range=[float(minfield),float(maxfield)], $
249	; format=format_colorbar
250	;endif
251	;goto, skip_new
252
253	if (flag_cb eq 'true') then begin
254	format_colorbar=format
255
256	;;; ndiv peut etre regle
257	if (n_elements(ndiv) eq 0) then ndiv=10
258
259	;if (n_elements(poscb) eq 0) then poscb=0.95
260	;if (n_elements(poscb) eq 0) then poscb=0.80
261	if (n_elements(poscb) eq 0) then poscb=0.70
262
263	if (isotropic eq 'true') then begin
264	;poscb=0.85 &
265	pos=[poscb, 0.12, poscb+0.03, 0.92]
266	;poscb=0.95 & pos=[poscb, 0.12, poscb+0.03, 0.92]
267	; poscb=0.96 & pos=[poscb, 0.12, poscb+0.03, 0.92]
268	; poscb=0.05 & pos=[poscb, 0.12, poscb+0.03, 0.92]
269	;colorbar = Obj_New("COLORBAR",$
270	; Range=[float(minfield),float(maxfield)],$
271	; /vertical, $
272	; Ncolors=255,$
273	; charsize=0.85,$
274	; format=format_colorbar,$
275	; major=ndiv,$
276	; ticklen=-0.15,$
277	; position=pos)
278	;colorbar->Draw
279	;Obj_Destroy, colorbar
280	colorbar, $
281	/vertical, $
282	position=pos, $
283	range=[float(minfield),float(maxfield)], $
284	divisions=ndiv, $
285	ncolors=255,$
286	; charsize=0.85,$
287	ticklen=-0.15,$
288	format=format_colorbar
289
290	xyouts, pos(2), 0.01, subtitle_user, charsize=0.8, alignment=1., /NORMAL
291
292	if (n_elements(nam2) ne 0) then begin
293	thedate=STRSPLIT(nam2, /EXTRACT)
294	thedate=STRSPLIT(thedate(2),"'",/EXTRACT)
295	thedate=STRSPLIT(thedate(0),"_",/EXTRACT)
296	thedate=STRSPLIT(thedate(0),"-",/EXTRACT)
297	xyouts, pos(2), 0.01, 'LMD Martian Mesoscale Model - d'+thedate(2)+'/m'+thedate(1), charsize=0.7, alignment=1., /NORMAL
298	endif else begin
299	thedate=''
300	endelse
301	;;
302	;;
303	;;
304	;xyouts, pos(2), 0.01, 'LMD Mars Mesoscale Model / '+subtitle_hunhun, charsize=0.95, alignment=1., /NORMAL
305	;;
306	;;
307	;;
308
309	endif else begin
310
311	pos=[0.15, 0.17, 0.95, 0.20]
312	if (subtitle_user eq '') then pos=[0.15, 0.20, 0.95, 0.23]
313
314	;colorbar = Obj_New("COLORBAR",$
315	; Range=[float(minfield),float(maxfield)],$
316	; Ncolors=255,$
317	; charsize=0.85,$
318	; format=format_colorbar,$
319	; major=ndiv,$
320	; ticklen=-0.15,$
321	; position=pos)
322	;colorbar->Draw
323	;Obj_Destroy, colorbar
324	colorbar, $
325	position=pos, $
326	range=[float(minfield),float(maxfield)], $
327	divisions=ndiv, $
328	ncolors=255,$
329	; charsize=0.85,$
330	ticklen=-0.15,$
331	format=format_colorbar
332
333	if (n_elements(nam2) ne 0) then begin
334	thedate=STRSPLIT(nam2, /EXTRACT)
335	thedate=STRSPLIT(thedate(2),"'",/EXTRACT)
336	thedate=STRSPLIT(thedate(0),"_",/EXTRACT)
337	thedate=STRSPLIT(thedate(0),"-",/EXTRACT)
338	xyouts, 0.95, pos(1)-0.08, 'LMD Martian Mesoscale Model - d'+thedate(2)+'/m'+thedate(1), charsize=0.7, alignment=1., /NORMAL
339	endif else begin
340	thedate=''
341	endelse
342
343	endelse
344	endif
345
346	skip_new:
347
348	;;;;;;;;;;;;;;;;;;;;;;;;
349	;;;;;;;;;;;;;;;;;;;;;;;;
350	;;; ajouter un point
351	;loadct, 0
352	;xyouts, 353.87-360, -1.88, '+', ALIGNMENT=0.5, CHARSIZE=2, CHARTHICK=2
353	;;;;;;;;;;;;;;;;;;;;;;;;
354	;;;;;;;;;;;;;;;;;;;;;;;;
355
356	;goto, nocontour
357	;print, 'contour'
358	;--------------------
359	; overplot contour
360	;--------------------
361	if (n_elements(overcontour) ne 1) then begin
362
363	w=where(abs(overcontour) lt missing_value)
364	if (w(0) ne -1) then begin
365	min_contour=min(overcontour[w])
366	max_contour=max(overcontour[w])
367	endif
368
369
370	; to avoid spurious blank zones
371	nlevels=colors/2
372	if (n_elements(lev) eq 0) then lev=min_contour+(findgen(nlevels)/float(nlevels-1))*(max_contour-min_contour)
373
374	;lev=round(lev)
375	;yeah=UNIQ(lev) & lev=lev[yeah]
376
377
378	if (isotropic eq 'true') then begin
379	loadct,0,/silent
380	;reverse_ct
381	contour, overcontour, $
382	lon,lat, $
383	/overplot, $
384	; nlevels=colors/2, $
385	levels=lev, $
386	C_LINESTYLE = (lev LT 0.0), $ ;; dotted lines if negative
387	max_value=max_contour, $
388	min_value=min_contour, $
389	c_labels=findgen(n_elements(lev))*0., $
390	; c_charsize=!P.charsize/2, $
391	color=75;0;50;75;100tropclair;150;0;, $
392	; /noerase, $
393	; xtitle=xtitleset, $
394	; xrange=windowx, $
395	; ytitle=ytitleset, $
396	; yrange=windowy, $
397	; position=pospos, $
398	; /isotropic, $
399	; xtickinterval=intervalx, $
400	; ytickinterval=intervaly
401	endif else begin
402	loadct,0,/silent
403	contour, overcontour, $
404	lon,lat, $
405	/overplot, $
406	; nlevels=colors/2, $
407	levels=lev, $
408	C_LINESTYLE = (lev LT 0.0), $ ;; dotted lines if negative
409	max_value=max_contour, $
410	min_value=min_contour, $
411	c_labels=findgen(n_elements(lev))*0., $
412	color=0;, $
413	; /noerase, $
414	; xtitle=xtitleset, $
415	; xrange=windowx, $
416	; ytitle=ytitleset, $
417	; yrange=windowy, $
418	;; position=[0.15, 0.15, 0.95, 0.75], $
419	; xtickinterval=intervalx, $
420	; ytickinterval=intervaly
421	endelse
422
423
424	endif
425	nocontour:
426
427	;--------------------
428	; overplot wind
429	;--------------------
430	if (n_elements(overvector_x)*n_elements(overvector_y) ne 1) then begin
431
432	w=where(abs(overvector_x) ge missing_value) ; u and v have the same missing values ...
433	if (w(0) ne -1) then begin
434	overvector_x[w]=0.
435	overvector_y[w]=0.
436	endif
437
438
439	len=2.0
440	if (n_elements(windex) eq 0) then windex=20.
441	ref_mag=round(windex)
442	mytext=STRTRIM(STRING(ref_mag),2)+'ms!E-1!N'
443
444	if (n_elements(stride) eq 0) then begin
445	;stride=3 ; STRIDE: Pick every nth vector for plotting. (1)
446	;stride=2
447	stride=5
448	endif else begin
449	stride=round(stride)
450	endelse
451
452	ref_pos=[0.10, 0.90]
453	ref_pos=[0.15, 0.10]
454	ref_pos=[0.20, 0.05]
455	;ref_pos=[0.20, 0.02]
456
457	s=size(lon)
458	if (s[0] eq 1) then begin
459	zeu=overvector_x
460	zev=overvector_y
461	longs=lon
462	lats=lat
463	endif else begin
464	minlat=min(lat) & maxlat=max(lat) & minlon=min(lon) & maxlon=max(lon)
465	npoints=2n_elements(lon(,0)) ;; trop de points, mais au moins on ne perd rien (2 et 3 donnent results similaires)
466	TRIANGULATE, lon, lat, tr
467	zeu = GRIDDATA( lon, lat, overvector_x, /LINEAR, triangles=tr, dimension=npoints, MISSING=!VALUES.F_NAN )
468	zev = GRIDDATA( lon, lat, overvector_y, /LINEAR, triangles=tr, dimension=npoints, MISSING=!VALUES.F_NAN )
469	longs = minlon + (maxlon - minlon)*findgen(npoints)/float(npoints-1)
470	lats = minlat + (maxlat - minlat)*findgen(npoints)/float(npoints-1)
471	endelse
472
473	wlon=where( (longs gt windowx(0)) and (longs lt windowx(1)) )
474	wlat=where( (lats gt windowy(0)) and (lats lt windowy(1)) )
475
476	if ((wlon(0) eq -1) or (wlat(0) eq -1)) then begin
477	print, 'wrong window settings'
478	print, windowx(0), windowx(1)
479	print, longs
480	print, windowy(0), windowy(1)
481	print, lats
482	stop
483	endif
484
485	myu=dblarr(n_elements(longs),n_elements(lats))
486	myv=dblarr(n_elements(longs),n_elements(lats))
487
488	myu(,)=!Values.F_NAN
489	myv(,)=!Values.F_NAN
490
491	for i=0,n_elements(wlon)-1 do begin
492	for j=0,n_elements(wlat)-1 do begin
493	myu(wlon(i),wlat(j))=zeu(wlon(i),wlat(j))
494	myv(wlon(i),wlat(j))=zev(wlon(i),wlat(j))
495	endfor
496	endfor
497
498	;; eliminate vector for colorbar clarity ?
499	;nxs=1
500	;nxe=n_elements(longs(*))-2
501	;nys=1
502	;nye=n_elements(lats(*))-2
503	;myu=overvector_x(nxs:nxe,nys:nye)
504	;myv=overvector_y(nxs:nxe,nys:nye)
505	;longs=lon(nxs:nxe)
506	;lats=lat(nys:nye)
507
508	loadct,0,/silent
509	VECTOR, myu, myv, longs, lats, LENGTH=len, REF_MAG=ref_mag,$
510	STRIDE=stride, TYPE=2, HEAD_LEN=0.3, ANGLE=30, $
511	REF_POS=ref_pos, REF_TEXT=mytext,$
512	ALIGN=0.5, COLOR=1, /OVERPLOT;, /NOERASE
513
514	yeah:
515	endif
516
517
518	endelse
519
520
521
522	goto, why_contour_end
523	;;;;;;;;;;;;;;;;;;;;;;
524	;;;;;;;;;;;;;;;;;;;;;;
525	;;;;;;;;;;;;;;;;;;;;;;
526
527	if (n_elements(overcontour) ne 1) then begin
528
529	w=where(abs(overcontour) lt missing_value)
530	if (w(0) ne -1) then begin
531	min_contour=min(overcontour[w])
532	max_contour=max(overcontour[w])
533	endif
534
535
536	; to avoid spurious blank zones
537	nlevels=colors/2
538	if (n_elements(lev) eq 0) then lev=min_contour+(findgen(nlevels)/float(nlevels-1))*(max_contour-min_contour)
539	;lev=round(lev)
540	;yeah=UNIQ(lev) & lev=lev[yeah]
541
542
543	if (isotropic eq 'true') then begin
544	loadct,0,/silent
545	contour, overcontour, $
546	lon,lat, $
547	; nlevels=colors/2, $
548	levels=lev, $
549	C_LINESTYLE = (lev LT 0.0), $ ;; dotted lines if negative
550	max_value=max_contour, $
551	min_value=min_contour, $
552	c_labels=findgen(nlevels), $ ;*0., $
553	color=0, $
554	/noerase, $
555	xtitle=xtitleset, $
556	xrange=windowx, $
557	ytitle=ytitleset, $
558	yrange=windowy, $
559	; position=[0.15, 0.15, 0.95, 0.75], $
560	/isotropic, $
561	xtickinterval=intervalx, $
562	ytickinterval=intervaly
563	endif else begin
564	loadct,0,/silent
565	contour, overcontour, $
566	lon,lat, $
567	; nlevels=colors/2, $
568	levels=lev, $
569	C_LINESTYLE = (lev LT 0.0), $ ;; dotted lines if negative
570	max_value=max_contour, $
571	min_value=min_contour, $
572	c_labels=findgen(nlevels), $
573	color=0, $
574	/noerase, $
575	xtitle=xtitleset, $
576	xrange=windowx, $
577	ytitle=ytitleset, $
578	yrange=windowy, $
579	position=[0.15, 0.15, 0.95, 0.75], $
580	xtickinterval=intervalx, $
581	ytickinterval=intervaly
582	endelse
583
584
585	endif
586
587	why_contour_end:
588
589	end

Note: See TracBrowser for help on using the repository browser.

Download in other formats:

Original Format