1 | ;+ |
---|
2 | ; |
---|
3 | ; @file_comments |
---|
4 | ; get the land/sea mask array from a NetCDF file |
---|
5 | ; |
---|
6 | ; @categories |
---|
7 | ; Read NetCDF file |
---|
8 | ; |
---|
9 | ; @param fileid {in}{required}{type=salar string or long} |
---|
10 | ; if fileid is a scalar string then it is the name of the file (with |
---|
11 | ; the full path) to be opened (in that case, the file will be opened |
---|
12 | ; and closed within ncdf_getmask). |
---|
13 | ; |
---|
14 | ; if fileid is a scalar then it is the id of the file return by a call |
---|
15 | ; to ncdf_open outside of ncdf_getmask (in that case, the file will |
---|
16 | ; NOT be opened and closed within ncdf_getmask) |
---|
17 | ; |
---|
18 | ; @keyword ADDSCL_BEFORE {default=0}{type=scalar: 0 or 1} |
---|
19 | ; put 1 to apply add_offset and scale factor on data before looking for |
---|
20 | ; missing values when using USEASMASK keyword |
---|
21 | ; |
---|
22 | ; @keyword INVMASK {default=0}{type=scalar: 0 or 1} |
---|
23 | ; Inverse the land/sea mask (that should have 0/1 values for land/sea): mask = 1-mask |
---|
24 | ; |
---|
25 | ; @keyword MASKNAME {type=string} |
---|
26 | ; A string giving the name of the variable in the file |
---|
27 | ; that contains the land/sea mask |
---|
28 | ; |
---|
29 | ; @keyword MISSING_VALUE {type=scalar} |
---|
30 | ; To define (or redefine if the attribute is already existing) the |
---|
31 | ; missing values used with USEASMASK keyword. Note that this value is |
---|
32 | ; not used if TESTOP keyword is given and contains 2 words. |
---|
33 | ; |
---|
34 | ; @keyword USEASMASK {type=scalar string} |
---|
35 | ; A string giving the name of the variable in the file |
---|
36 | ; that will be used to build the land/sea mask. In this case the |
---|
37 | ; mask is based on the first record (if record dimension |
---|
38 | ; exists). The mask is build according to operator defined by TESTOP |
---|
39 | ; keyword (default NE) and the testing values defined as |
---|
40 | ; 1) the second word of TESTOP if existing |
---|
41 | ; 2) MISSING_VALUE keyword |
---|
42 | ; 3) attribute missing_value or _fillvalue of the variable USEASMASK |
---|
43 | ; 4) !Values.f_nan (can be used only with NE and EQ operators) |
---|
44 | ; |
---|
45 | ; @keyword TESTOP {default='NE'} {type=scalar string, for example 'GT 0.5'} |
---|
46 | ; a string describing the type of test that will be done to define the |
---|
47 | ; mask. The test is performed on the variable specified by USEASMASK |
---|
48 | ; keyword. |
---|
49 | ; |
---|
50 | ; TESTOP can contain 1 or 2 words. The first word is the operator |
---|
51 | ; definition: "EQ" "NE" "GE" "GT" "LE" "LT" (default is NE). The |
---|
52 | ; second word define the testing value. If TESTOP contains only 1 |
---|
53 | ; word, then the test value is denifed by |
---|
54 | ; 1) MISSING_VALUE keyword |
---|
55 | ; 2) attribute missing_value or _fillvalue of the variable USEASMASK |
---|
56 | ; 3) !Values.f_nan (can be used only with NE and EQ operators) |
---|
57 | ; |
---|
58 | ; @keyword XMINMESH {default=0L}{type=scalar} |
---|
59 | ; Define common (cm_4mesh) variables ixminmesh used to define the localization |
---|
60 | ; of the first point of the grid along the x direction in a zoom of the original grid |
---|
61 | ; |
---|
62 | ; @keyword YMINMESH {default=0L}{type=scalar} |
---|
63 | ; Define common (cm_4mesh) variables iyminmesh used to define the localization |
---|
64 | ; of the first point of the grid along the y direction in a zoom of the original grid |
---|
65 | ; |
---|
66 | ; @keyword ZMINMESH {default=0L}{type=scalar} |
---|
67 | ; Define common (cm_4mesh) variables izminmesh used to define the localization |
---|
68 | ; of the first point of the grid along the z direction in a zoom of the original grid |
---|
69 | ; |
---|
70 | ; @keyword XMAXMESH {default=jpiglo-1}{type=scalar} |
---|
71 | ; Define common (cm_4mesh) variables ixmaxmesh used to define the localization |
---|
72 | ; of the last point of the grid along the x direction in a zoom of the original grid |
---|
73 | ; Note that if XMAXMESH < 0 then ixmaxmesh is defined as ixmaxmesh = jpiglo -1 + xmaxmesh |
---|
74 | ; |
---|
75 | ; @keyword YMAXMESH {default=jpjglo-1}{type=scalar} |
---|
76 | ; Define common (cm_4mesh) variables iymaxmesh used to define the localization |
---|
77 | ; of the last point of the grid along the y direction in a zoom of the original grid |
---|
78 | ; Note that if YMAXMESH < 0 then iymaxmesh is defined as iymaxmesh = jpjglo -1 + ymaxmesh |
---|
79 | ; |
---|
80 | ; @keyword ZMAXMESH {default=jpkglo-1}{type=scalar} |
---|
81 | ; Define common (cm_4mesh) variables izmaxmesh used to define the localization |
---|
82 | ; of the last point of the grid along the z direction in a zoom of the original grid |
---|
83 | ; Note that if ZMAXMESH < 0 then izmaxmesh is defined as izmaxmesh = jpkglo -1 + maxmesh |
---|
84 | ; |
---|
85 | ; @keyword |
---|
86 | ; _EXTRA to be able to call ncdf_getmask with _extra keyword |
---|
87 | ; |
---|
88 | ; @returns |
---|
89 | ; the land/sea mask 2D or 3D array or -1 in case of error or mask absence |
---|
90 | ; |
---|
91 | ; @examples |
---|
92 | ; |
---|
93 | ; IDL> mask = ncdf_getmask('HadISST1_1m_187001_200702_sst_reg1m.nc',useasmask = 'sst', missing_value = -1.00000e+30) |
---|
94 | ; |
---|
95 | ; IDL> mask = ncdf_getmask('meshmaskORCA2.nc', maskname = 'tmask') |
---|
96 | ; |
---|
97 | ; IDL> mask = ncdf_getmask('t106.nc', useasmask = 'SLM', testop = 'le 0.5') |
---|
98 | ; |
---|
99 | ; @history |
---|
100 | ; August 2007: Sebastien Masson (smasson\@lodyc.jussieu.fr) |
---|
101 | ; |
---|
102 | ; @version |
---|
103 | ; $Id$ |
---|
104 | ; |
---|
105 | ;- |
---|
106 | FUNCTION ncdf_getmask, fileid, ADDSCL_BEFORE = addscl_before $ |
---|
107 | , MASKNAME = maskname, USEASMASK = useasmask $ |
---|
108 | , MISSING_VALUE = missing_value, INVMASK = invmask $ |
---|
109 | , XMINMESH = xminmesh, XMAXMESH = xmaxmesh $ |
---|
110 | , YMINMESH = yminmesh, YMAXMESH = ymaxmesh $ |
---|
111 | , ZMINMESH = zminmesh, ZMAXMESH = zmaxmesh $ |
---|
112 | , TESTOP = testop, _EXTRA = ex |
---|
113 | ; |
---|
114 | compile_opt idl2, strictarrsubs |
---|
115 | ; |
---|
116 | @cm_general ; needed for iodir |
---|
117 | ; |
---|
118 | IF n_elements(xminmesh) NE 0 THEN ixminmesh = 0L > long(xminmesh[0]) ELSE ixminmesh = 0l |
---|
119 | IF n_elements(yminmesh) NE 0 THEN iyminmesh = 0L > long(yminmesh[0]) ELSE iyminmesh = 0l |
---|
120 | IF n_elements(zminmesh) NE 0 THEN izminmesh = 0L > long(zminmesh[0]) ELSE izminmesh = 0l |
---|
121 | ; |
---|
122 | IF NOT (keyword_set(maskname) OR keyword_set(useasmask)) AND keyword_set(romsgrid) THEN maskname = 'mask_rho' |
---|
123 | IF NOT (keyword_set(maskname) OR keyword_set(useasmask)) THEN return, -1 |
---|
124 | ;---------------------------------------------------------- |
---|
125 | ; should we open the file? |
---|
126 | IF size(fileid, /type) EQ 7 THEN $ |
---|
127 | cdfid = ncdf_open(isafile(fileid, title = 'which file must be open by ncdf_getmask?', IODIRECTORY = iodir, _extra = ex)) $ |
---|
128 | ELSE cdfid = fileid |
---|
129 | ; what is inside the file |
---|
130 | inq = ncdf_inquire(cdfid) |
---|
131 | ; name of the variables |
---|
132 | namevar = strarr(inq.nvars) |
---|
133 | for varid = 0, inq.nvars-1 do begin |
---|
134 | invar = ncdf_varinq(cdfid, varid) |
---|
135 | namevar[varid] = strlowcase(invar.name) |
---|
136 | ENDFOR |
---|
137 | ;---------------------------------------------------------- |
---|
138 | CASE 1 OF |
---|
139 | keyword_set(maskname):mskid = (where(namevar EQ strlowcase(maskname)))[0] |
---|
140 | keyword_set(useasmask):mskid = (where(namevar EQ strlowcase(useasmask)))[0] |
---|
141 | ENDCASE |
---|
142 | ; |
---|
143 | if mskid NE -1 THEN BEGIN |
---|
144 | mskinq = ncdf_varinq(cdfid, mskid) |
---|
145 | ; is the mask variable containing the record dimension? |
---|
146 | withrcd = (where(mskinq.dim EQ inq.recdim))[0] |
---|
147 | ; in order to read only the first record |
---|
148 | ; we need to get the size of each dimension |
---|
149 | offset = lonarr(mskinq.ndims) |
---|
150 | count = replicate(1L, mskinq.ndims) |
---|
151 | FOR d = 0, mskinq.ndims -1 DO BEGIN |
---|
152 | IF d NE withrcd THEN BEGIN |
---|
153 | ncdf_diminq, cdfid, mskinq.dim[d], name, size |
---|
154 | count[d] = size |
---|
155 | ENDIF |
---|
156 | ENDFOR |
---|
157 | |
---|
158 | IF n_elements(xmaxmesh) NE 0 THEN BEGIN |
---|
159 | IF xmaxmesh GE 0 THEN ixmaxmesh = long(xmaxmesh[0]) ELSE ixmaxmesh = count[0] - 1L + long(xmaxmesh[0]) |
---|
160 | ENDIF ELSE ixmaxmesh = count[0] - 1L |
---|
161 | ixmaxmesh = 0 > ixmaxmesh < (count[0] - 1L) ; make sure ixmaxmesh is not too big |
---|
162 | jpiglo = count[0] |
---|
163 | offset[0] = ixminmesh |
---|
164 | count[0] = ixmaxmesh - ixminmesh + 1L |
---|
165 | jpi = count[0] |
---|
166 | |
---|
167 | IF n_elements(ymaxmesh) NE 0 THEN BEGIN |
---|
168 | IF ymaxmesh GE 0 THEN iymaxmesh = long(ymaxmesh[0]) ELSE iymaxmesh = count[1] - 1L + long(ymaxmesh[0]) |
---|
169 | ENDIF ELSE iymaxmesh = count[1] - 1L |
---|
170 | iymaxmesh = 0 > iymaxmesh < (count[1] - 1L) ; make sure ixmaxmesh is not too big |
---|
171 | jpjglo = count[1] |
---|
172 | offset[1] = iyminmesh |
---|
173 | count[1] = iymaxmesh - iyminmesh + 1L |
---|
174 | jpj = count[1] |
---|
175 | |
---|
176 | IF (mskinq.ndims EQ 3 AND withrcd EQ -1) OR (mskinq.ndims EQ 4) THEN BEGIN |
---|
177 | IF n_elements(zmaxmesh) NE 0 THEN BEGIN |
---|
178 | IF zmaxmesh GE 0 THEN izmaxmesh = long(zmaxmesh[0]) ELSE izmaxmesh = count[2] - 1L + long(zmaxmesh[0]) |
---|
179 | ENDIF ELSE izmaxmesh = count[2] - 1L |
---|
180 | izmaxmesh = 0 > izmaxmesh < (count[2] - 1L) ; make sure ixmaxmesh is not too big |
---|
181 | jpkglo = count[2] |
---|
182 | offset[2] = izminmesh |
---|
183 | count[2] = izmaxmesh - izminmesh + 1L |
---|
184 | jpk = count[2] |
---|
185 | ENDIF |
---|
186 | |
---|
187 | ; read the variable for the first record |
---|
188 | ncdf_varget, cdfid, mskid, mask, count = count, offset = offset |
---|
189 | |
---|
190 | ; get add_offset, scale factor and missing value attributes |
---|
191 | ncdf_getatt, cdfid, mskid, add_offset = add, scale_factor = scl, missing_value = miss |
---|
192 | ; do we apply add_offset and scale factor ? |
---|
193 | IF keyword_set(addscl_before) THEN BEGIN |
---|
194 | IF scl NE 1 THEN mask = mask * scl |
---|
195 | IF add NE 0 THEN mask = mask + add |
---|
196 | ENDIF |
---|
197 | |
---|
198 | IF keyword_set(useasmask) THEN BEGIN |
---|
199 | |
---|
200 | IF n_elements(missing_value) NE 0 THEN miss = missing_value |
---|
201 | IF size(miss, /type) EQ 7 THEN miss = !values.f_nan |
---|
202 | |
---|
203 | IF NOT keyword_set(testop) THEN testop = 'NE' |
---|
204 | tmp = strsplit(testop, ' ', /extract) |
---|
205 | op = strupcase(tmp[0]) |
---|
206 | IF op EQ 'EQ' THEN BEGIN |
---|
207 | op = 'NE' |
---|
208 | invmask = 1b - keyword_set(invmask) |
---|
209 | ENDIF |
---|
210 | IF n_elements(tmp) EQ 1 THEN testval = miss ELSE testval = float(tmp[1]) |
---|
211 | IF finite(testval) EQ 0 THEN BEGIN |
---|
212 | IF op NE 'NE' THEN mask = report('NaN test value can be used only with EQ or NE operator') ELSE mask = finite(mask) |
---|
213 | ENDIF ELSE BEGIN |
---|
214 | CASE op OF |
---|
215 | 'GE':mask = mask GE testval |
---|
216 | 'GT':mask = mask GT testval |
---|
217 | 'LE':mask = mask LE testval |
---|
218 | 'LT':mask = mask LT testval |
---|
219 | 'NE':BEGIN |
---|
220 | ; we have to take care of the float accuracy |
---|
221 | CASE 1 OF |
---|
222 | testval GE 1.e6:mask = mask LT (testval - 10) |
---|
223 | testval LE -1.e6:mask = mask GT (testval + 10) |
---|
224 | abs(testval) LE 1.e-6:mask = abs(mask) GT 1.e-6 |
---|
225 | ELSE:mask = mask NE testval |
---|
226 | ENDCASE |
---|
227 | END |
---|
228 | ENDCASE |
---|
229 | ENDELSE |
---|
230 | |
---|
231 | ENDIF |
---|
232 | |
---|
233 | IF mask[0] NE -1 THEN BEGIN |
---|
234 | mask = byte(round(mask)) |
---|
235 | if keyword_set(invmask) then mask = 1b - mask |
---|
236 | ENDIF |
---|
237 | |
---|
238 | ENDIF ELSE mask = -1 |
---|
239 | |
---|
240 | IF size(fileid, /type) EQ 7 THEN ncdf_close, cdfid |
---|
241 | |
---|
242 | return, mask |
---|
243 | END |
---|