Context Navigation

read_ncdf_varget.pro @ 212

Last change on this file since 212 was 212, checked in by smasson, 17 years ago
include [yz]reverse in read_ncdf_varget
Property svn:eol-style set to `native` Property svn:keywords set to `Id`
File size: 23.6 KB

Line
1	;+
2	;
3	; @todo seb
4	;
5	;-
6	; le format netcdf et IDL nous permet d''aller lire uniquement la partie de
7	; tableau qui nous interesse.
8	;
9	; attention, qd il y a un shift il faut faire attention a ce que
10	; l''on fait et en particulier arriver a positionner la boite
11	; specifiee par firstx:lastx par rapport au tableau contenu dans
12	; le fichier NetCdf.
13	; pour s''y reperer voici un petit dessin
14	;
15	; dans le repere 0:jpi, avec le zoom definit par firstx:lastx
16	;
17	; 0 j key-1 key i jpi-1
18	; \|......\|........\| \|-----------\|----------------\|
19	;
20	; si on veut savoir comment c''etait avant le shift:
21	;
22	; 0 i-key jpi-key-1 jpi-key jpi-key+j jpi-1
23	; \|--------\|--------------\| \|........\|..........\|
24	;
25	; avant que l''on ajuste au zoom definit par ixminmesh,ixmaxmesh c''etait:
26	;
27	; 0 ixmin i-key jpi-key-1 jpi-key jpi-key+j jpi-1 jpidta
28	; +ixmin +ixmin +ixmin +ixmin +ixmin
29	; \|,,,,\|--------\|--------------\| \|........\|..........\|,,,,,,,,,,,\|
30	;
31	; apres il suffit de remplacer i et j par firstx ou lastx qd on
32	; ne fait pas de stride.
33	;
34	; When using stride, things are getting more complicated....
35	;
36	; 1) when must replace firstx by firstx*key_stride[0] and
37	; lastx by lastx*key_stride[0]
38	; 2) we must use jpitotal instead of jpi
39	; 3) finding the position of the offset in the right side part is
40	; tricky. look at this illustration of stride = 4 ...
41	; 0 j key-1 key i jpi-1
42	; +...+..\|+...+...+ \|--+---+---+\|--+---+---+---+-\|
43	;
44	; 0 i-key jpi-key-1 jpi-key jpi-key+j jpi-1
45	; \|--+---+---+\|--+---+---+---+-\| +...+..\|+...+...+
46	; in the ........ part, it is easy the first point is number 0 but in
47	; the --------- part the first point is number 3... Ihe position of
48	; the first point in the is given by:
49	; (key_stride[0]-1)-((key-1) MOD key_stride[0])
50	; 4) last point...by default, the number of element read by IDL when
51	; using stride is given by n_elements/stride. However, we must read:
52	; (n_elements+stride-1)/stride or (n_elements-1)/key_stride+1
53	; for example: for n_elements between 10 and 12, with a stride of 3,
54	; we must read 4 points instead of 3...: +--+--+--+--
55	; This problem as an easy solution by using the keyword count with
56	; the appropiate value!
57	;
58	ixmin = ixminmesh-ixmindta
59	iymin = iyminmesh-iymindta
60	izmin = izminmesh-izmindta
61	jpitotal = long(ixmaxmesh-ixminmesh+1)
62	key = long(key_shift MOD jpitotal)
63	if key LT 0 then key = key+jpitotal
64	;
65	case 1 OF
66	;......................................................................
67	;......................................................................
68	varcontient.ndims eq 2:BEGIN ;xy array
69	;......................................................................
70	;......................................................................
71	case 1 OF
72	;,,,,,,,,,,,,,,,,,,,,,,,,
73	keyword_set(key_shift) EQ 0 AND total(key_stride) EQ 3:BEGIN
74	; on ne shift pas et il n''y a pas de stride
75	;,,,,,,,,,,,,,,,,,,,,,,,,
76	ncdf_varget,cdfid,name,res,offset=[firstx+ixmin,firsty+iymin] $
77	,count=[nx,ny]
78	END
79	;,,,,,,,,,,,,,,,,,,,,,,,,
80	keyword_set(key_shift) EQ 0 AND total(key_stride) NE 3:BEGIN
81	; on ne shift pas mais il y a un stride
82	;,,,,,,,,,,,,,,,,,,,,,,,,
83	ncdf_varget,cdfid,name,res,offset=[firstx*key_stride[0]+ixmin $
84	,firsty*key_stride[1]+iymin] $
85	,count=[nx,ny], stride = key_stride[0:1]
86	END
87	;,,,,,,,,,,,,,,,,,,,,,,,,
88	keyword_set(key_shift) NE 0 AND total(key_stride) EQ 3:BEGIN
89	; on shift mais il n''y a pas de stride
90	;,,,,,,,,,,,,,,,,,,,,,,,,
91	case 1 of
92	; --------- part, we can directly extract the array in one piece
93	firstx GE key:BEGIN
94	ncdf_varget,cdfid,name,res,offset=[ixmin-key+firstx $
95	,firsty+iymin] $
96	,count=[nx,ny]
97	END
98	; ......... part, we can directly extract the array in one piece
99	lastx LE key-1:BEGIN
100	ncdf_varget,cdfid,name,res,offset=[jpi-key+ixmin+firstx $
101	,firsty+iymin] $
102	,count=[nx,ny]
103	END
104	ELSE:BEGIN ; we have to extract the array in 2 pieces...
105	; ......... part, first part...
106	ncdf_varget,cdfid,name,tab1,offset=[jpi-key+ixmin+firstx $
107	,firsty+iymin] $
108	,count=[key-1-firstx+1,ny]
109	; --------- part, second part...
110	ncdf_varget,cdfid,name,tab2,offset=[ixmin,firsty+iymin] $
111	,count=[lastx-key+1,ny]
112	res = [temporary(tab1), temporary(tab2)]
113	END
114	ENDCASE
115	END
116	;,,,,,,,,,,,,,,,,,,,,,,,,
117	keyword_set(key_shift) NE 0 AND total(key_stride) NE 3:BEGIN
118	; on shift et il y a un stride
119	;,,,,,,,,,,,,,,,,,,,,,,,,
120	case 1 OF ; case sur la facon de fire le champ
121	; --------- part, we can directly extract the array in one piece
122	firstx*key_stride[0] GE key:BEGIN
123	ncdf_varget,cdfid,name,res,offset=[ixmin $
124	+(key_stride[0]-1)-((key-1) MOD key_stride[0]) $
125	+(firstx-((key-1)/key_stride[0]+1))*key_stride[0] $
126	,firsty*key_stride[1]+iymin] $
127	,count=[nx,ny], stride = key_stride[0:1]
128	END
129	; ......... part, we can directly extract the array in one piece
130	lastx*key_stride[0] LE key-1:BEGIN
131	ncdf_varget,cdfid,name,res,offset=[jpitotal-key+ixmin+firstx*key_stride[0] $
132	,firsty*key_stride[1]+iymin] $
133	,count=[nx,ny], stride = key_stride[0:1]
134	END
135	ELSE:BEGIN ; we have to extract the array in 2 pieces...
136	; ......... part, first part...
137	ncdf_varget,cdfid,name,tab1,offset=[jpitotal-key+ixmin+firstx*key_stride[0] $
138	,firsty*key_stride[1]+iymin] $
139	,count=[(key-firstx*key_stride[0]-1)/key_stride[0]+1,ny] $
140	, stride = key_stride[0:1]
141	; --------- part, second part...
142	ncdf_varget,cdfid,name,tab2,offset=[ixmin $
143	+(key_stride[0]-1)-((key-1) MOD key_stride[0]) $
144	, firsty*key_stride[1]+iymin] $
145	,count=[nx-(key-firstx*key_stride[0]-1)/key_stride[0]-1,ny] $
146	, stride = key_stride[0:1]
147	; in one unique array...
148	res = [temporary(tab1), temporary(tab2)]
149	END
150	ENDCASE ; case sur la facon de fire le champ
151	END
152	ENDCASE ; differentes possibilites de key_shift et key_performance
153	END
154	;......................................................................
155	;......................................................................
156	varcontient.ndims eq 3 AND (where(varcontient.dim EQ contient.recdim))[0] EQ -1:BEGIN ;xyz array
157	;......................................................................
158	;......................................................................
159	case 1 of
160	;,,,,,,,,,,,,,,,,,,,,,,,,
161	keyword_set(key_shift) EQ 0 AND total(key_stride) EQ 3:BEGIN
162	; on ne shift pas et il n''y a pas de stride
163	;,,,,,,,,,,,,,,,,,,,,,,,,
164	ncdf_varget,cdfid,name,res,offset=[firstx+ixmin,firsty+iymin,firstz+izmin] $
165	,count=[nx,ny,nz]
166	END
167	;,,,,,,,,,,,,,,,,,,,,,,,,
168	keyword_set(key_shift) EQ 0 AND total(key_stride) NE 3:BEGIN
169	; on ne shift pas mais il y a un stride
170	;,,,,,,,,,,,,,,,,,,,,,,,,
171	ncdf_varget,cdfid,name,res,offset=[firstx*key_stride[0]+ixmin $
172	,firsty*key_stride[1]+iymin $
173	,firstz*key_stride[2]+izmin] $
174	,count=[nx,ny,nz], stride = key_stride
175	END
176	;,,,,,,,,,,,,,,,,,,,,,,,,
177	keyword_set(key_shift) NE 0 AND total(key_stride) EQ 3:BEGIN
178	; on shift mais il n''y a pas de stride
179	;,,,,,,,,,,,,,,,,,,,,,,,,
180	case 1 of
181	firstx GE key:BEGIN ; on peut tout couper d''un coup
182	ncdf_varget,cdfid,name,res,offset=[ixmin-key+firstx $
183	,firsty+iymin $
184	,firstz+izmin] $
185	,count=[nx,ny,nz]
186	END
187	lastx LE key-1:BEGIN ; on peut tout couper d''un coup
188	ncdf_varget,cdfid,name,res,offset=[jpi-key+ixmin+firstx $
189	,firsty+iymin $
190	,firstz+izmin] $
191	,count=[nx,ny,nz]
192	END
193	ELSE:BEGIN ; Le tableau est separe en 2 morceaux...
194	ncdf_varget,cdfid,name,tab1,offset=[jpi-key+ixmin+firstx $
195	,firsty+iymin $
196	,firstz+izmin] $
197	,count=[key-1-firstx+1,ny,nz]
198	ncdf_varget,cdfid,name,tab2,offset=[ixmin $
199	,firsty+iymin $
200	,firstz+izmin] $
201	,count=[lastx-key+1,ny,nz]
202	res = [temporary(tab1), temporary(tab2)]
203	END
204	ENDCASE
205	END
206	;,,,,,,,,,,,,,,,,,,,,,,,,
207	keyword_set(key_shift) NE 0 AND total(key_stride) NE 3:BEGIN
208	; on shift et il y a un stride
209	;,,,,,,,,,,,,,,,,,,,,,,,,
210	case 1 OF ; case sur la facon de fire le champ
211	firstx*key_stride[0] GE key:BEGIN ; on peut tout couper d''un coup
212	ncdf_varget,cdfid,name,res,offset=[ixmin $
213	+(key_stride[0]-1)-((key-1) MOD key_stride[0]) $
214	+(firstx-((key-1)/key_stride[0]+1))*key_stride[0] $
215	,firsty*key_stride[1]+iymin $
216	,firstz*key_stride[2]+izmin] $
217	,count=[nx,ny,nz], stride = key_stride
218	END
219	lastx*key_stride[0] LE key-1:BEGIN ; on peut tout couper d''un coup
220	ncdf_varget,cdfid,name,res,offset=[jpitotal-key+ixmin+firstx*key_stride[0] $
221	,firsty*key_stride[1]+iymin $
222	,firstz*key_stride[2]+izmin] $
223	,count=[nx,ny,nz], stride = key_stride
224	END
225	ELSE:BEGIN ; le tableau est separe en 2 morceaux...
226	ncdf_varget,cdfid,name,tab1,offset=[jpitotal-key+ixmin+firstx*key_stride[0] $
227	,firsty*key_stride[1]+iymin $
228	,firstz*key_stride[2]+izmin] $
229	,count=[(key-firstx*key_stride[0]-1)/key_stride[0]+1,ny,nz] $
230	, stride = key_stride
231	; 2eme bout...
232	ncdf_varget,cdfid,name,tab2,offset=[ixmin $
233	+(key_stride[0]-1)-((key-1) MOD key_stride[0]) $
234	,firsty*key_stride[1]+iymin $
235	,firstz*key_stride[2]+izmin] $
236	,count=[nx-(key-firstx*key_stride[0]-1)/key_stride[0]-1,ny,nz] $
237	, stride = key_stride
238	; on recolle le tout
239	res = [temporary(tab1), temporary(tab2)]
240	END
241	ENDCASE ; case sur la facon de fire le champ
242	END
243	ENDCASE ; differentes possibilites de key_shift et key_performance
244	END
245	;......................................................................
246	;......................................................................
247	varcontient.ndims eq 3 AND (where(varcontient.dim EQ contient.recdim))[0] NE -1:BEGIN ;xyt array
248	;......................................................................
249	;......................................................................
250	case 1 of
251	;,,,,,,,,,,,,,,,,,,,,,,,,
252	keyword_set(key_shift) EQ 0 AND total(key_stride) EQ 3:BEGIN
253	; on ne shift pas et il n''y a pas de stride
254	;,,,,,,,,,,,,,,,,,,,,,,,,
255	ncdf_varget,cdfid,name,res,offset=[firstx+ixmin $
256	,firsty+iymin $
257	,firsttps] $
258	,count=[nx,ny,jpt]
259	END
260	;,,,,,,,,,,,,,,,,,,,,,,,,
261	keyword_set(key_shift) EQ 0 AND total(key_stride) NE 3:BEGIN
262	; on ne shift pas mais il y a un stride
263	;,,,,,,,,,,,,,,,,,,,,,,,,
264	ncdf_varget,cdfid,name,res,offset=[firstx*key_stride[0]+ixmin $
265	,firsty*key_stride[1]+iymin $
266	,firsttps] $
267	,count=[nx,ny,jpt], stride = [key_stride[0:1], 1]
268	END
269	;,,,,,,,,,,,,,,,,,,,,,,,,
270	keyword_set(key_shift) NE 0 AND total(key_stride) EQ 3:BEGIN
271	; on shift mais il n''y a pas de stride
272	;,,,,,,,,,,,,,,,,,,,,,,,,
273	case 1 of
274	; --------- part, we can directly extract the array in one piece
275	firstx GE key:BEGIN ; on peut tout couper d''un coup
276	ncdf_varget,cdfid,name,res,offset=[ixmin-key+firstx $
277	,firsty+iymin $
278	,firsttps] $
279	,count=[nx,ny,jpt]
280	END
281	; ......... part, we can directly extract the array in one piece
282	lastx LE key-1:BEGIN ; on peut tout couper d''un coup
283	ncdf_varget,cdfid,name,res,offset=[jpi-key+ixmin+firstx $
284	,firsty+iymin $
285	,firsttps] $
286	,count=[nx,ny,jpt]
287	END
288	ELSE:BEGIN ; Le tableau est separe en 2 morceaux...
289	; ......... part, first part...
290	ncdf_varget,cdfid,name,tab1,offset=[jpi-key+ixmin+firstx $
291	,firsty+iymin $
292	,firsttps] $
293	,count=[key-1-firstx+1,ny,jpt]
294	; --------- part, second part...
295	ncdf_varget,cdfid,name,tab2,offset=[ixmin $
296	,firsty+iymin $
297	,firsttps] $
298	,count=[lastx-key+1,ny,jpt]
299	res = [temporary(tab1), temporary(tab2)]
300	END
301	ENDCASE
302	END
303	;,,,,,,,,,,,,,,,,,,,,,,,,
304	keyword_set(key_shift) NE 0 AND total(key_stride) NE 3:BEGIN
305	; on shift et il y a un stride
306	;,,,,,,,,,,,,,,,,,,,,,,,,
307	; --------- part, we can directly extract the array in one piece
308	case 1 OF ; case sur la facon de fire le champ
309	firstx*key_stride[0] GE key:BEGIN
310	ncdf_varget,cdfid,name,res,offset=[ixmin $
311	+(key_stride[0]-1)-((key-1) MOD key_stride[0]) $
312	+( firstx-((key-1)/key_stride[0]+1) )*key_stride[0] $
313	,firsty*key_stride[1]+iymin $
314	,firsttps] $
315	,count=[nx,ny,jpt], stride = [key_stride[0:1], 1]
316	END
317	; ......... part, we can directly extract the array in one piece
318	lastx*key_stride[0] LE key-1:BEGIN
319	ncdf_varget,cdfid,name,res,offset=[jpitotal-key+ixmin+firstx*key_stride[0] $
320	,firsty*key_stride[1]+iymin $
321	,firsttps] $
322	,count=[nx,ny,jpt], stride = [key_stride[0:1], 1]
323	END
324	ELSE:BEGIN
325	; ......... part, first part...
326	ncdf_varget,cdfid,name,tab1,offset=[jpitotal-key+ixmin+firstx*key_stride[0] $
327	,firsty*key_stride[1]+iymin $
328	,firsttps] $
329	,count=[(key-firstx*key_stride[0]-1)/key_stride[0]+1,ny,jpt] $
330	, stride = [key_stride[0:1], 1]
331	; --------- part, second part...
332	ncdf_varget,cdfid,name,tab2,offset=[ixmin $
333	+(key_stride[0]-1)-((key-1) MOD key_stride[0]) $
334	, firsty*key_stride[1]+iymin $
335	,firsttps] $
336	,count=[nx-(key-firstx*key_stride[0]-1)/key_stride[0]-1,ny,jpt] $
337	, stride = [key_stride[0:1], 1]
338	; on recolle le tout
339	res = [temporary(tab1), temporary(tab2)]
340	END
341	ENDCASE ; case sur la facon de fire le champ
342	END
343	endcase
344	END
345	;......................................................................
346	;......................................................................
347	varcontient.ndims eq 4:BEGIN ;xyzt array
348	;......................................................................
349	;......................................................................
350	case 1 of
351	;,,,,,,,,,,,,,,,,,,,,,,,,
352	keyword_set(key_shift) EQ 0 AND total(key_stride) EQ 3:BEGIN
353	; on ne shift pas et il n''y a pas de stride
354	;,,,,,,,,,,,,,,,,,,,,,,,,
355	ncdf_varget,cdfid,name,res,offset=[firstx+ixmin $
356	,firsty+iymin $
357	,firstz+izmin $
358	,firsttps] $
359	,count=[nx,ny,nz,jpt]
360	END
361	;,,,,,,,,,,,,,,,,,,,,,,,,
362	keyword_set(key_shift) EQ 0 AND total(key_stride) NE 3:BEGIN
363	; on ne shift pas mais il y a un stride
364	;,,,,,,,,,,,,,,,,,,,,,,,,
365	ncdf_varget,cdfid,name,res,offset=[firstx*key_stride[0]+ixmin $
366	,firsty*key_stride[1]+iymin $
367	,firstz*key_stride[2]+izmin $
368	,firsttps] $
369	,count=[nx,ny,nz,jpt], stride = [key_stride, 1]
370	END
371	;,,,,,,,,,,,,,,,,,,,,,,,,
372	keyword_set(key_shift) NE 0 AND total(key_stride) EQ 3:BEGIN
373	; on shift mais il n''y a pas de stride
374	;,,,,,,,,,,,,,,,,,,,,,,,,
375	case 1 of
376	firstx GE key:BEGIN ; on peut tout couper d''un coup
377	ncdf_varget,cdfid,name,res,offset=[ixmin-key+firstx $
378	,firsty+iymin $
379	,firstz+izmin $
380	,firsttps] $
381	,count=[nx,ny,nz,jpt]
382	END
383	lastx LE key-1:BEGIN ; on peut tout couper d''un coup
384	ncdf_varget,cdfid,name,res,offset=[jpi-key+ixmin+firstx $
385	,firsty+iymin $
386	,firstz+izmin $
387	,firsttps] $
388	,count=[nx,ny,nz,jpt]
389	END
390	ELSE:BEGIN ; Le tableau est separe en 2 morceaux...
391	ncdf_varget,cdfid,name,tab1,offset=[jpi-key+ixmin+firstx $
392	,firsty+iymin $
393	,firstz+izmin $
394	,firsttps] $
395	,count=[key-1-firstx+1,ny,nz,jpt]
396	ncdf_varget,cdfid,name,tab2,offset=[ixmin $
397	,firsty+iymin $
398	,firstz+izmin $
399	,firsttps] $
400	,count=[lastx-key+1,ny,nz,jpt]
401	res = [temporary(tab1), temporary(tab2)]
402	END
403	ENDCASE
404	END
405	;,,,,,,,,,,,,,,,,,,,,,,,,
406	keyword_set(key_shift) NE 0 AND total(key_stride) NE 3:BEGIN
407	; on shift et il y a un stride
408	;,,,,,,,,,,,,,,,,,,,,,,,,
409	case 1 OF ; case sur la facon de fire le champ
410	firstx*key_stride[0] GE key:BEGIN ; on peut tout coupe d''un coup
411	ncdf_varget,cdfid,name,res,offset=[ixmin $
412	+(key_stride[0]-1)-((key-1) MOD key_stride[0]) $
413	+(firstx-((key-1)/key_stride[0]+1))*key_stride[0] $
414	,firsty*key_stride[1]+iymin $
415	,firstz*key_stride[2]+izmin $
416	,firsttps] $
417	,count=[nx,ny,nz,jpt], stride = [key_stride, 1]
418	END
419	lastx*key_stride[0] LE key-1:BEGIN ; on peut tout coupe d''un coup
420	ncdf_varget,cdfid,name,res,offset=[jpitotal-key+ixmin+firstx*key_stride[0] $
421	,firsty*key_stride[1]+iymin $
422	,firstz*key_stride[2]+izmin $
423	,firsttps] $
424	,count=[nx,ny,nz,jpt], stride = [key_stride, 1]
425	END
426	ELSE:BEGIN ; le tableau est separe en 2 morceaux...
427	ncdf_varget,cdfid,name,tab1,offset=[jpitotal-key+ixmin+firstx*key_stride[0] $
428	,firsty*key_stride[1]+iymin $
429	,firstz*key_stride[2]+izmin $
430	,firsttps] $
431	,count=[(key-firstx*key_stride[0]-1)/key_stride[0]+1,ny,nz,jpt] $
432	, stride = [key_stride, 1]
433	; 2eme bout...
434	ncdf_varget,cdfid,name,tab2,offset=[ixmin $
435	+(key_stride[0]-1)-((key-1) MOD key_stride[0]) $
436	, firsty*key_stride[1]+iymin $
437	, firstz*key_stride[2]+izmin $
438	, firsttps] $
439	,count=[nx-(key-firstx*key_stride[0]-1)/key_stride[0]-1,ny,nz,jpt] $
440	, stride = [key_stride, 1]
441	; on recolle le tout
442	res = [temporary(tab1), temporary(tab2)]
443	END
444	ENDCASE ; case sur la facon de fire le champ
445	END
446	endcase
447	END
448	ENDCASE
449
450	; we apply reverse
451	IF keyword_set(key_yreverse) AND ny NE 1 THEN BEGIN
452	IF varcontient.ndims - ((where(varcontient.dim EQ contient.recdim))[0] NE -1) EQ 2 THEN $
453	res = reverse(reform(res, nx, ny, jpt, /overwrite), 2) $
454	ELSE res = reverse(reform(res, nx, ny, nz, jpt, /overwrite), 2)
455	ENDIF
456	if keyword_set(key_zreverse) AND nz NE 1 $
457	AND varcontient.ndims - ((where(varcontient.dim EQ contient.recdim))[0] NE -1) EQ 3 THEN $
458	res = reverse(reform(res, nx, ny, nz, jpt, /overwrite), 3)

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

Download in other formats:

Original Format