source: trunk/step2_diff.pro @ 48

Last change on this file since 48 was 41, checked in by pinsard, 15 years ago

small headers improvements
  • Property svn:keyword set to ID
File size: 7.1 KB
Line 
1;+
2;
3; @file_comments
4; compute delta between *_5d_yyyy_grid_T_orcares.pro Netcdf files with same dimension
5;
6; @param file1 {in}{required} {type=string}
7; name of the first file to be read
8;
9; @param file2 {in}{required} {type=string}
10; name of the second file to be read
11;
12; @param file3 {in}{required} {type=string}
13; file where differences between variables are written
14;
15; @examples
16; test (for header and delta)
17; to compute difference between Sigma in Sigma_5d_1993_grid_T_ORCA2.nc and
18; itself and write it in ginette.nc :
19; $ cd ${GEOMAG}
20; $ idl
21; IDL> file1=getenv('GEOMAG_OD') + '/Sigma_5d_1993_grid_T_ORCA2.nc'
22; IDL> file2=getenv('GEOMAG_OD') + '/Sigma_5d_1993_grid_T_ORCA2.nc'
23; IDL> step2_diff, file1, file2, 'ginette.nc'
24; values of delta must be 0 everywhere
25; $ rm ginette.nc
26;
27; real life
28; to compute difference between Sigma in
29; /usr/work/sur/fvi/OPA/ORCA2/Sigma_5d_1994_grid_T.nc and
30; ${GEOMGA_OD}/Sigma_5d_1994_grid_T_ORCA2.nc and write it in
31; Sigma_5d_1994_grid_T_ORCA2_diff.nc :
32; $ cd ${GEOMAG}
33; $ idl
34; IDL> year=1994
35; IDL> file1= '/usr/work/sur/fvi/OPA/ORCA2/Sigma_5d_1994_grid_T.nc'
36; IDL> file2=getenv('GEOMAG_OD') + '/Sigma_5d_1994_grid_T_ORCA2.nc'
37; IDL> step2_diff, file1, file2, 'Sigma_5d_1994_grid_T_ORCA2_diff.nc'
38;
39; to see the difference, for example
40; IDL> xxx,'Sigma_5d_1994_grid_T_ORCA2_diff.nc'
41; select xy on plt wigdet
42;
43; see also step2_diff.sh for quick start
44;
45; @history
46; fplod 2007-08-22T09:39:45Z aedon.locean-ipsl.upmc.fr (Darwin)
47; correction
48; fplod 2007-07-30T10:58:12Z aedon.locean-ipsl.upmc.fr (Darwin)
49; like <progeomag>step1_diff</progeomag>
50;
51; @todo
52; check difference between headers
53;
54; @version
55; $Id$
56;
57PRO step2_diff, file1, file2, file3
58;
59ncverbose=1
60;
61; check if file1 exists
62fullfile1 = isafile(file1, NEW=0,/MUST_EXIST)
63IF fullfile1[0] EQ '' THEN BEGIN
64   msg = 'eee : the file ' + fullfile1 + ' was not found.'
65   ras = report(msg)
66   RETURN
67ENDIF
68;
69; protection
70IF (FILE_TEST(fullfile1[0], /READ) EQ 0) THEN BEGIN
71   msg = 'eee : the file ' + fullfile1[0] + ' is not readable.'
72   ras = report(msg)
73   RETURN
74ENDIF
75;
76; check if file2 exists
77fullfile2 = isafile(file2, NEW=0,/MUST_EXIST)
78IF fullfile2[0] EQ '' THEN BEGIN
79   msg = 'eee : the file ' + fullfile2 + ' was not found.'
80   ras = report(msg)
81   RETURN
82ENDIF
83;
84; protection
85IF (FILE_TEST(fullfile2[0], /READ) EQ 0) THEN BEGIN
86   msg = 'eee : the file ' + fullfile2[0] + ' is not readable.'
87   ras = report(msg)
88   RETURN
89ENDIF
90;
91; read the first file
92ncdf_read,  fullfile1[0], file1info, file1dinfo, file1vinfo, $
93            file1gatts, file1vatts, file1data
94;
95; read the second file
96ncdf_read,  fullfile2[0], file2info, file2dinfo, file2vinfo, $
97            file2gatts, file2vatts, file2data
98;
99; compute delta
100delta = *file1data[4] - *file2data[4]
101minarr = min(delta, max = maxarr)
102;
103; check for number of non zero in delta
104checknonzero=where(delta NE 0.,count)
105IF count EQ 0 THEN BEGIN
106   msg = 'iii : delta is zero everywhere'
107   ras = report(msg)
108ENDIF ELSE BEGIN
109   msg = 'iii : delta is not zero ' + STRING(count)  + ' times'
110   ras = report(msg)
111ENDELSE
112;
113; write output
114; in order to avoid unexpected overwritten
115IF (FILE_TEST(file3) EQ 1) THEN BEGIN
116   msg = 'eee : the file ' + file3  + ' already exists.'
117   ras = report(msg)
118   RETURN
119ENDIF
120;
121;---------------------------
122; Creation of the NetCdf output file
123;---------------------------
124;
125  netcdf_id = NCDF_CREATE(file3, /clobber)
126  NCDF_CONTROL, netcdf_id, /NOFILL
127;
128; dimension
129  dimidx = NCDF_DIMDEF(netcdf_id, file1dinfo[0].name ,  file1dinfo[0].size)
130  dimidy = NCDF_DIMDEF(netcdf_id, file1dinfo[1].name ,  file1dinfo[1].size)
131  dimidz = NCDF_DIMDEF(netcdf_id, file1dinfo[2].name ,  file1dinfo[2].size)
132  dimidt = NCDF_DIMDEF(netcdf_id, file1dinfo[3].name, /UNLIMITED)
133  dimidz2 = NCDF_DIMDEF(netcdf_id, file1dinfo[4].name ,  file1dinfo[4].size)
134;
135; global attributes
136  NCDF_ATTPUT, netcdf_id, 'Conventions', 'GDT 1.2', /GLOBAL ; ++ conformite
137  NCDF_ATTPUT, netcdf_id, 'file_name'  , file3, /GLOBAL
138  NCDF_ATTPUT, netcdf_id, 'Title'      , $
139               'delta between '+ file1 + ' and ' + file2, /GLOBAL
140;
141; declaration of variables
142  varid = lonarr(5)
143;
144  varid[0] = NCDF_VARDEF(netcdf_id, file1vinfo[0].name  , [dimidx, dimidy], /FLOAT)
145  NCDF_ATTPUT, netcdf_id, varid[0], $
146               (*file1vatts[0])[0].name , $
147               STRING(*(*file1vatts[0])[0].values)
148  NCDF_ATTPUT, netcdf_id, varid[0], $
149               (*file1vatts[0])[1].name, $
150               STRING(*(*file1vatts[0])[1].values)
151  NCDF_ATTPUT, netcdf_id, varid[0], $
152               (*file1vatts[0])[2].name, $
153               STRING(*(*file1vatts[0])[2].values)
154;
155  varid[1] = NCDF_VARDEF(netcdf_id, file1vinfo[1].name  , [dimidx, dimidy], /FLOAT)
156  NCDF_ATTPUT, netcdf_id, varid[1], $
157               (*file1vatts[1])[0].name, $
158               STRING(*(*file1vatts[1])[0].values)
159  NCDF_ATTPUT, netcdf_id, varid[1], $
160               (*file1vatts[1])[1].name, $
161               STRING(*(*file1vatts[1])[1].values)
162  NCDF_ATTPUT, netcdf_id, varid[1], $
163               (*file1vatts[1])[2].name, $
164               STRING(*(*file1vatts[1])[2].values)
165;
166  varid[2] = NCDF_VARDEF(netcdf_id, file1vinfo[2].name, [dimidz2], /FLOAT);
167  NCDF_ATTPUT, netcdf_id, varid[2], $
168               (*file1vatts[2])[0].name, $
169               STRING(*(*file1vatts[2])[0].values)
170  NCDF_ATTPUT, netcdf_id, varid[2], $
171               (*file1vatts[2])[1].name, $
172               STRING(*(*file1vatts[2])[1].values)
173  NCDF_ATTPUT, netcdf_id, varid[2], $
174               (*file1vatts[2])[2].name, $
175               STRING(*(*file1vatts[2])[2].values)
176;
177  varid[3] = NCDF_VARDEF(netcdf_id, file1vinfo[3].name, [dimidt], /FLOAT);
178  NCDF_ATTPUT, netcdf_id, varid[3], $
179               (*file1vatts[3])[0].name, $
180               STRING(*(*file1vatts[3])[0].values)
181  NCDF_ATTPUT, netcdf_id, varid[3], $
182               (*file1vatts[3])[1].name, $
183               STRING(*(*file1vatts[3])[1].values)
184  NCDF_ATTPUT, netcdf_id, varid[3], $
185               (*file1vatts[3])[2].name, $
186               STRING(*(*file1vatts[3])[2].values)
187  NCDF_ATTPUT, netcdf_id, varid[3], $
188               (*file1vatts[3])[3].name, $
189               STRING(*(*file1vatts[3])[3].values)
190  NCDF_ATTPUT, netcdf_id, varid[3], $
191               (*file1vatts[3])[4].name, $
192               STRING(*(*file1vatts[3])[4].values)
193help,delta
194
195  varid[4] =  NCDF_VARDEF(netcdf_id, 'delta',   [dimidx, dimidy,dimidz, dimidt], /DOUBLE)
196  NCDF_ATTPUT, netcdf_id, varid[4], 'long_name', 'delta'
197  NCDF_ATTPUT, netcdf_id, varid[4], 'valid_min', minarr ,/DOUBLE
198  NCDF_ATTPUT, netcdf_id, varid[4], 'valid_max', maxarr ,/DOUBLE
199;
200;---------------------------
201; end of header definition, writing of the NetCdf files
202;---------------------------
203;
204  NCDF_CONTROL, netcdf_id, /ENDEF
205;
206  NCDF_VARPUT, netcdf_id, file1vinfo[0].name, (*file1data[0])
207  NCDF_VARPUT, netcdf_id, file1vinfo[1].name, (*file1data[1])
208  NCDF_VARPUT, netcdf_id, file1vinfo[2].name, (*file1data[2])
209  NCDF_VARPUT, netcdf_id, file1vinfo[3].name, (*file1data[3])
210;
211; write the data
212  NCDF_VARPUT, netcdf_id, 'delta', delta
213;
214;---------------------------
215; close the netcdf files
216  NCDF_CLOSE, netcdf_id
217;
218  msg = 'iii : ' + file3 + ' created'
219  ras = report(msg)
220;
221END
Note: See TracBrowser for help on using the repository browser.