Context Navigation

source: trunk/condmag_on_orca.pro @ 10

Last change on this file since 10 was 10, checked in by pinsard, 17 years ago
passage sur rhodes et init adaptation ORCA025
Property svn:keyword set to `Id`
File size: 15.9 KB

Line
1	;+
2	;
3	; @file_comments
4	; interpolate condmag.nc file (sediment and magnetic fields) on ORCA grid
5	; and produce cond_sed_<orcares>.nc and Br_<orcares>.nc
6	;
7	; @categories
8	; interpolation, orca grid
9	;
10	; @param orcares {in}{required}{type=string}
11	; code of ORCA grid to use for output results
12	; must be 'ORCA2' or 'ORCA025'
13	;
14	; @param method {in}{required}{type=string}
15	; code for interpolation method
16	; must be 'bilinear' or 'imoms3'
17	; in fact bilinear is used in this geomag project
18	; ++ est-ce que ce sera tjs vrai ?
19	;
20	; @param gridtype {in}{required}{type=string}
21	; to specify on which grid we do the interpolation T, U, V
22	; must belong to T, U or V
23	;
24	; @restrictions
25	; - Requires SAXO tools
26	; - must have condmag.nc in ${GEOMAG_ID}/
27	; - must have meshmask for ORCA2 or ORCA025 in ${GEOMAG_ID}/
28	; - must not have cond_sed_*.nc in ${GEOMAG_OD}/
29	; - must not have Br_*.nc in ${GEOMAG_OD}/
30	;
31	; @todo
32	; provide tools to plot output files
33	; produce a NetCDF GDT or CF compliant
34	;
35	; @pre
36	; be sure to have datafile condmag.nc in the directory defined in
37	; ${GEOMAG_ID}/ see geomag_env.sh
38	; be sure to have meshmask of ORCA grid you choose in the directory defined in
39	; ${GEOMAG_ID}/
40	; for ORCA2 filename is meshmask_bab.nc
41	; ++ au pif entre mesh_hgr, mesh_zgr et mask
42	; for ORCA0252 filename is ORCA025-G42_mesh_hgr.nc
43	;
44	; @post
45	; cond_sed_<I>orcasres</I>.nc is now present in ${GEOMAG_OD}/
46	; Br_<I>orcasres</I>.nc is now present in ${GEOMAG_OD}/
47	; see geomag_env.sh
48	;
49	; @examples
50	; IDL> .run condmag_on_orca
51	; IDL> condmag_on_orca, 'ORCA2', 'bilinear','T'
52	;
53	; @history
54	; reee522 2006-11-23T13:34:31Z rhodes (IRIX64)
55	; add ORCA025 (beginning)
56	; reee522 2006-11-23T13:28:07Z rhodes (IRIX64)
57	; information about the opening of meshmask file was missing ... because
58	; the open is hidden in the call of get_gridparam
59	; fplod 2006-11-23T08:57:10Z aedon.locean-ipsl.upmc.fr (Darwin)
60	; rename to condmag_on_orca.pro
61	; fplod 2006-11-22T15:10:31Z aedon.locean-ipsl.upmc.fr (Darwin)
62	; add creation of Br_<I>orcasres</I>.nc
63	; ++ not very beautifuly implemented
64	; fplod 2006-11-22T10:38:51Z aedon.locean-ipsl.upmc.fr (Darwin)
65	; always use msg = 'iii : ...' or msg = 'eee : ...'
66	; information on open io
67	; fplod 2006-11-20T10:32:02Z aedon.locean-ipsl.upmc.fr (Darwin)
68	; cleaning before introduction into svn repository
69	; fplod 2006-03-23T13:05:39Z aedon.lodyc.jussieu.fr (Darwin)
70	; cond_sed_ORCA2 presque ok manque attributes min/max de cond_sed
71	; fplod 2006-03-22T09:35:42Z aedon.lodyc.jussieu.fr (Darwin)
72	; created from
73	; /Users/fplod/incas/seb/DORAEMON/wind/idl/quikscat_to_orca_scalar.pro written
74	; by Sebastien Masson
75	; to reproduce /usr/work/sur/fvi/OPA/geomag/cond_sed_ORCA2.nc
76	; main differences : no yyyy parameter, no time loop , no mask and no missing
77	; values (++ to be checked ) in data input (condmag.nc), no scale factor,
78	; no OFFSET, no save of weight and addresses
79	;
80	; @version
81	; $Id$
82	;
83	;-
84	PRO condmag_on_orca, orcares, method, gridtype
85	;
86	compile_opt idl2, strictarrsubs
87	;
88	;----
89	; check input parameters
90	;----
91	;
92	; check orcares definition
93	;
94	CASE orcares OF
95	'ORCA2': BEGIN
96	msg = 'iii : valid orcares parameter = ' + orcares
97	PRINT, msg
98	filename_oce = 'meshmask_bab.nc'
99	END
100	'ORCA025': BEGIN
101	msg = 'iii : valid orcares parameter = ' + orcares
102	PRINT, msg
103	filename_oce = orcares + '-G42_mesh_hgr.nc'
104	END
105	ELSE : BEGIN
106	msg = 'eee : invalid orcares parameter = ' + orcares
107	PRINT, msg
108	msg = 'eee : orcares must be ORCA2 or ORCA025'
109	PRINT, msg
110	RETURN
111	END
112	ENDCASE
113	;
114	; check method definition
115	CASE method OF
116	'bilinear': BEGIN
117	msg = 'iii : valid method parameter = ' + method
118	PRINT, msg
119	END
120	'imoms3' : BEGIN
121	msg = 'iii : valid method parameter = ' + method
122	PRINT, msg
123	END
124	ELSE : BEGIN
125	msg = 'eee : invalid method parameter = ' + method
126	PRINT, msg
127	msg = 'eee : method must be bilinear or imoms3'
128	PRINT, msg
129	RETURN
130	END
131	ENDCASE
132	;
133	; check gridtype definition
134	gridtype = strupcase(gridtype)
135	CASE gridtype OF
136	'T' : BEGIN
137	msg = 'iii : valid gridtype parameter = ' + gridtype
138	PRINT, msg
139	END
140	'U' : BEGIN
141	msg = 'iii : valid gridtype parameter = ' + gridtype
142	PRINT, msg
143	END
144	'V' : BEGIN
145	msg = 'iii : valid gridtype parameter = ' + gridtype
146	PRINT, msg
147	END
148	ELSE: BEGIN
149	msg = 'eee : invalid gridtype parameter = ' + gridtype
150	PRINT, msg
151	msg = 'eee : gridtype must be T, U or V'
152	PRINT, msg
153	RETURN
154	END
155	ENDCASE
156	;
157	; check for input files
158	;
159	; test if ${GEOMAG_ID} defined
160	geomag_id_env=GETENV('GEOMAG_ID')
161	CASE geomag_id_env OF
162	'' : BEGIN
163	msg = 'eee : ${GEOMAG_ID} is not defined'
164	PRINT, msg
165	RETURN
166	END
167	ELSE: BEGIN
168	msg = 'iii : ${GEOMAG_ID} is ' + geomag_id_env
169	PRINT, msg
170	END
171	ENDCASE
172	;
173	iodirin = isadirectory(geomag_id_env)
174	;
175	; existence and protection of ${GEOMAG_ID}
176	IF (FILE_TEST(iodirin, /DIRECTORY,/EXECUTABLE , /READ) EQ 0) THEN BEGIN
177	msg = 'eee : the directory' + iodirin + ' is not accessible.'
178	PRINT, msg
179	RETURN
180	ENDIF
181	;
182	; conductivity and magnetic field
183	filename_condmag = 'condmag.nc'
184	;
185	; check if this file exists
186	fullfilename_condmag = isafile(iodirin + filename_condmag, NEW=0,/MUST_EXIST)
187	IF fullfilename_condmag[0] EQ '' THEN BEGIN
188	msg = 'eee : the file ' + fullfilename_condmag + ' was not found.'
189	PRINT, msg
190	RETURN
191	ENDIF
192	;
193	; protection
194	IF (FILE_TEST(fullfilename_condmag[0], /READ) EQ 0) THEN BEGIN
195	msg = 'eee : the file ' + fullfilename_condmag[0] + ' is not readable.'
196	PRINT, msg
197	RETURN
198	ENDIF
199	;
200	; mesh mask
201	; check if this file exists
202	fullfilename_oce = isafile(iodirin + filename_oce, NEW=0,/MUST_EXIST, RECURSIVE=0)
203	IF fullfilename_oce[0] EQ '' THEN BEGIN
204	msg = 'eee : the file ' + fullfilename_oce + ' was not found.'
205	PRINT, msg
206	RETURN
207	ENDIF
208	;
209	; protection
210	IF (FILE_TEST(fullfilename_oce[0], /READ) EQ 0) THEN BEGIN
211	msg = 'eee : the file ' + fullfileoce_condmag[0] + ' is not readable.'
212	PRINT, msg
213	RETURN
214	ENDIF
215	;
216	; test if ${GEOMAG_OD} defined
217	geomag_od_env=GETENV('GEOMAG_OD')
218	CASE geomag_od_env OF
219	'' : BEGIN
220	msg = 'eee : ${GEOMAG_OD} is not defined'
221	PRINT, msg
222	RETURN
223	END
224	ELSE: BEGIN
225	msg = 'iii : ${GEOMAG_OD} is ' + geomag_od_env
226	PRINT, msg
227	END
228	ENDCASE
229	;
230	; check if output data will be possible
231	iodirout = isadirectory(geomag_od_env)
232	;
233	; existence and protection
234	IF (FILE_TEST(iodirout, /DIRECTORY,/WRITE) EQ 0) THEN BEGIN
235	msg = 'eee : the directory' + iodirout + ' was not found.'
236	PRINT, msg
237	RETURN
238	ENDIF
239	;
240	; build output filenames
241	filename_cond_sed = 'cond_sed' + '_' + orcares +'.nc'
242	fullfilename_cond_sed = iodirout + filename_cond_sed
243	;
244	; in order to avoid unexpected overwritten
245	IF (FILE_TEST(fullfilename_cond_sed) EQ 1) THEN BEGIN
246	msg = 'eee : the file ' + fullfilename_cond_sed + ' already exists.'
247	PRINT, msg
248	RETURN
249	ENDIF
250	;
251	filename_br = 'Br' + '_' + orcares +'.nc'
252	fullfilename_br = iodirout + filename_br
253	;
254	; in order to avoid unexpected overwritten
255	IF (FILE_TEST(fullfilename_br) EQ 1) THEN BEGIN
256	msg = 'eee : the file ' + fullfilename_br + ' already exists.'
257	PRINT, msg
258	RETURN
259	ENDIF
260	;
261	; d'après ncdump -h /usr/work/sur/fvi/OPA/geomag/condmag.nc
262	condmaglonname = 'lo'
263	condmaglatname = 'la'
264	varname_cond_sed = 'cond_sed'
265	varname_br = 'Br'
266	;
267	;----
268	; conductivity and magnetic field grid parameters
269	;----
270	;
271	get_gridparams, fullfilename_condmag[0], $
272	condmaglonname, condmaglatname, $
273	condmaglon, condmaglat, jpia, jpja, 1,/DOUBLE
274	;
275	;----
276	; Oceanic grid parameters
277	;----
278	;
279	olonname = 'glam' + STRLOWCASE(gridtype)
280	olatname = 'gphi' + STRLOWCASE(gridtype)
281	get_gridparams, fullfilename_oce[0], $
282	olonname, olatname, $
283	olon, olat, jpio, jpjo, 2,/DOUBLE
284	msg = 'iii : ' + fullfilename_oce[0] + ' opened for read'
285	PRINT, msg
286	;
287	;---------------
288	; Compute weight and address
289	;---------------
290	;
291	CASE method OF
292	'bilinear': compute_fromreg_bilinear_weigaddr, $
293	condmaglon, condmaglat, olon, olat, weig, addr
294	'imoms3' : compute_fromreg_imoms3_weigaddr, $
295	condmaglon, condmaglat, olon, olat, weig, addr
296	ENDCASE
297	;
298	; reading condmag file
299	netcdf_id_condmag = NCDF_OPEN(fullfilename_condmag[0], /NOWRITE)
300	msg = 'iii : ' + fullfilename_condmag[0] + ' opened for read'
301	PRINT, msg
302
303	varinq_cond_sed = NCDF_VARINQ(netcdf_id_condmag, varname_cond_sed)
304	varinq_br = NCDF_VARINQ(netcdf_id_condmag, varname_br)
305	;
306	;---------------------------
307	; Creation of the NetCdf file for cond_sed and Br
308	;---------------------------
309	;
310	netcdf_id_cond_sed = NCDF_CREATE(fullfilename_cond_sed, /clobber)
311	NCDF_CONTROL, netcdf_id_cond_sed, /NOFILL
312	netcdf_id_br = NCDF_CREATE(fullfilename_br, /clobber)
313	NCDF_CONTROL, netcdf_id_br, /NOFILL
314	;
315	; dimension
316	dimidx = NCDF_DIMDEF(netcdf_id_cond_sed, 'x' , jpio)
317	dimidy = NCDF_DIMDEF(netcdf_id_cond_sed, 'y' , jpjo)
318	dimidt = NCDF_DIMDEF(netcdf_id_cond_sed, 'lo', /UNLIMITED)
319	dimidx = NCDF_DIMDEF(netcdf_id_br, 'x' , jpio)
320	dimidy = NCDF_DIMDEF(netcdf_id_br, 'y' , jpjo)
321	dimidt = NCDF_DIMDEF(netcdf_id_br, 'lo', /UNLIMITED)
322	;
323	; global attributes
324	NCDF_ATTPUT, netcdf_id_cond_sed, 'Conventions', 'GDT 1.2', /GLOBAL ; ++ conformité
325	NCDF_ATTPUT, netcdf_id_cond_sed, 'file_name' , filename_cond_sed, /GLOBAL
326	NCDF_ATTPUT, netcdf_id_cond_sed, 'Title' , 'Conductance', /GLOBAL
327	NCDF_ATTPUT, netcdf_id_br, 'Conventions', 'GDT 1.2', /GLOBAL ; ++ conformité
328	NCDF_ATTPUT, netcdf_id_br, 'file_name' , filename_Br, /GLOBAL
329	NCDF_ATTPUT, netcdf_id_br, 'Title' , 'Magnetic field', /GLOBAL
330	;
331	; declaration of variables
332	; 4 common variables for the two files to produce
333	varid = lonarr(3)
334	;
335	varid[0] = NCDF_VARDEF(netcdf_id_cond_sed, 'nav_lon' , [dimidx, dimidy], /FLOAT)
336	NCDF_ATTPUT, netcdf_id_cond_sed, varid[0], 'units' , 'degrees_east'
337	NCDF_ATTPUT, netcdf_id_cond_sed, varid[0], 'valid_min', min(olon, max = omax),/FLOAT
338	NCDF_ATTPUT, netcdf_id_cond_sed, varid[0], 'valid_max', omax,/FLOAT
339	NCDF_ATTPUT, netcdf_id_cond_sed, varid[0], 'long_name', 'Longitude at t-point'
340	varid[0] = NCDF_VARDEF(netcdf_id_br, 'nav_lon' , [dimidx, dimidy], /FLOAT)
341	NCDF_ATTPUT, netcdf_id_br, varid[0], 'units' , 'degrees_east'
342	NCDF_ATTPUT, netcdf_id_br, varid[0], 'valid_min', min(olon, max = omax),/FLOAT
343	NCDF_ATTPUT, netcdf_id_br, varid[0], 'valid_max', omax,/FLOAT
344	NCDF_ATTPUT, netcdf_id_br, varid[0], 'long_name', 'Longitude at t-point'
345	;
346	varid[1] = NCDF_VARDEF(netcdf_id_cond_sed, 'nav_lat' , [dimidx, dimidy], /FLOAT)
347	NCDF_ATTPUT, netcdf_id_cond_sed, varid[1], 'units' , 'degrees_north'
348	NCDF_ATTPUT, netcdf_id_cond_sed, varid[1], 'valid_min', min(olat, max = omax),/FLOAT
349	NCDF_ATTPUT, netcdf_id_cond_sed, varid[1], 'valid_max', omax,/FLOAT
350	NCDF_ATTPUT, netcdf_id_cond_sed, varid[1], 'long_name', 'Latitude at t-point'
351	varid[1] = NCDF_VARDEF(netcdf_id_br, 'nav_lat' , [dimidx, dimidy], /FLOAT)
352	NCDF_ATTPUT, netcdf_id_br, varid[1], 'units' , 'degrees_north'
353	NCDF_ATTPUT, netcdf_id_br, varid[1], 'valid_min', min(olat, max = omax),/FLOAT
354	NCDF_ATTPUT, netcdf_id_br, varid[1], 'valid_max', omax,/FLOAT
355	NCDF_ATTPUT, netcdf_id_br, varid[1], 'long_name', 'Latitude at t-point'
356	;
357	varid[2] = NCDF_VARDEF(netcdf_id_cond_sed, 'time' , [dimidt], /FLOAT)
358	varid[2] = NCDF_VARDEF(netcdf_id_br, 'time' , [dimidt], /FLOAT)
359	;
360	; each of the two files to produce contains a specific variable
361	;
362	varid_cond_sed = lonarr(1)
363	varid_cond_sed[0] = NCDF_VARDEF(netcdf_id_cond_sed, varname_cond_sed, [dimidx, dimidy, dimidt], /FLOAT)
364	;
365	; variable 3
366	NCDF_ATTPUT, netcdf_id_cond_sed, varid_cond_sed[0], 'long_name', 'Conductance' ; ++ non cf
367	NCDF_ATTPUT, netcdf_id_cond_sed, varid_cond_sed[0], 'units', 'siemens' ; ++ récupérer l'unite de cond_sed dans condmag.nc
368	; pour min et max, il faut avoir lu la variable ... cf. plus bas ++
369	; donc pour l'instant on les met à valeur manquante
370	NCDF_ATTPUT, netcdf_id_cond_sed, varid_cond_sed[0], 'valid_min', !VALUES.F_NAN ,/FLOAT
371	NCDF_ATTPUT, netcdf_id_cond_sed, varid_cond_sed[0], 'valid_max', !VALUES.F_NAN ,/FLOAT
372	;
373	varid_br = lonarr(1)
374	varid_br[0] = NCDF_VARDEF(netcdf_id_br, varname_br, [dimidx, dimidy, dimidt], /FLOAT)
375	;
376	; variable 3
377	NCDF_ATTPUT, netcdf_id_br, varid_br[0], 'long_name', 'magnetic field' ; ++ non cf
378	NCDF_ATTPUT, netcdf_id_br, varid_br[0], 'units', 'tesla' ; ++ récupérer l'unité de Br dans condmag.nc
379	; pour min et max, il faut avoir lu la variable ... cf. plus bas ++
380	; donc pour l'instant on les met à valeur manquante
381	NCDF_ATTPUT, netcdf_id_br, varid_br[0], 'valid_min', !VALUES.F_NAN ,/FLOAT
382	NCDF_ATTPUT, netcdf_id_cond_sed, varid_cond_sed[0], 'valid_max', !VALUES.F_NAN ,/FLOAT
383	;
384	;---------------------------
385	; end of header definition, writing of the NetCdf files
386	;---------------------------
387	NCDF_CONTROL, netcdf_id_cond_sed, /ENDEF
388	NCDF_CONTROL, netcdf_id_br, /ENDEF
389	NCDF_CLOSE, netcdf_id_condmag ; ++ si le close n'est pas ici ça merde mais je ne sais pas pourquoi
390	;
391	; grid
392	NCDF_VARPUT, netcdf_id_cond_sed, 'nav_lon', olon
393	NCDF_VARPUT, netcdf_id_cond_sed, 'nav_lat', olat
394	NCDF_VARPUT, netcdf_id_cond_sed, varid[2], FLOAT(0.5) ; ++ c'est quoi cette valeur
395	NCDF_VARPUT, netcdf_id_br, 'nav_lon', olon
396	NCDF_VARPUT, netcdf_id_br, 'nav_lat', olat
397	NCDF_VARPUT, netcdf_id_br, varid[2], FLOAT(0.5) ; ++ c'est quoi cette valeur
398
399	;---------------------------
400	; réouverture du fichier ... mais pourquoi on le relit !!++
401	netcdf_id_condmag = NCDF_OPEN(fullfilename_condmag[0],/NOWRITE)
402	msg = 'iii : ' + fullfilename_condmag[0] + ' reopened for read'
403	PRINT, msg
404	;
405	NCDF_VARGET, netcdf_id_condmag, varname_cond_sed, varin_cond_sed ; , COUNT = [jpia, jpja, 1], OFFSET = [0, 0, 0]
406	NCDF_VARGET, netcdf_id_condmag, varname_br, varin_br ; , COUNT = [jpia, jpja, 1], OFFSET = [0, 0, 0]
407	;
408	; do the interpolation
409	varin_cond_sed = TOTAL(weig*varin_cond_sed[addr], 1)
410	varin_cond_sed = REFORM(varin_cond_sed, jpio, jpjo, /OVER)
411	varin_br = TOTAL(weig*varin_br[addr], 1)
412	varin_br = REFORM(varin_br, jpio, jpjo, /OVER)
413	;
414	NCDF_CLOSE, netcdf_id_condmag
415	varout_cond_sed = TEMPORARY(varin_cond_sed)
416	varout_br = TEMPORARY(varin_br)
417	; compute min max
418	minarr_cond_sed = min(varout_cond_sed, max = maxarr_cond_sed)
419	minarr_br = min(varout_br, max = maxarr_br)
420	;
421	; put back the masked value
422	;++ IF bad[0] NE -1 THEN varout_cond_sed[TEMPORARY(bad)] = 32767
423	;++ IF bad[0] NE -1 THEN varout_br[TEMPORARY(bad)] = 32767
424	;
425	; write the data
426	NCDF_VARPUT, netcdf_id_cond_sed, varname_cond_sed, varout_cond_sed, COUNT = [jpio, jpjo, 1], OFFSET = [0, 0, 0]
427	NCDF_VARPUT, netcdf_id_br, varname_br, varout_br, COUNT = [jpio, jpjo, 1], OFFSET = [0, 0, 0]
428	;
429	; update min max attributes
430	NCDF_CONTROL, netcdf_id_cond_sed, /REDEF
431	NCDF_ATTPUT, netcdf_id_cond_sed, varid_cond_sed[0], 'valid_min', minarr_cond_sed,/FLOAT
432	NCDF_ATTPUT, netcdf_id_cond_sed, varid_cond_sed[0], 'valid_max', maxarr_cond_sed,/FLOAT
433	NCDF_CONTROL, netcdf_id_cond_sed, /ENDEF
434	NCDF_CONTROL, netcdf_id_br, /REDEF
435	NCDF_ATTPUT, netcdf_id_br, varid_br[0], 'valid_min', minarr_br,/FLOAT
436	NCDF_ATTPUT, netcdf_id_br, varid_br[0], 'valid_max', maxarr_br,/FLOAT
437	NCDF_CONTROL, netcdf_id_br, /ENDEF
438	;
439	;---------------------------
440	; close the netcdf files
441	NCDF_CLOSE, netcdf_id_cond_sed
442	NCDF_CLOSE, netcdf_id_br
443	;
444	msg = 'iii : ' + fullfilename_cond_sed + ' created'
445	PRINT, msg
446	msg = 'iii : ' + fullfilename_br + ' created'
447	PRINT, msg
448	;
449	END

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

Download in other formats:

Original Format