[2] | 1 | ;+ |
---|
| 2 | ; |
---|
[226] | 3 | ; @file_comments |
---|
[172] | 4 | ; compute the grid parameters (cm_4mesh) common |
---|
[103] | 5 | ; |
---|
[172] | 6 | ; domains sizes: |
---|
[192] | 7 | ; --------------- |
---|
[172] | 8 | ; jpi, jpj, jpk, jpiglo, jpjglo, jpkglo, jpidta, jpjdta, jpkdta |
---|
[118] | 9 | ; |
---|
[172] | 10 | ; domains positions regarding to the original grid: |
---|
[192] | 11 | ; -------------------------------------------------- |
---|
[172] | 12 | ; ixminmesh, ixmaxmesh, iyminmesh, iymaxmesh, izminmesh, izmaxmesh |
---|
| 13 | ; ixmindta, ixmaxdta, iymindta, iymaxdta, izmindta, izmaxdta |
---|
[124] | 14 | ; |
---|
[172] | 15 | ; horizontal parameters: |
---|
| 16 | ; ---------------------- |
---|
| 17 | ; glamt, glamf, gphit, gphit, e1t, e2t |
---|
[2] | 18 | ; |
---|
[172] | 19 | ; additional horizontal parameters if FULLCGRID keyword is defined: |
---|
| 20 | ; ----------------------------------------------------------------- |
---|
| 21 | ; glamu, glamv, gphiu, gphiv, e1u, e1v, e1f, e2u, e2v, e2f |
---|
| 22 | ; |
---|
| 23 | ; verticals parameters: |
---|
| 24 | ; --------------------- |
---|
| 25 | ; gdept, gdepw, e3t, e3w |
---|
| 26 | ; |
---|
| 27 | ; masks: |
---|
| 28 | ; ------ |
---|
[124] | 29 | ; tmask |
---|
[13] | 30 | ; |
---|
[172] | 31 | ; additional masks if FULLCGRID keyword is defined: |
---|
| 32 | ; ------------------------------------------------- |
---|
| 33 | ; umaskred, vmaskred, fmaskredx, fmaskredy |
---|
[124] | 34 | ; |
---|
[172] | 35 | ; triangles_list: |
---|
| 36 | ; --------------- |
---|
[118] | 37 | ; triangulation |
---|
[13] | 38 | ; |
---|
[172] | 39 | ; key_* variables: |
---|
| 40 | ; ---------------- |
---|
| 41 | ; key_onearth, key_periodic, key_shift, key_stride, key_partialstep, |
---|
| 42 | ; key_yreverse, key_zreverse, key_gridtype |
---|
| 43 | ; |
---|
| 44 | ; xxx related variables: |
---|
| 45 | ; ---------------------- |
---|
| 46 | ; ccmeshparameters, ccreadparameters |
---|
| 47 | ; |
---|
[226] | 48 | ; @categories |
---|
[157] | 49 | ; Grid |
---|
[13] | 50 | ; |
---|
[172] | 51 | ; @param startx {in}{optional}{type=scalar} |
---|
| 52 | ; x starting point, optional if [XY]AXIS keyword is used |
---|
[2] | 53 | ; |
---|
[172] | 54 | ; @param starty {in}{optional}{type=scalar} |
---|
| 55 | ; y starting point, optional if [XY]AXIS keyword is used |
---|
[13] | 56 | ; |
---|
[172] | 57 | ; @param stepxin {in}{optional}{type=scalar or vector} |
---|
| 58 | ; x direction step, optional if [XY]AXIS keyword is used, must be > 0 |
---|
| 59 | ; if stepxin is a vector nx is not used |
---|
[13] | 60 | ; |
---|
[172] | 61 | ; @param stepyin {in}{optional}{type=scalar or vector} |
---|
| 62 | ; y direction step, optional if [XY]AXIS keyword is used, |
---|
| 63 | ; could be > 0 (south to north) or < 0 (north to south) |
---|
| 64 | ; if stepyin is a vector ny is not used |
---|
[13] | 65 | ; |
---|
[172] | 66 | ; @param nxin {in}{optional}{type=scalar} |
---|
| 67 | ; number of points in x direction, |
---|
| 68 | ; optional if [XY]AXIS keyword is used or stepxin is a vector |
---|
[13] | 69 | ; |
---|
[226] | 70 | ; @param nyin {in}{optional}{type=scalar} |
---|
[172] | 71 | ; number of points in y direction, |
---|
| 72 | ; optional if [XY]AXIS keyword is used or stepyin is a vector |
---|
[13] | 73 | ; |
---|
[172] | 74 | ; @keyword FULLCGRID {default=0}{type=scalar: 0 or 1} |
---|
| 75 | ; Activate to specify that you want to compute all the C grid parameters: |
---|
| 76 | ; definition of glam[uv], gphi[uv], e1[uvf], e2[uvf], [uv]maskred and |
---|
| 77 | ; fmaskred[xy] will be add to the default computations |
---|
| 78 | ; |
---|
| 79 | ; @keyword GLAMBOUNDARY {default=those defined in the file}{type=2 elements vector} |
---|
| 80 | ; Longitude boundaries that should be used to visualize the data. |
---|
| 81 | ; lon2 > lon1 |
---|
| 82 | ; lon2 - lon1 le 360 |
---|
| 83 | ; By default, the common (cm_4mesh) variable key_shift will be automatically |
---|
| 84 | ; defined according to GLAMBOUNDARY. |
---|
| 85 | ; |
---|
| 86 | ; @keyword MASK {default=array of 1}{type=2D or 3D array} |
---|
| 87 | ; Specify the land(0)/sea(1) mask |
---|
| 88 | ; |
---|
| 89 | ; @keyword ONEARTH {default=1}{type=scalar: 0 or 1} |
---|
| 90 | ; Force the manual definition of data localization on the earth or not |
---|
| 91 | ; 0) if the data are not on the earth |
---|
| 92 | ; 1) if the data are on earth (in that case we can for example use |
---|
| 93 | ; the labels 'longitude', 'latitude' in plots). |
---|
| 94 | ; The resulting value will be stored in the common (cm_4mesh) variable key_onearth |
---|
[327] | 95 | ; ONEARTH=0 forces PERIODIC=0, SHIFT=0 and is cancelling GLAMBOUNDARY |
---|
[172] | 96 | ; |
---|
| 97 | ; @keyword PERIODIC {default=computed by using the first line of glamt}{type=scalar: 0 or 1} |
---|
| 98 | ; Force the manual definition of the grid zonal periodicity. |
---|
| 99 | ; The resulting value will be stored in the common (cm_4mesh) variable key_periodic |
---|
[327] | 100 | ; PERIODIC=0 forces SHIFT=0 |
---|
[172] | 101 | ; |
---|
| 102 | ; @keyword PLAIN {default=0}{type=scalar: 0 or 1} |
---|
[327] | 103 | ; Force YREVERSE=0, ZREVERSE=0, PERIODIC=0, SHIFT=0, STRIDE=[1, 1, 1] and |
---|
[252] | 104 | ; suppress the automatic redefinition of the domain in case of x periodicity overlap, |
---|
| 105 | ; y periodicity overlap (ORCA type only) and mask border to 0. |
---|
[69] | 106 | ; |
---|
[172] | 107 | ; @keyword SHIFT {default=computed according to glamboundary}{type=scalar} |
---|
[226] | 108 | ; Force the manual definition of the zonal shift that must be apply to the data. |
---|
[172] | 109 | ; The resulting value will be stored in the common (cm_4mesh) variable key_shift |
---|
[327] | 110 | ; Note that if key_periodic=0 then in any case key_shift=0. |
---|
[13] | 111 | ; |
---|
[163] | 112 | ; @keyword STRCALLING {type=string} |
---|
| 113 | ; a string containing the calling command used to |
---|
[231] | 114 | ; call computegrid (this is used by <pro>xxx</pro>) |
---|
[69] | 115 | ; |
---|
[172] | 116 | ; @keyword STRIDE {default=[1, 1, 1]}{type=3 elements vector} |
---|
| 117 | ; Specify the stride in x, y and z direction. The resulting |
---|
| 118 | ; value will be stored in the common (cm_4mesh) variable key_stride |
---|
[13] | 119 | ; |
---|
[172] | 120 | ; @keyword XAXIS {type=1D or 2D array} |
---|
| 121 | ; Specify longitudes in this case startx, stepx and nx are not used but |
---|
| 122 | ; could be necessary if the y axis is not defined with yaxis. It must be |
---|
| 123 | ; possible to sort the first line of xaxis in the increasing order by |
---|
| 124 | ; shifting its elements. |
---|
[13] | 125 | ; |
---|
[172] | 126 | ; @keyword YAXIS {type=1D or 2D array} |
---|
| 127 | ; Specify latitudes in this case starty, stepy and ny are not used but |
---|
| 128 | ; starty and stepy could be necessary if the x axis is not defined with |
---|
[238] | 129 | ; xaxis. It must be sorted in the increasing or decreasing order (along each column if 2d array). |
---|
[13] | 130 | ; |
---|
[124] | 131 | ; @keyword XYINDEX activate to specify that the horizontal grid should |
---|
[13] | 132 | ; be simply defined by using the index of the points |
---|
| 133 | ; (xaxis = findgen(nx) and yaxis = findgen(ny)) |
---|
| 134 | ; using this keyword forces key_onearth=0 |
---|
| 135 | ; |
---|
[172] | 136 | ; @keyword XMINMESH {default=0L}{type=scalar} |
---|
| 137 | ; Define common (cm_4mesh) variables ixminmesh used to define the localization |
---|
| 138 | ; of the first point of the grid along the x direction in a zoom of the original grid |
---|
[13] | 139 | ; |
---|
[172] | 140 | ; @keyword YMINMESH {default=0L}{type=scalar} |
---|
| 141 | ; Define common (cm_4mesh) variables iyminmesh used to define the localization |
---|
| 142 | ; of the first point of the grid along the y direction in a zoom of the original grid |
---|
[13] | 143 | ; |
---|
[172] | 144 | ; @keyword ZMINMESH {default=0L}{type=scalar} |
---|
| 145 | ; Define common (cm_4mesh) variables izminmesh used to define the localization |
---|
| 146 | ; of the first point of the grid along the z direction in a zoom of the original grid |
---|
| 147 | ; |
---|
| 148 | ; @keyword XMAXMESH {default=jpiglo-1}{type=scalar} |
---|
| 149 | ; Define common (cm_4mesh) variables ixmaxmesh used to define the localization |
---|
| 150 | ; of the last point of the grid along the x direction in a zoom of the original grid |
---|
[192] | 151 | ; Note that if XMAXMESH < 0 then ixmaxmesh is defined as ixmaxmesh = jpiglo -1 + xmaxmesh |
---|
[172] | 152 | ; |
---|
| 153 | ; @keyword YMAXMESH {default=jpjglo-1}{type=scalar} |
---|
| 154 | ; Define common (cm_4mesh) variables iymaxmesh used to define the localization |
---|
| 155 | ; of the last point of the grid along the y direction in a zoom of the original grid |
---|
[192] | 156 | ; Note that if YMAXMESH < 0 then iymaxmesh is defined as iymaxmesh = jpjglo -1 + ymaxmesh |
---|
[172] | 157 | ; |
---|
| 158 | ; @keyword ZMAXMESH {default=jpkglo-1}{type=scalar} |
---|
| 159 | ; Define common (cm_4mesh) variables izmaxmesh used to define the localization |
---|
| 160 | ; of the last point of the grid along the z direction in a zoom of the original grid |
---|
[192] | 161 | ; Note that if ZMAXMESH < 0 then izmaxmesh is defined as izmaxmesh = jpkglo -1 + maxmesh |
---|
[172] | 162 | ; |
---|
[118] | 163 | ; @keyword FBASE2TBASE |
---|
[172] | 164 | ; Activate when the model is a C grid based on a F point |
---|
| 165 | ; (with a F point at the bottom-left corner and a T point at the |
---|
[226] | 166 | ; upper-right corner). In this case, we ignore |
---|
[209] | 167 | ; - the first line of F and V points |
---|
| 168 | ; - the last line of T and U points |
---|
| 169 | ; - if the grid is not x-periodic, the first column of F and U points |
---|
| 170 | ; - if the grid is not x-periodic, the last column of T and V points. |
---|
| 171 | ; => we are back to a C grid based on T point as for OPA model. |
---|
| 172 | ; Note that in that case, key_gridtype = 'c_f' and not 'c' (-> used in read_ncdf) |
---|
| 173 | ; Note that activate FBASE2TBASE forces FULLCGRID=1 |
---|
[118] | 174 | ; |
---|
[491] | 175 | ; @keyword FFBASE2TBASE |
---|
| 176 | ; Activate when the model is a C grid based on a F point |
---|
| 177 | ; (with a F point at the bottom-left corner and again a F point at the |
---|
| 178 | ; upper-right corner). In this case, we ignore |
---|
| 179 | ; - the first line of F and V points |
---|
| 180 | ; - if the grid is not x-periodic, the first column of F and U points |
---|
| 181 | ; => we are back to a C grid based on T point as for OPA model. |
---|
| 182 | ; Note that in that case, key_gridtype = 'c_ff' and not 'c' (-> used in read_ncdf) |
---|
| 183 | ; Note that activate FFBASE2TBASE forces FULLCGRID=1 |
---|
| 184 | ; |
---|
[209] | 185 | ; @keyword UBASE2TBASE |
---|
| 186 | ; Activate when the model is a C grid based on a U point |
---|
| 187 | ; (with a U point at the bottom-left corner and a T point at the |
---|
[226] | 188 | ; upper-right corner). In this case, we ignore |
---|
[209] | 189 | ; - if the grid is not x-periodic, the first column of F and U points |
---|
| 190 | ; - if the grid is not x-periodic, the last column of T and V points. |
---|
| 191 | ; => we are back to a C grid based on T point as for OPA model. |
---|
| 192 | ; Note that in that case, key_gridtype = 'c_u' and not 'c' (-> used in read_ncdf) |
---|
| 193 | ; Note that activate UBASE2TBASE forces FULLCGRID=1 |
---|
| 194 | ; |
---|
| 195 | ; @keyword VBASE2TBASE |
---|
| 196 | ; Activate when the model is a C grid based on a V point |
---|
| 197 | ; (with a V point at the bottom-left corner and a T point at the |
---|
[226] | 198 | ; upper-right corner). In this case, we ignore |
---|
[209] | 199 | ; - the first line of F and V points |
---|
| 200 | ; - the last line of T and U points |
---|
| 201 | ; => we are back to a C grid based on T point as for OPA model. |
---|
| 202 | ; Note that in that case, key_gridtype = 'c_v' and not 'c' (-> used in read_ncdf) |
---|
| 203 | ; Note that activate VBASE2TBASE forces FULLCGRID=1 |
---|
| 204 | ; |
---|
[172] | 205 | ; @keyword ROMSH {type=2D array} |
---|
| 206 | ; This array is the final bathymetry at RHO-points. It is stored in the common |
---|
| 207 | ; variable (cm_4mesh) romszinfos.h |
---|
| 208 | ; Used when the model is a ROMS C-grid with one more point |
---|
| 209 | ; in longitude for T and V grid and one more point in latitude |
---|
[226] | 210 | ; for T and U grid. In this case, we ignore |
---|
[172] | 211 | ; - the last line of T and U points |
---|
| 212 | ; - the last column of T and V points. |
---|
| 213 | ; => we are back to a C grid based on T point as for OPA model. |
---|
| 214 | ; Note that activate ROMSH forces FULLCGRID=1 |
---|
[118] | 215 | ; |
---|
[172] | 216 | ; @keyword STRCALLING {type=scalar string} |
---|
| 217 | ; Used by xxx... |
---|
[13] | 218 | ; |
---|
[252] | 219 | ; @keyword YREVERSE {default=computed according to gphit[0, 1:jpj-1] LT gphit[0, 0:jpj-2]}{type=scalar} |
---|
| 220 | ; Force the manual definition of the y reverse that must be apply to the data. |
---|
| 221 | ; The resulting value will be stored in the common (cm_4mesh) variable key_yreverse |
---|
| 222 | ; |
---|
[172] | 223 | ; @keyword ZAXIS {type=1D} |
---|
[238] | 224 | ; Specify the vertical axis. Must be sorted in the increasing or decreasing order |
---|
[118] | 225 | ; |
---|
[498] | 226 | ; @keyword ZZAXIS {type=1D} |
---|
| 227 | ; same a zaxis (overwrite zaxis definition). Needed when |
---|
| 228 | ; computegrid is called by inincdf because inincdf as alreday a keyword with a |
---|
| 229 | ; name sarting by zaxis: zaxisname |
---|
| 230 | ; |
---|
[252] | 231 | ; @keyword ZREVERSE {default=computed according to gdept[0] GT gdept[1]}{type=scalar} |
---|
| 232 | ; Force the manual definition of the z reverse that must be apply to the data. |
---|
| 233 | ; The resulting value will be stored in the common (cm_4mesh) variable key_zreverse |
---|
| 234 | ; |
---|
[226] | 235 | ; @keyword _EXTRA |
---|
[237] | 236 | ; not used in the present case ... |
---|
[172] | 237 | ; |
---|
[238] | 238 | ; @uses |
---|
[371] | 239 | ; <pro>cm_4mesh</pro> |
---|
| 240 | ; <pro>cm_4data</pro> |
---|
| 241 | ; <pro>cm_4cal</pro> |
---|
[2] | 242 | ; |
---|
[252] | 243 | |
---|
[238] | 244 | ; @restrictions |
---|
[237] | 245 | ; if the grid has x/y periodicity overlap and/or if |
---|
[371] | 246 | ; the mask has 0 everywhere at the border (like a closed sea) and |
---|
| 247 | ; if (we did not activate /plain and xminmesh, xmaxmesh, yminmesh, |
---|
| 248 | ; ymaxmesh keywords are defined to their default values), we redefine |
---|
[421] | 249 | ; xminmesh, xmaxmesh, yminmesh, ymaxmesh in order to remove the |
---|
[371] | 250 | ; overlapping part and/or to open the domain (avoid it be forced |
---|
| 251 | ; to use cell_fill = 1). |
---|
[2] | 252 | ; |
---|
[237] | 253 | ; FUV points definition is not exact if the grid is irregular |
---|
[2] | 254 | ; |
---|
[226] | 255 | ; @history |
---|
| 256 | ; Sebastien Masson (smasson\@lodyc.jussieu.fr) |
---|
[2] | 257 | ; 2000-04-20 |
---|
[371] | 258 | ; Sept 2004, several bug fix to suit C grid type... |
---|
[493] | 259 | ; Aug 2005, rewrite almost everything... |
---|
[118] | 260 | ; |
---|
[226] | 261 | ; @version |
---|
| 262 | ; $Id$ |
---|
[118] | 263 | ; |
---|
[2] | 264 | ;- |
---|
[13] | 265 | PRO computegrid, startx, starty, stepxin, stepyin, nxin, nyin $ |
---|
[501] | 266 | , XAXIS = xaxis, YAXIS = yaxis, ZAXIS = zaxis, ZZAXIS = zzaxis $ |
---|
| 267 | , MASK = mask, GLAMBOUNDARY = glamboundary $ |
---|
| 268 | , XMINMESH = xminmesh, XMAXMESH = xmaxmesh $ |
---|
| 269 | , YMINMESH = yminmesh, YMAXMESH = ymaxmesh $ |
---|
| 270 | , ZMINMESH = zminmesh, ZMAXMESH = zmaxmesh $ |
---|
| 271 | , ONEARTH = onearth, PERIODIC = periodic $ |
---|
| 272 | , PLAIN = plain, SHIFT = shift, STRIDE = stride $ |
---|
| 273 | , YREVERSE = yreverse, ZREVERSE = zreverse $ |
---|
| 274 | , FULLCGRID = fullcgrid, XYINDEX = xyindex $ |
---|
| 275 | , UBASE2TBASE = ubase2tbase, VBASE2TBASE = vbase2tbase $ |
---|
| 276 | , FBASE2TBASE = fbase2tbase, FFBASE2TBASE = ffbase2tbase $ |
---|
| 277 | , STRCALLING = strcalling, ROMSH = romsh, _EXTRA = ex |
---|
[114] | 278 | ; |
---|
| 279 | compile_opt idl2, strictarrsubs |
---|
| 280 | ; |
---|
[13] | 281 | @cm_4mesh |
---|
| 282 | @cm_4data |
---|
| 283 | @cm_4cal |
---|
| 284 | IF NOT keyword_set(key_forgetold) THEN BEGIN |
---|
| 285 | @updatenew |
---|
| 286 | @updatekwd |
---|
| 287 | ENDIF |
---|
| 288 | ;--------------------------------------------------------- |
---|
[2] | 289 | ;------------------------------------------------------------ |
---|
[13] | 290 | time1 = systime(1) ; for key_performance |
---|
| 291 | ;------------------------------------------------------------ |
---|
[2] | 292 | ; |
---|
[13] | 293 | ;==================================================== |
---|
| 294 | ; Check input parameters |
---|
| 295 | ;==================================================== |
---|
[2] | 296 | ; |
---|
[13] | 297 | ; xaxis related parameters |
---|
[2] | 298 | ; |
---|
[124] | 299 | if n_elements(xaxis) NE 0 then BEGIN |
---|
[13] | 300 | CASE (size(xaxis))[0] OF |
---|
| 301 | 0:nx = 1L |
---|
| 302 | 1:nx = (size(xaxis))[1] |
---|
| 303 | 2:nx = (size(xaxis))[1] |
---|
| 304 | ENDCASE |
---|
[124] | 305 | ENDIF ELSE BEGIN |
---|
| 306 | IF n_elements(startx) EQ 0 THEN BEGIN |
---|
[13] | 307 | dummy = report('If xaxis is not given, startx must be defined') |
---|
| 308 | return |
---|
| 309 | ENDIF |
---|
| 310 | CASE n_elements(stepxin) OF |
---|
[124] | 311 | 0:BEGIN |
---|
[13] | 312 | dummy = report('If xaxis is not given, stepxin must be defined') |
---|
| 313 | return |
---|
| 314 | END |
---|
[124] | 315 | 1:BEGIN |
---|
| 316 | IF n_elements(nxin) EQ 0 THEN BEGIN |
---|
[13] | 317 | dummy = report('If xaxis is not given and stepxin has only one element, nx must be defined') |
---|
| 318 | return |
---|
| 319 | ENDIF ELSE nx = nxin |
---|
| 320 | END |
---|
[124] | 321 | ELSE:nx = n_elements(stepxin) |
---|
[13] | 322 | ENDCASE |
---|
[124] | 323 | ENDELSE |
---|
[2] | 324 | ; |
---|
[13] | 325 | ; yaxis related parameters |
---|
[2] | 326 | ; |
---|
[124] | 327 | if n_elements(yaxis) NE 0 then BEGIN |
---|
[13] | 328 | CASE (size(yaxis))[0] OF |
---|
| 329 | 0:ny = 1L |
---|
| 330 | 1:ny = (size(yaxis))[1] |
---|
| 331 | 2:ny = (size(yaxis))[2] |
---|
| 332 | ENDCASE |
---|
[124] | 333 | ENDIF ELSE BEGIN |
---|
| 334 | IF n_elements(starty) EQ 0 THEN BEGIN |
---|
[13] | 335 | dummy = report('If yaxis is not given, starty must be defined') |
---|
| 336 | return |
---|
| 337 | ENDIF |
---|
| 338 | CASE n_elements(stepyin) OF |
---|
[124] | 339 | 0:BEGIN |
---|
[13] | 340 | dummy = report('If yaxis is not given, stepyin must be defined') |
---|
| 341 | return |
---|
| 342 | END |
---|
[124] | 343 | 1:BEGIN |
---|
| 344 | IF n_elements(nyin) EQ 0 THEN BEGIN |
---|
[13] | 345 | dummy = report('If yaxis is not given and stepyin has only one element, ny must be defined') |
---|
| 346 | return |
---|
| 347 | ENDIF ELSE ny = nyin |
---|
| 348 | END |
---|
[124] | 349 | ELSE:ny = n_elements(stepyin) |
---|
[13] | 350 | ENDCASE |
---|
| 351 | ENDELSE |
---|
[2] | 352 | ; |
---|
[13] | 353 | ; zaxis related parameters |
---|
[2] | 354 | ; |
---|
[499] | 355 | if n_elements(zzaxis) NE 0 then zaxis = zzaxis |
---|
[124] | 356 | if n_elements(zaxis) NE 0 then BEGIN |
---|
[13] | 357 | CASE (size(zaxis))[0] OF |
---|
| 358 | 0:nz = 1L |
---|
| 359 | 1:nz = (size(zaxis))[1] |
---|
| 360 | ELSE:BEGIN |
---|
[236] | 361 | ras = report( 'not coded') |
---|
[13] | 362 | stop |
---|
| 363 | END |
---|
| 364 | ENDCASE |
---|
| 365 | ENDIF ELSE nz = 1L |
---|
[2] | 366 | ; |
---|
[13] | 367 | ;==================================================== |
---|
| 368 | ; Others automatic definitions... |
---|
| 369 | ;==================================================== |
---|
| 370 | ; |
---|
[399] | 371 | IF keyword_set(romsh) THEN fullcgrid = 1 |
---|
| 372 | CASE 1 OF |
---|
[491] | 373 | keyword_set(ffbase2tbase):key_gridtype = 'c_ff' |
---|
[399] | 374 | keyword_set(fbase2tbase):key_gridtype = 'c_f' |
---|
| 375 | keyword_set(ubase2tbase):key_gridtype = 'c_u' |
---|
| 376 | keyword_set(vbase2tbase):key_gridtype = 'c_v' |
---|
| 377 | else:key_gridtype = 'c' |
---|
| 378 | ENDCASE |
---|
[491] | 379 | IF strlen(key_gridtype) NE 1 THEN fullcgrid = 1 |
---|
[399] | 380 | ; |
---|
[421] | 381 | IF n_elements(xminmesh) NE 0 OR n_elements(xmaxmesh) NE 0 THEN BEGIN |
---|
| 382 | IF n_elements(xminmesh) EQ 0 THEN xminmesh = ixminmesh |
---|
| 383 | IF n_elements(xmaxmesh) EQ 0 THEN xmaxmesh = ixmaxmesh |
---|
[399] | 384 | IF nx EQ jpi AND xminmesh EQ ixminmesh AND xmaxmesh EQ ixmaxmesh THEN xalreadycut = 1 |
---|
[172] | 385 | ENDIF |
---|
[399] | 386 | IF keyword_set(xalreadycut) THEN BEGIN |
---|
[495] | 387 | xmin = 0 |
---|
| 388 | xmax = jpi - 1 |
---|
[399] | 389 | nxx = jpi |
---|
| 390 | ENDIF ELSE BEGIN |
---|
| 391 | jpiglo = long(nx) |
---|
| 392 | IF n_elements(xminmesh) NE 0 THEN ixminmesh = long(xminmesh[0]) ELSE ixminmesh = 0l |
---|
| 393 | IF n_elements(xmaxmesh) NE 0 THEN ixmaxmesh = long(xmaxmesh[0]) ELSE ixmaxmesh = jpiglo-1 |
---|
| 394 | IF ixmaxmesh LT 0 THEN ixmaxmesh = jpiglo -1 + ixmaxmesh |
---|
| 395 | ixmaxmesh = 0 > ixmaxmesh < (jpiglo-1) |
---|
| 396 | ixminmesh = 0 > ixminmesh < ixmaxmesh |
---|
| 397 | jpi = ixmaxmesh-ixminmesh+1 |
---|
[495] | 398 | xmin = ixminmesh |
---|
[399] | 399 | xmax = ixmaxmesh |
---|
| 400 | nxx = jpiglo |
---|
| 401 | ENDELSE |
---|
| 402 | |
---|
[421] | 403 | IF n_elements(yminmesh) NE 0 OR n_elements(ymaxmesh) NE 0 THEN BEGIN |
---|
| 404 | IF n_elements(yminmesh) EQ 0 THEN yminmesh = iyminmesh |
---|
| 405 | IF n_elements(ymaxmesh) EQ 0 THEN ymaxmesh = iymaxmesh |
---|
[399] | 406 | IF ny EQ jpj AND yminmesh EQ iyminmesh AND ymaxmesh EQ iymaxmesh THEN yalreadycut = 1 |
---|
| 407 | ENDIF |
---|
| 408 | IF keyword_set(yalreadycut) THEN BEGIN |
---|
[495] | 409 | ymin = 0 |
---|
| 410 | ymax = jpj - 1 |
---|
[399] | 411 | nyy = jpj |
---|
[495] | 412 | ENDIF ELSE BEGIN |
---|
[399] | 413 | jpjglo = long(ny) |
---|
| 414 | IF n_elements(yminmesh) NE 0 THEN iyminmesh = long(yminmesh[0]) ELSE iyminmesh = 0l |
---|
| 415 | IF n_elements(ymaxmesh) NE 0 THEN iymaxmesh = long(ymaxmesh[0]) ELSE iymaxmesh = jpjglo-1 |
---|
| 416 | IF key_gridtype EQ 'c_v' OR key_gridtype EQ 'c_f' THEN iymaxmesh = iymaxmesh-1 |
---|
| 417 | IF iymaxmesh LT 0 THEN iymaxmesh = jpjglo -1 + iymaxmesh |
---|
| 418 | iymaxmesh = 0 > iymaxmesh < (jpjglo-1) |
---|
| 419 | iyminmesh = 0 > iyminmesh < iymaxmesh |
---|
| 420 | jpj = iymaxmesh-iyminmesh+1 |
---|
[495] | 421 | ymin = iyminmesh |
---|
[399] | 422 | ymax = iymaxmesh |
---|
| 423 | nyy = jpjglo |
---|
| 424 | ENDELSE |
---|
| 425 | |
---|
| 426 | IF n_elements(zminmesh) NE 0 AND n_elements(zmaxmesh) NE 0 THEN BEGIN |
---|
| 427 | IF nz EQ jpk AND zminmesh EQ izminmesh AND zmaxmesh EQ izmaxmesh THEN zalreadycut = 1 |
---|
| 428 | ENDIF |
---|
| 429 | IF keyword_set(zalreadycut) THEN BEGIN |
---|
[495] | 430 | zmin = 0 |
---|
| 431 | zmax = jpk - 1 |
---|
[399] | 432 | nzz = jpk |
---|
[495] | 433 | ENDIF ELSE BEGIN |
---|
[501] | 434 | jpkglo = long(nz) |
---|
[399] | 435 | IF n_elements(zminmesh) NE 0 THEN izminmesh = long(zminmesh[0]) ELSE izminmesh = 0l |
---|
| 436 | IF n_elements(zmaxmesh) NE 0 THEN izmaxmesh = long(zmaxmesh[0]) ELSE izmaxmesh = jpkglo-1 |
---|
| 437 | IF izmaxmesh LT 0 THEN izmaxmesh = jpkglo -1 + izmaxmesh |
---|
| 438 | izmaxmesh = 0 > izmaxmesh < (jpkglo-1) |
---|
| 439 | izminmesh = 0 > izminmesh < izmaxmesh |
---|
| 440 | jpk = izmaxmesh-izminmesh+1 |
---|
[495] | 441 | zmin = izminmesh |
---|
[399] | 442 | zmax = izmaxmesh |
---|
| 443 | nzz = jpkglo |
---|
| 444 | ENDELSE |
---|
[13] | 445 | ; |
---|
[69] | 446 | ; impact of plain keyword: |
---|
| 447 | ; |
---|
| 448 | IF keyword_set(plain) THEN BEGIN |
---|
[252] | 449 | yreverse = 0 |
---|
| 450 | zreverse = 0 |
---|
[69] | 451 | periodic = 0 |
---|
| 452 | shift = 0 |
---|
| 453 | stride = [1, 1, 1] |
---|
| 454 | ENDIF |
---|
| 455 | ; |
---|
[13] | 456 | ; avoid basics errors... |
---|
| 457 | ; |
---|
| 458 | jpidta = jpiglo |
---|
| 459 | jpjdta = jpjglo |
---|
| 460 | jpkdta = jpkglo |
---|
| 461 | ixmindta = 0 |
---|
| 462 | ixmaxdta = jpidta-1 |
---|
| 463 | iymindta = 0 |
---|
| 464 | iymaxdta = jpjdta-1 |
---|
| 465 | izmindta = 0 |
---|
| 466 | izmaxdta = jpkdta-1 |
---|
| 467 | ; |
---|
| 468 | key_partialstep = 0 |
---|
| 469 | if n_elements(stride) eq 3 then key_stride = stride $ |
---|
| 470 | ELSE key_stride = [1, 1, 1] |
---|
| 471 | ; |
---|
| 472 | ; check xyindex and its consequences |
---|
| 473 | ; |
---|
| 474 | if keyword_set(xyindex) then onearth = 0 |
---|
| 475 | ; |
---|
| 476 | ; check onearth and its consequences |
---|
[124] | 477 | ; |
---|
[13] | 478 | IF n_elements(onearth) EQ 0 THEN key_onearth = 1b $ |
---|
[124] | 479 | ELSE key_onearth = keyword_set(onearth) |
---|
[13] | 480 | IF NOT key_onearth THEN BEGIN |
---|
| 481 | periodic = 0 |
---|
| 482 | shift = 0 |
---|
[69] | 483 | ENDIF |
---|
[2] | 484 | |
---|
[13] | 485 | r = 6371000. |
---|
| 486 | ; |
---|
| 487 | ;==================================================== |
---|
| 488 | ; X direction : glamt |
---|
| 489 | ;==================================================== |
---|
| 490 | ; |
---|
| 491 | ; def of glamt |
---|
| 492 | ; |
---|
[124] | 493 | if n_elements(xaxis) NE 0 then BEGIN |
---|
[13] | 494 | if keyword_set(xyindex) THEN glamt = findgen(jpiglo) ELSE glamt = xaxis |
---|
| 495 | ENDIF ELSE BEGIN |
---|
| 496 | if keyword_set(xyindex) THEN stepx = 1. ELSE stepx = stepxin |
---|
| 497 | CASE 1 OF |
---|
| 498 | n_elements(stepx):glamt = startx + findgen(jpiglo)*stepx |
---|
| 499 | size(stepx, /n_dimensions):glamt = startx + total(stepx, /cumulative) |
---|
[124] | 500 | ELSE:BEGIN |
---|
[13] | 501 | dummy = report('Wrong definition of stepx...') |
---|
| 502 | return |
---|
[124] | 503 | END |
---|
[13] | 504 | ENDCASE |
---|
| 505 | ENDELSE |
---|
| 506 | ; |
---|
| 507 | ; apply glamboundary |
---|
| 508 | ; |
---|
| 509 | IF keyword_set(glamboundary) AND key_onearth THEN BEGIN |
---|
| 510 | IF glamboundary[0] GE glamboundary[1] THEN stop |
---|
| 511 | IF glamboundary[1]-glamboundary[0] GT 360 THEN stop |
---|
| 512 | glamt = glamt MOD 360 |
---|
| 513 | smaller = where(glamt LT glamboundary[0]) |
---|
| 514 | if smaller[0] NE -1 then glamt[smaller] = glamt[smaller]+360. |
---|
| 515 | bigger = where(glamt GE glamboundary[1]) |
---|
| 516 | if bigger[0] NE -1 then glamt[bigger] = glamt[bigger]-360. |
---|
| 517 | ENDIF |
---|
| 518 | ; |
---|
| 519 | ; force glamt to have 2 dimensions |
---|
| 520 | ; |
---|
[399] | 521 | IF n_elements(glamt) EQ nxx*nyy THEN glamt = reform(glamt, nxx, nyy, /over) $ |
---|
| 522 | ELSE glamt = reform(glamt, nxx, /over) |
---|
| 523 | CASE size(glamt, /n_dimensions) OF |
---|
[495] | 524 | 1:BEGIN |
---|
[399] | 525 | IF n_elements(glamt) EQ 1 THEN glamt = replicate(glamt[0], jpi, jpj) $ |
---|
| 526 | ELSE glamt = glamt[xmin:xmax]#replicate(1, jpj) |
---|
| 527 | END |
---|
| 528 | 2:glamt = glamt[xmin:xmax, ymin:ymax] |
---|
[13] | 529 | ENDCASE |
---|
[121] | 530 | ; keep 2d array even with degenerated dimension |
---|
[13] | 531 | IF jpj EQ 1 THEN glamt = reform(glamt, jpi, jpj, /over) |
---|
| 532 | ; |
---|
| 533 | ;==================================================== |
---|
| 534 | ; Y direction : gphit |
---|
| 535 | ;==================================================== |
---|
| 536 | ; |
---|
| 537 | ; def of gphit |
---|
| 538 | ; |
---|
| 539 | if n_elements(yaxis) NE 0 THEN BEGIN |
---|
| 540 | if keyword_set(xyindex) THEN gphit = findgen(jpjglo) ELSE gphit = yaxis |
---|
| 541 | ENDIF ELSE BEGIN |
---|
| 542 | if keyword_set(xyindex) THEN stepy = 1. ELSE stepy = stepyin |
---|
| 543 | CASE 1 OF |
---|
| 544 | n_elements(stepy):gphit = starty + findgen(jpjglo)*stepy |
---|
| 545 | size(stepy, /n_dimensions):gphit = starty + total(stepy, /cumulative) |
---|
[124] | 546 | ELSE:BEGIN |
---|
[13] | 547 | dummy = report('Wrong definition of stepy...') |
---|
| 548 | return |
---|
[124] | 549 | END |
---|
[13] | 550 | ENDCASE |
---|
| 551 | ENDELSE |
---|
| 552 | ; |
---|
| 553 | ; force gphit to have 2 dimensions |
---|
| 554 | ; |
---|
[399] | 555 | IF n_elements(gphit) EQ nxx*nyy THEN gphit = reform(gphit, nxx, nyy, /over) $ |
---|
| 556 | ELSE gphit = reform(gphit, nyy, /over) |
---|
| 557 | CASE size(gphit, /n_dimensions) OF |
---|
[495] | 558 | 1:BEGIN |
---|
[399] | 559 | IF n_elements(gphit) EQ 1 THEN gphit = replicate(gphit[0], jpi, jpj) $ |
---|
| 560 | ELSE gphit = replicate(1, jpi)#gphit[ymin:ymax] |
---|
| 561 | END |
---|
| 562 | 2:gphit = gphit[xmin:xmax, ymin:ymax] |
---|
[13] | 563 | ENDCASE |
---|
[121] | 564 | ; keep 2d array even with degenerated dimension |
---|
[124] | 565 | IF jpj EQ 1 THEN gphit = reform(gphit, jpi, jpj, /over) |
---|
[13] | 566 | ; |
---|
| 567 | ;==================================================== |
---|
[69] | 568 | ; check y periodicity... Only according to ORCA grid |
---|
[124] | 569 | ;==================================================== |
---|
[297] | 570 | ; check the periodicity if iyminmesh and iymaxmesh have the default definitions... |
---|
[252] | 571 | IF NOT keyword_set(plain) AND key_onearth EQ 1 $ |
---|
[501] | 572 | AND iyminmesh EQ 0l AND iymaxmesh eq jpjglo-1 AND jpj GE 3 AND jpi GE 2 THEN BEGIN |
---|
[69] | 573 | |
---|
| 574 | CASE 1 OF |
---|
| 575 | ixminmesh EQ 0l AND ixmaxmesh eq jpiglo-1 $ |
---|
[501] | 576 | AND array_equal(gphit[1:*, jpj-1], reverse(gphit[1:*, jpj-3])) EQ 1:BEGIN |
---|
[69] | 577 | ; T pivot |
---|
| 578 | ymaxmesh = -1 |
---|
| 579 | recall = 1 |
---|
[124] | 580 | END |
---|
[69] | 581 | ixminmesh EQ 1l AND ixmaxmesh eq jpiglo-2 $ |
---|
[124] | 582 | AND array_equal(gphit[*, jpj-1], reverse(shift(gphit[*, jpj-3], -1))) EQ 1:BEGIN |
---|
[69] | 583 | ; T pivot |
---|
| 584 | ymaxmesh = -1 |
---|
| 585 | recall = 1 |
---|
[124] | 586 | END |
---|
[69] | 587 | ixminmesh EQ 0l AND ixmaxmesh eq jpiglo-1 $ |
---|
[501] | 588 | AND array_equal(gphit[*, jpj-1], reverse(gphit[*, jpj-2])) EQ 1:BEGIN |
---|
[69] | 589 | ; F pivot |
---|
| 590 | ymaxmesh = -1 |
---|
| 591 | recall = 1 |
---|
[124] | 592 | END |
---|
[69] | 593 | ixminmesh EQ 1l AND ixmaxmesh eq jpiglo-2 $ |
---|
[124] | 594 | AND array_equal(gphit[*, jpj-1], reverse(gphit[*, jpj-2])) EQ 1:BEGIN |
---|
[69] | 595 | ; F pivot |
---|
| 596 | ymaxmesh = -1 |
---|
| 597 | recall = 1 |
---|
[124] | 598 | END |
---|
[69] | 599 | ELSE: |
---|
| 600 | ENDCASE |
---|
[124] | 601 | ENDIF |
---|
[69] | 602 | ; |
---|
| 603 | ;==================================================== |
---|
| 604 | ; check x periodicity... |
---|
[124] | 605 | ;==================================================== |
---|
[501] | 606 | IF n_elements(periodic) NE 0 THEN forcenoperio = 1 - keyword_set(periodic) |
---|
[297] | 607 | ; check the periodicity if ixminmesh and ixmaxmesh have the default definitions... |
---|
[69] | 608 | IF NOT keyword_set(plain) AND NOT keyword_set(forcenoperio) AND key_onearth EQ 1 $ |
---|
[252] | 609 | AND ixminmesh EQ 0l AND ixmaxmesh eq jpiglo-1 AND jpi GE 3 THEN BEGIN |
---|
[69] | 610 | CASE 0 OF |
---|
| 611 | total((glamt[0, *] - glamt[jpi-2, *]) MOD 360) $ |
---|
[501] | 612 | + total((glamt[1, *] - glamt[jpi-1, *]) MOD 360):BEGIN |
---|
[69] | 613 | xminmesh = 1 |
---|
| 614 | xmaxmesh = -1 |
---|
| 615 | recall = 1 |
---|
[124] | 616 | END |
---|
| 617 | total((glamt[0, *] - glamt[jpi-2, *]) MOD 360):BEGIN |
---|
[69] | 618 | xminmesh = 1 |
---|
| 619 | recall = 1 |
---|
[124] | 620 | END |
---|
| 621 | total((glamt[1, *] - glamt[jpi-1, *]) MOD 360):BEGIN |
---|
[69] | 622 | xmaxmesh = -1 |
---|
| 623 | recall = 1 |
---|
| 624 | END |
---|
[488] | 625 | total((glamt[0, *] - glamt[jpi-1, *]) MOD 360):BEGIN |
---|
| 626 | xmaxmesh = -1 |
---|
| 627 | recall = 1 |
---|
| 628 | END |
---|
[69] | 629 | ELSE: |
---|
| 630 | ENDCASE |
---|
[124] | 631 | ENDIF |
---|
[69] | 632 | ;==================================================== |
---|
| 633 | ; recall computegrid if needed... |
---|
[124] | 634 | ;==================================================== |
---|
[69] | 635 | IF keyword_set(recall) THEN BEGIN |
---|
| 636 | computegrid, XAXIS = glamt, YAXIS = gphit, ZAXIS = zaxis $ |
---|
| 637 | , MASK = mask, GLAMBOUNDARY = glamboundary $ |
---|
| 638 | , XMINMESH = xminmesh, XMAXMESH = xmaxmesh $ |
---|
| 639 | , YMINMESH = yminmesh, YMAXMESH = ymaxmesh $ |
---|
| 640 | , ZMINMESH = zminmesh, ZMAXMESH = zmaxmesh $ |
---|
| 641 | , PERIODIC = periodic, SHIFT = shift, STRIDE = stride $ |
---|
| 642 | , FULLCGRID = fullcgrid, XYINDEX = xyindex $ |
---|
[209] | 643 | , STRCALLING = strcalling $ |
---|
[172] | 644 | , ROMSH = romsh, _extra = ex |
---|
[69] | 645 | return |
---|
| 646 | ENDIF |
---|
[209] | 647 | ; |
---|
[69] | 648 | ;==================================================== |
---|
[209] | 649 | ; def key_yreverse |
---|
| 650 | ;==================================================== |
---|
| 651 | ; |
---|
[495] | 652 | IF n_elements(yreverse) EQ 0 THEN BEGIN |
---|
[252] | 653 | IF jpj GT 1 THEN BEGIN |
---|
| 654 | IF total(gphit[0, 1:jpj-1] LT gphit[0, 0:jpj-2]) GT jpj/2 THEN key_yreverse = 1 ELSE key_yreverse = 0 |
---|
[209] | 655 | ENDIF ELSE key_yreverse = 0 |
---|
[252] | 656 | ENDIF ELSE key_yreverse = yreverse |
---|
| 657 | IF keyword_set(key_yreverse) THEN BEGIN |
---|
| 658 | gphit = reverse(gphit, 2) |
---|
| 659 | glamt = reverse(glamt, 2) |
---|
| 660 | ENDIF |
---|
[209] | 661 | ; |
---|
| 662 | ;==================================================== |
---|
[13] | 663 | ; def of key_shift |
---|
| 664 | ;==================================================== |
---|
| 665 | ; |
---|
[121] | 666 | ; definition of key_shift by shifting the array to have the min |
---|
[13] | 667 | ; values of glamt[*, 0] in glamt[0, 0] |
---|
| 668 | ; |
---|
| 669 | IF n_elements(shift) EQ 0 THEN BEGIN |
---|
| 670 | IF jpi GT 1 then BEGIN |
---|
| 671 | xtest = glamt[*, 0] |
---|
| 672 | key_shift = (where(xtest EQ min(xtest)))[0] |
---|
| 673 | IF key_shift NE 0 THEN key_shift = jpi - key_shift |
---|
| 674 | ENDIF ELSE key_shift = 0 |
---|
| 675 | ENDIF ELSE key_shift = shift |
---|
| 676 | ; |
---|
| 677 | ;==================================================== |
---|
| 678 | ; def of key_periodic |
---|
| 679 | ;==================================================== |
---|
| 680 | ; |
---|
| 681 | IF n_elements(periodic) EQ 0 THEN BEGIN |
---|
| 682 | IF jpi GT 1 THEN BEGIN |
---|
| 683 | xtest = shift(glamt[*, 0], key_shift) |
---|
| 684 | ; check that xtest is now sorted in the increasing order |
---|
| 685 | IF array_equal(sort(xtest), lindgen(jpi)) EQ 0 THEN BEGIN |
---|
[240] | 686 | ras = report(['WARNING: we cannot sort the xaxis with a simple shift...', $ |
---|
[501] | 687 | 'we force key_periodic = 0 and key_shift = 0', $ |
---|
| 688 | 'only horizontal plot may be ok...']) |
---|
[13] | 689 | key_periodic = 0 |
---|
| 690 | xnotsorted = 1 |
---|
[124] | 691 | ENDIF ELSE BEGIN |
---|
[13] | 692 | key_periodic = (xtest[jpi-1]+2*(xtest[jpi-1]-xtest[jpi-2])) $ |
---|
| 693 | GE (xtest[0]+360) |
---|
| 694 | ENDELSE |
---|
| 695 | ENDIF ELSE key_periodic = 0 |
---|
[124] | 696 | ENDIF ELSE key_periodic = keyword_set(periodic) |
---|
[13] | 697 | ; |
---|
| 698 | ; update key_shift |
---|
| 699 | ; |
---|
| 700 | key_shift = key_shift * (key_periodic EQ 1) |
---|
| 701 | ; |
---|
[209] | 702 | IF (key_gridtype EQ 'c_u' OR key_gridtype EQ 'c_f') AND NOT keyword_set(key_periodic) THEN BEGIN |
---|
[13] | 703 | ixmaxmesh = ixmaxmesh-1 |
---|
| 704 | jpi = jpi-1 |
---|
| 705 | ENDIF |
---|
| 706 | ; |
---|
| 707 | ;==================================================== |
---|
| 708 | ; apply key_shift |
---|
| 709 | ;==================================================== |
---|
| 710 | ; |
---|
[124] | 711 | if keyword_set(key_shift) then BEGIN |
---|
[13] | 712 | glamt = shift(glamt, key_shift, 0) |
---|
| 713 | gphit = shift(gphit, key_shift, 0) |
---|
[124] | 714 | IF jpj EQ 1 THEN BEGIN |
---|
[13] | 715 | glamt = reform(glamt, jpi, jpj, /over) |
---|
| 716 | gphit = reform(gphit, jpi, jpj, /over) |
---|
| 717 | ENDIF |
---|
| 718 | ENDIF |
---|
| 719 | ; |
---|
| 720 | ;==================================================== |
---|
| 721 | ; Are we using a "regular" grid (that can be described |
---|
| 722 | ; with x vector and y vector)? |
---|
| 723 | ;==================================================== |
---|
| 724 | ; |
---|
| 725 | ; to get faster, we first test the most basic cases before |
---|
| 726 | ; testing the full array. |
---|
| 727 | ; |
---|
| 728 | CASE 1 OF |
---|
| 729 | keyword_set(xyindex):key_irregular = 0b |
---|
| 730 | jpi EQ 1 OR jpj EQ 1:key_irregular = 0b |
---|
| 731 | n_elements(xaxis) EQ 0 AND n_elements(yaxis) EQ 0:key_irregular = 0b |
---|
| 732 | size(reform(xaxis), /n_dimensions) EQ 1 AND size(reform(xaxis), /n_dimensions) EQ 1:key_irregular = 0b |
---|
| 733 | n_elements(xaxis) EQ 0 AND size(reform(yaxis), /n_dimensions) EQ 1:key_irregular = 0b |
---|
| 734 | n_elements(yaxis) EQ 0 AND size(reform(xaxis), /n_dimensions) EQ 1:key_irregular = 0b |
---|
| 735 | array_equal(glamt[*, 0], glamt[*, jpj-1]) EQ 0:key_irregular = 1b |
---|
| 736 | array_equal(gphit[0, *], gphit[jpi-1, *]) EQ 0:key_irregular = 1b |
---|
| 737 | array_equal(glamt, glamt[*, 0]#replicate(1, jpj)) EQ 0:key_irregular = 1b |
---|
| 738 | array_equal(gphit, replicate(1, jpi)#(gphit[0, *])[*]) EQ 0:key_irregular = 1b |
---|
| 739 | ELSE:key_irregular = 0b |
---|
| 740 | ENDCASE |
---|
| 741 | ; |
---|
| 742 | ;==================================================== |
---|
| 743 | ; def of glamf: defined as the middle of T(i,j) T(i+1,j+1) |
---|
| 744 | ;==================================================== |
---|
| 745 | ; |
---|
[124] | 746 | IF jpi GT 1 THEN BEGIN |
---|
[13] | 747 | ; we must compute stepxf: x distance between T(i,j) T(i+1,j+1) |
---|
| 748 | CASE 1 OF |
---|
| 749 | n_elements(stepx):stepxf = stepx |
---|
| 750 | size(stepx, /n_dimensions):stepxf = stepx#replicate(1, jpj) |
---|
| 751 | ELSE:BEGIN |
---|
| 752 | if (keyword_set(key_onearth) AND keyword_set(xnotsorted)) $ |
---|
[501] | 753 | OR (keyword_set(key_periodic) AND key_irregular) then BEGIN |
---|
[13] | 754 | stepxf = (glamt + 720) MOD 360 |
---|
| 755 | IF jpj EQ 1 THEN stepxf = reform(stepxf, jpi, jpj, /over) |
---|
| 756 | stepxf = shift(stepxf, -1, -1) - stepxf |
---|
| 757 | stepxf = [ [[stepxf]], [[stepxf + 360]], [[stepxf - 360]] ] |
---|
| 758 | stepxf = min(abs(stepxf), dimension = 3) |
---|
| 759 | IF NOT keyword_set(key_periodic) THEN $ |
---|
[501] | 760 | stepxf[jpi-1, *] = stepxf[jpi-2, *] |
---|
[13] | 761 | ENDIF ELSE BEGIN |
---|
| 762 | stepxf = shift(glamt, -1, -1) - glamt |
---|
| 763 | IF keyword_set(key_periodic) THEN $ |
---|
[501] | 764 | stepxf[jpi-1, *] = 360 + stepxf[jpi-1, *] $ |
---|
[13] | 765 | ELSE stepxf[jpi-1, *] = stepxf[jpi-2, *] |
---|
| 766 | ENDELSE |
---|
[124] | 767 | IF jpj GT 1 THEN BEGIN |
---|
[13] | 768 | stepxf[*, jpj-1] = stepxf[*, jpj-2] |
---|
| 769 | stepxf[jpi-1, jpj-1] = stepxf[jpi-2, jpj-2] |
---|
| 770 | ENDIF |
---|
| 771 | END |
---|
| 772 | ENDCASE |
---|
| 773 | glamf = glamt + 0.5 * stepxf |
---|
| 774 | ENDIF ELSE glamf = glamt + 0.5 |
---|
| 775 | ; |
---|
[69] | 776 | IF keyword_set(key_periodic) AND (max(glamf)-min(glamt)) GE 360 THEN BEGIN |
---|
[124] | 777 | IF NOT keyword_set(glamboundary) THEN BEGIN |
---|
[69] | 778 | bigger = where(glamf GE min(glamt)+360) |
---|
| 779 | glamf[bigger] = glamf[bigger]-360. |
---|
| 780 | ENDIF ELSE glamf = glamboundary[0] > temporary(glamf) < glamboundary[1] |
---|
| 781 | ENDIF |
---|
| 782 | ; |
---|
[74] | 783 | IF jpj EQ 1 THEN glamf = reform(glamf, jpi, jpj, /over) |
---|
| 784 | ; |
---|
[13] | 785 | ;==================================================== |
---|
| 786 | ; def of gphif: defined as the middle of T(i,j) T(i+1,j+1) |
---|
| 787 | ;==================================================== |
---|
| 788 | ; |
---|
[124] | 789 | IF jpj GT 1 THEN BEGIN |
---|
[13] | 790 | ; we must compute stepyf: y distance between T(i,j) T(i+1,j+1) |
---|
| 791 | CASE 1 OF |
---|
| 792 | n_elements(stepy):stepyf = stepy |
---|
| 793 | size(stepy, /n_dimensions):stepyf = replicate(1, jpi)#stepy |
---|
| 794 | ELSE:BEGIN |
---|
| 795 | stepyf = shift(gphit, -1, -1) - gphit |
---|
| 796 | stepyf[*, jpj-1] = stepyf[*, jpj-2] |
---|
| 797 | IF jpi GT 1 THEN BEGIN |
---|
| 798 | if NOT keyword_set(key_periodic) THEN $ |
---|
[501] | 799 | stepyf[jpi-1, *] = stepyf[jpi-2, *] |
---|
[13] | 800 | stepyf[jpi-1, jpj-1] = stepyf[jpi-2, jpj-2] |
---|
[124] | 801 | ENDIF |
---|
| 802 | END |
---|
[13] | 803 | ENDCASE |
---|
| 804 | gphif = gphit + 0.5 * stepyf |
---|
| 805 | ENDIF ELSE gphif = gphit + 0.5 |
---|
| 806 | IF key_onearth THEN gphif = -90. > gphif < 90. |
---|
| 807 | ; |
---|
[74] | 808 | IF jpj EQ 1 THEN gphif = reform(gphif, jpi, jpj, /over) |
---|
| 809 | ; |
---|
[13] | 810 | ;==================================================== |
---|
| 811 | ; e1t: x distance between U(i-1,j) and U(i,j) |
---|
| 812 | ;==================================================== |
---|
| 813 | ; |
---|
| 814 | ; *-|-*---|---*---| |
---|
| 815 | ; |
---|
| 816 | IF jpi GT 1 THEN BEGIN |
---|
| 817 | IF n_elements(stepx) NE 1 THEN BEGIN |
---|
| 818 | IF keyword_set(irregular) THEN BEGIN |
---|
| 819 | ; we must compute stepxu: x distance between T(i,j) T(i+1,j) |
---|
| 820 | IF keyword_set(key_periodic) THEN BEGIN |
---|
| 821 | stepxu = (glamt + 720) MOD 360 |
---|
| 822 | stepxu = shift(stepxu, -1, 0) - stepxu |
---|
| 823 | stepxu = [ [[stepxu]], [[stepxu + 360]], [[stepxu - 360]] ] |
---|
| 824 | stepxu = min(abs(stepxu), dimension = 3) |
---|
[124] | 825 | ENDIF ELSE BEGIN |
---|
[13] | 826 | stepxu = shift(glamt, -1, 0) - glamt |
---|
| 827 | stepxu[jpi-1, *] = stepxf[jpi-2, *] |
---|
[124] | 828 | ENDELSE |
---|
[13] | 829 | ENDIF ELSE stepxu = stepxf |
---|
| 830 | IF jpj EQ 1 THEN stepxu = reform(stepxu, jpi, jpj, /over) |
---|
| 831 | e1t = 0.5*(stepxu+shift(stepxu, 1, 0)) |
---|
| 832 | IF NOT keyword_set(key_periodic) THEN $ |
---|
[501] | 833 | e1t[0, *] = e1t[1, *] |
---|
[13] | 834 | ENDIF ELSE e1t = replicate(stepx, jpi, jpj) |
---|
| 835 | ENDIF ELSE e1t = replicate(1b, jpi, jpj) |
---|
| 836 | ; |
---|
[74] | 837 | IF jpj EQ 1 THEN e1t = reform(e1t, jpi, jpj, /over) |
---|
| 838 | ; |
---|
[13] | 839 | ;==================================================== |
---|
| 840 | ; e2t: y distance between V(i,j-1) and V(i,j) |
---|
| 841 | ;==================================================== |
---|
| 842 | ; |
---|
| 843 | IF jpj GT 1 THEN BEGIN |
---|
| 844 | ; we must compute stepyv: y distance between T(i,j) T(i,j+1) |
---|
| 845 | IF n_elements(stepy) NE 1 THEN BEGIN |
---|
| 846 | IF keyword_set(key_irregular) THEN BEGIN |
---|
| 847 | stepyv = shift(gphit, 0, -1) - gphit |
---|
| 848 | stepyv[*, jpj-1] = stepyv[*, jpj-2] |
---|
| 849 | ENDIF ELSE stepyv = stepyf |
---|
| 850 | e2t = 0.5*(stepyv+shift(stepyv, 0, 1)) |
---|
| 851 | e2t[*, 0] = e2t[*, 1] |
---|
| 852 | ENDIF ELSE e2t = replicate(stepy, jpi, jpj) |
---|
| 853 | ENDIF ELSE e2t = replicate(1b, jpi, jpj) |
---|
| 854 | ; |
---|
| 855 | IF key_onearth THEN e2t = r * !pi/180. * temporary(e2t) |
---|
| 856 | ; |
---|
[74] | 857 | IF jpj EQ 1 THEN e2t = reform(e2t, jpi, jpj, /over) |
---|
| 858 | ; |
---|
[13] | 859 | ;==================================================== |
---|
| 860 | ; def of glamu: defined as the middle of T(i,j) T(i+1,j) |
---|
| 861 | ;==================================================== |
---|
| 862 | ; |
---|
[501] | 863 | IF keyword_set(irregular) THEN BEGIN |
---|
| 864 | glamu = glamt + 0.5 * stepxu |
---|
| 865 | IF keyword_set(glamboundary) AND key_onearth THEN $ |
---|
| 866 | glamu = glamboundary[0] > temporary(glamu) < glamboundary[1] |
---|
| 867 | ENDIF ELSE glamu = glamf |
---|
[13] | 868 | ; |
---|
[501] | 869 | IF jpj EQ 1 THEN glamu = reform(glamu, jpi, jpj, /over) |
---|
[74] | 870 | ; |
---|
[13] | 871 | ;==================================================== |
---|
| 872 | ; def of gphiu: defined as the middle of T(i,j) T(i+1,j) |
---|
| 873 | ;==================================================== |
---|
| 874 | ; |
---|
[501] | 875 | IF jpi GT 1 THEN BEGIN |
---|
| 876 | ; we must compute stepyu: y distance between T(i+1,j) T(i,j) |
---|
| 877 | IF keyword_set(key_irregular) THEN BEGIN |
---|
| 878 | stepyu = shift(gphit, -1, 0) - gphit |
---|
| 879 | IF NOT keyword_set(key_periodic) THEN $ |
---|
| 880 | stepyu[jpi-1, *] = stepyu[jpi-2, *] |
---|
| 881 | gphiu = gphit + 0.5 * stepyu |
---|
[13] | 882 | ENDIF ELSE gphiu = gphit |
---|
[501] | 883 | ENDIF ELSE gphiu = gphit |
---|
[13] | 884 | IF key_onearth THEN gphiu = -90. > gphiu < 90. |
---|
| 885 | ; |
---|
[74] | 886 | IF jpj EQ 1 THEN gphiu = reform(gphiu, jpi, jpj, /over) |
---|
| 887 | ; |
---|
[13] | 888 | ;==================================================== |
---|
| 889 | ; def of glamv: defined as the middle of T(i,j) T(i,j+1) |
---|
| 890 | ;==================================================== |
---|
| 891 | ; |
---|
[501] | 892 | IF jpj GT 1 THEN BEGIN |
---|
| 893 | ; we must compute stepxv: x distance between T(i,j) T(i,j+1) |
---|
| 894 | IF keyword_set(irregular) THEN BEGIN |
---|
| 895 | IF keyword_set(key_periodic) THEN BEGIN |
---|
| 896 | stepxv = (glamt + 720) MOD 360 |
---|
| 897 | stepxv = shift(stepxv, 0, -1) - stepxv |
---|
| 898 | stepxv = [ [[stepxv]], [[stepxv + 360]], [[stepxv - 360]] ] |
---|
| 899 | stepxv = min(abs(stepxv), dimension = 3) |
---|
| 900 | ENDIF ELSE stepxv = shift(glamt, 0, -1) - glamt |
---|
| 901 | stepxv[*, jpj-1] = stepxv[*, jpj-2] |
---|
| 902 | glamv = glamt + 0.5 * stepxv |
---|
| 903 | IF keyword_set(glamboundary) AND key_onearth THEN $ |
---|
| 904 | glamv = glamboundary[0] > temporary(glamv) < glamboundary[1] |
---|
[13] | 905 | ENDIF ELSE glamv = glamt |
---|
[501] | 906 | ENDIF ELSE glamv = glamt |
---|
[13] | 907 | ; |
---|
| 908 | ;==================================================== |
---|
| 909 | ; def of gphiv: defined as the middle of T(i,j) T(i,j+1) |
---|
| 910 | ;==================================================== |
---|
| 911 | ; |
---|
[501] | 912 | IF keyword_set(key_irregular) THEN $ |
---|
| 913 | gphiv = gphit + 0.5 * stepyv $ |
---|
| 914 | ELSE gphiv = gphif |
---|
| 915 | IF key_onearth THEN gphiv = -90. > gphiv < 90. |
---|
[13] | 916 | ; |
---|
[501] | 917 | IF jpj EQ 1 THEN gphiv = reform(gphiv, jpi, jpj, /over) |
---|
| 918 | |
---|
| 919 | ;~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ |
---|
| 920 | IF keyword_set(fullcgrid) THEN BEGIN |
---|
[74] | 921 | ; |
---|
[13] | 922 | ;==================================================== |
---|
| 923 | ; e1u: x distance between T(i,j) and T(i+1,j) |
---|
| 924 | ;==================================================== |
---|
| 925 | ; |
---|
| 926 | IF jpi GT 1 AND n_elements(stepx) NE 1 THEN $ |
---|
[501] | 927 | e1u = stepxu ELSE e1u = e1t |
---|
[13] | 928 | ; |
---|
| 929 | ;==================================================== |
---|
| 930 | ; e2u: y distance between F(i,j-1) and F(i,j) |
---|
| 931 | ;==================================================== |
---|
| 932 | ; |
---|
| 933 | IF keyword_set(key_irregular) THEN BEGIN |
---|
[124] | 934 | e2u = gphif - shift(gphif, 0, 1) |
---|
[13] | 935 | e2u[*, 0] = e2u[*, 1] |
---|
| 936 | IF key_onearth THEN e2u = r * !pi/180. * temporary(e2u) |
---|
| 937 | ENDIF ELSE e2u = e2t |
---|
| 938 | ; |
---|
[74] | 939 | IF jpj EQ 1 THEN e2u = reform(e2u, jpi, jpj, /over) |
---|
| 940 | ; |
---|
[13] | 941 | ;==================================================== |
---|
| 942 | ; e1v: x distance between F(i-1,j) and F(i,j) |
---|
| 943 | ;==================================================== |
---|
| 944 | ; |
---|
| 945 | IF keyword_set(irregular) THEN BEGIN |
---|
| 946 | IF keyword_set(key_periodic) THEN BEGIN |
---|
| 947 | e1v = (glamf + 720) MOD 360 |
---|
[124] | 948 | e1v = e1v - shift(e1v, 1, 0) |
---|
[13] | 949 | e1v = [ [[e1v]], [[e1v + 360]], [[e1v - 360]] ] |
---|
| 950 | e1v = min(abs(e1v), dimension = 3) |
---|
[124] | 951 | ENDIF ELSE BEGIN |
---|
[13] | 952 | e1v = glamf - shift(glamf, 1, 0) |
---|
| 953 | e1v[0, *] = stepxf[1, *] |
---|
[124] | 954 | ENDELSE |
---|
| 955 | ENDIF ELSE e1v = e1t |
---|
[13] | 956 | ; |
---|
[74] | 957 | IF jpj EQ 1 THEN e1v = reform(e1v, jpi, jpj, /over) |
---|
| 958 | ; |
---|
[13] | 959 | ;==================================================== |
---|
| 960 | ; e2v: y distance between T(i,j) and T(i+1,j) |
---|
| 961 | ;==================================================== |
---|
| 962 | ; |
---|
| 963 | IF jpj GT 1 and n_elements(stepy) NE 1 THEN BEGIN |
---|
| 964 | e2v = stepyv |
---|
| 965 | IF key_onearth THEN e2v = r * !pi/180. * temporary(e2v) |
---|
| 966 | ENDIF ELSE e2v = e2t |
---|
| 967 | ; |
---|
| 968 | ;==================================================== |
---|
| 969 | ; e1f: x distance between V(i,j) and V(i+1,j) |
---|
| 970 | ;==================================================== |
---|
| 971 | ; |
---|
| 972 | IF keyword_set(irregular) THEN BEGIN |
---|
| 973 | IF keyword_set(key_periodic) THEN BEGIN |
---|
| 974 | e1f = (glamv + 720) MOD 360 |
---|
| 975 | e1f = shift(e1f, -1, 0) - e1f |
---|
| 976 | e1f = [ [[e1f]], [[e1f + 360]], [[e1f - 360]] ] |
---|
| 977 | e1f = min(abs(e1f), dimension = 3) |
---|
[124] | 978 | ENDIF ELSE BEGIN |
---|
[13] | 979 | e1f = shift(glamv, -1, 0) - glamt |
---|
| 980 | e1f[jpi-1, *] = stepxf[jpi-2, *] |
---|
[124] | 981 | ENDELSE |
---|
[13] | 982 | ENDIF ELSE e1f = e1u |
---|
| 983 | ; |
---|
[74] | 984 | IF jpj EQ 1 THEN e1f = reform(e1f, jpi, jpj, /over) |
---|
| 985 | ; |
---|
[13] | 986 | ;==================================================== |
---|
| 987 | ; e2f: y distance between U(i,j) and U(i,j+1) |
---|
| 988 | ;==================================================== |
---|
| 989 | ; |
---|
| 990 | IF keyword_set(key_irregular) THEN BEGIN |
---|
| 991 | e2f = shift(gphiu, 0, -1) - gphiu |
---|
| 992 | e2f[*, jpj-1] = e2f[*, jpj-2] |
---|
| 993 | IF key_onearth THEN e2f = r * !pi/180. * temporary(e2f) |
---|
| 994 | ENDIF ELSE e2f = e2v |
---|
| 995 | ; |
---|
[74] | 996 | IF jpj EQ 1 THEN e2f = reform(e2f, jpi, jpj, /over) |
---|
| 997 | ; |
---|
[13] | 998 | ENDIF |
---|
| 999 | ;~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ |
---|
| 1000 | ; |
---|
| 1001 | ;==================================================== |
---|
| 1002 | ; e1[tuvf] from degree to meters |
---|
| 1003 | ;==================================================== |
---|
| 1004 | ; |
---|
| 1005 | IF keyword_set(key_onearth) THEN BEGIN |
---|
[501] | 1006 | IF keyword_set(key_irregular) THEN BEGIN |
---|
| 1007 | ; T |
---|
| 1008 | e1t = map_npoints(glamu, gphiu, shift(glamu, 1, 0), shift(gphiu, 1, 0), /two_by_two) |
---|
| 1009 | e1t = reform(e1t, jpi, jpj, /over) |
---|
| 1010 | IF NOT keyword_set(key_periodic) THEN e1t[0, *] = e1t[1, *] |
---|
| 1011 | e2t = map_npoints(glamv, gphiv, shift(glamv, 0, 1), shift(gphiv, 0, 1), /two_by_two) |
---|
| 1012 | e2t = reform(e2t, jpi, jpj, /over) |
---|
| 1013 | e2t[*, 0] = e2t[*, 1] |
---|
| 1014 | IF keyword_set(fullcgrid) THEN BEGIN |
---|
| 1015 | ; U |
---|
| 1016 | e1u = map_npoints(glamt, gphit, shift(glamt, -1, 0), shift(gphit, -1, 0), /two_by_two) |
---|
| 1017 | e1u = reform(e1u, jpi, jpj, /over) |
---|
| 1018 | IF NOT keyword_set(key_periodic) THEN e1u[jpi-1, *] = e1u[jpi-2, *] |
---|
| 1019 | e2u = map_npoints(glamf, gphif, shift(glamf, 0, 1), shift(gphif, 0, 1), /two_by_two) |
---|
| 1020 | e2u = reform(e2u, jpi, jpj, /over) |
---|
| 1021 | e2u[*, 0] = e2u[*, 1] |
---|
| 1022 | ; V |
---|
| 1023 | e1v = map_npoints(glamf, gphif, shift(glamf, 1, 0), shift(gphif, 1, 0), /two_by_two) |
---|
| 1024 | e1v = reform(e1v, jpi, jpj, /over) |
---|
| 1025 | IF NOT keyword_set(key_periodic) THEN e1v[0, *] = e1v[1, *] |
---|
| 1026 | e2v = map_npoints(glamt, gphit, shift(glamt, 0, -1), shift(gphit, 0, -1), /two_by_two) |
---|
| 1027 | e2v = reform(e2v, jpi, jpj, /over) |
---|
| 1028 | e2v[*, jpj-1] = e2v[*, jpj-2] |
---|
| 1029 | ; F |
---|
| 1030 | e1f = map_npoints(glamv, gphiv, shift(glamv, -1, 0), shift(gphiv, -1, 0), /two_by_two) |
---|
| 1031 | e1f = reform(e1f, jpi, jpj, /over) |
---|
| 1032 | IF NOT keyword_set(key_periodic) THEN e1f[jpi-1, *] = e1f[jpi-2, *] |
---|
| 1033 | e2f = map_npoints(glamu, gphiu, shift(glamu, 0, -1), shift(gphiu, 0, -1), /two_by_two) |
---|
| 1034 | e2f = reform(e2f, jpi, jpj, /over) |
---|
| 1035 | e2f[*, jpj-1] = e2f[*, jpj-2] |
---|
| 1036 | ENDIF |
---|
| 1037 | |
---|
| 1038 | ENDIF ELSE BEGIN |
---|
| 1039 | e1t = r * !pi/180. * temporary(e1t) * cos(!pi/180.*gphit) |
---|
| 1040 | IF keyword_set(fullcgrid) THEN BEGIN |
---|
| 1041 | e1u = r * !pi/180. * temporary(e1u) * cos(!pi/180.*gphiu) |
---|
| 1042 | e1v = r * !pi/180. * temporary(e1v) * cos(!pi/180.*gphiv) |
---|
| 1043 | e1f = r * !pi/180. * temporary(e1f) * cos(!pi/180.*gphif) |
---|
| 1044 | ENDIF |
---|
| 1045 | ENDELSE |
---|
[13] | 1046 | ENDIF |
---|
| 1047 | ; |
---|
[124] | 1048 | IF jpj EQ 1 THEN BEGIN |
---|
[74] | 1049 | e1t = reform(e1t, jpi, jpj, /over) |
---|
| 1050 | IF keyword_set(fullcgrid) THEN BEGIN |
---|
| 1051 | e1u = reform(e1u, jpi, jpj, /over) |
---|
| 1052 | e1v = reform(e1v, jpi, jpj, /over) |
---|
| 1053 | e1f = reform(e1f, jpi, jpj, /over) |
---|
[124] | 1054 | ENDIF |
---|
| 1055 | ENDIF |
---|
[74] | 1056 | ; |
---|
[13] | 1057 | ;==================================================== |
---|
| 1058 | ; if not fullcgrid: make sure we don't use glam[uv], gphi[uv], e[12][uvf] |
---|
| 1059 | ;==================================================== |
---|
| 1060 | ; |
---|
| 1061 | IF NOT keyword_set(fullcgrid) THEN BEGIN |
---|
| 1062 | glamu = !values.f_nan & glamv = !values.f_nan |
---|
| 1063 | gphiu = !values.f_nan & gphiv = !values.f_nan |
---|
| 1064 | e1u = !values.f_nan & e1v = !values.f_nan & e1f = !values.f_nan |
---|
| 1065 | e2u = !values.f_nan & e2v = !values.f_nan & e2f = !values.f_nan |
---|
[124] | 1066 | firstxu = !values.f_nan & lastxu = !values.f_nan & nxu = !values.f_nan |
---|
| 1067 | firstyu = !values.f_nan & lastyu = !values.f_nan & nyu = !values.f_nan |
---|
| 1068 | firstxv = !values.f_nan & lastxv = !values.f_nan & nxv = !values.f_nan |
---|
| 1069 | firstyv = !values.f_nan & lastyv = !values.f_nan & nyv = !values.f_nan |
---|
[13] | 1070 | ENDIF |
---|
| 1071 | ; |
---|
| 1072 | ;==================================================== |
---|
| 1073 | ; Z direction |
---|
| 1074 | ;==================================================== |
---|
| 1075 | ; |
---|
[124] | 1076 | ; z axis |
---|
[13] | 1077 | ; |
---|
[501] | 1078 | CASE n_elements(zaxis) OF |
---|
| 1079 | 0:BEGIN |
---|
| 1080 | gdept = 0. |
---|
| 1081 | key_zreverse = 0 |
---|
| 1082 | END |
---|
| 1083 | 1:BEGIN |
---|
| 1084 | gdept = zaxis |
---|
| 1085 | key_zreverse = 0 |
---|
| 1086 | END |
---|
| 1087 | ELSE:BEGIN |
---|
| 1088 | gdept = zaxis[zmin:zmax] |
---|
| 1089 | IF n_elements(zreverse) EQ 0 THEN BEGIN |
---|
| 1090 | IF jpk GT 1 THEN BEGIN |
---|
| 1091 | if gdept[0] GT gdept[1] then key_zreverse = 1 ELSE key_zreverse = 0 |
---|
| 1092 | ENDIF ELSE key_zreverse = 0 |
---|
| 1093 | ENDIF ELSE key_zreverse = zreverse |
---|
| 1094 | IF keyword_set(key_zreverse) THEN gdept = reverse(gdept) |
---|
| 1095 | END |
---|
| 1096 | ENDCASE |
---|
[13] | 1097 | ; |
---|
[501] | 1098 | if n_elements(gdept) GT 1 then BEGIN |
---|
| 1099 | stepz = shift(gdept, -1)-gdept |
---|
| 1100 | stepz[jpk-1] = stepz[jpk-2] |
---|
| 1101 | gdepw = 0. > (gdept-stepz/2.) |
---|
| 1102 | ENDIF ELSE BEGIN |
---|
| 1103 | stepz = 1. |
---|
| 1104 | gdepw = gdept |
---|
| 1105 | ENDELSE |
---|
| 1106 | IF keyword_set(romsh) THEN gdepw = gdept |
---|
[13] | 1107 | ; |
---|
| 1108 | ;==================================================== |
---|
| 1109 | ; e3[tw]: |
---|
| 1110 | ;==================================================== |
---|
| 1111 | ; |
---|
[501] | 1112 | e3t = stepz |
---|
| 1113 | IF n_elements(stepz) GT 1 THEN BEGIN |
---|
| 1114 | e3w = 0.5*(stepz+shift(stepz, 1)) |
---|
| 1115 | e3w[0] = 0.5*e3t[0] |
---|
| 1116 | ENDIF ELSE e3w = e3t |
---|
[13] | 1117 | ; |
---|
| 1118 | ;==================================================== |
---|
| 1119 | ; Mask |
---|
| 1120 | ;==================================================== |
---|
| 1121 | ; |
---|
[121] | 1122 | ; default mask eq 1 |
---|
[172] | 1123 | if NOT keyword_set(mask) then tmask = -1 ELSE tmask = mask |
---|
[13] | 1124 | ; |
---|
[172] | 1125 | if tmask[0] NE -1 then BEGIN |
---|
[399] | 1126 | tmask = byte(temporary(tmask)) |
---|
| 1127 | IF n_elements(mask) EQ nxx*nyy AND nzz GT 1 THEN BEGIN |
---|
| 1128 | tmask = tmask[*]#replicate(1b, nzz) |
---|
| 1129 | tmask = reform(tmask, nxx, nyy, nzz, /overwrite) |
---|
[226] | 1130 | ENDIF |
---|
[399] | 1131 | IF nxx EQ 1 OR nyy EQ 1 OR nzz EQ 1 THEN tmask = reform(tmask, nxx, nyy, nzz, /overwrite) |
---|
| 1132 | tmask = tmask[xmin:xmax, ymin:ymax, zmin:zmax] |
---|
| 1133 | IF jpi EQ 1 OR jpj EQ 1 OR jpk EQ 1 THEN tmask = reform(tmask, jpi, jpj, jpk, /over) |
---|
[13] | 1134 | if key_shift NE 0 then tmask = shift(tmask, key_shift, 0, 0) |
---|
| 1135 | ; because tmask = reverse(tmask, 2) is not working if the 3rd |
---|
| 1136 | ; dimension of tmask = 1, we call reform. |
---|
| 1137 | IF jpk EQ 1 THEN tmask = reform(tmask, /over) |
---|
| 1138 | IF key_yreverse EQ 1 THEN tmask = reverse(tmask, 2) |
---|
| 1139 | IF jpk EQ 1 THEN tmask = reform(tmask, jpi, jpj, jpk, /over) |
---|
| 1140 | IF key_zreverse EQ 1 THEN tmask = reverse(tmask, 3) |
---|
| 1141 | IF jpk EQ 1 THEN tmask = reform(tmask, jpi, jpj, jpk, /over) |
---|
| 1142 | IF keyword_set(fullcgrid) THEN BEGIN |
---|
| 1143 | IF keyword_set(key_periodic) THEN BEGIN |
---|
| 1144 | msk = tmask*shift(tmask, -1, 0, 0) |
---|
| 1145 | umaskred = msk[jpi-1, *, *] |
---|
| 1146 | ENDIF ELSE umaskred = tmask[jpi-1, *, *] |
---|
| 1147 | vmaskred = tmask[*, jpj-1, *] |
---|
[124] | 1148 | ENDIF |
---|
[315] | 1149 | fmaskredy = tmask[jpi-1, *, *] |
---|
| 1150 | fmaskredx = tmask[*, jpj-1, *] |
---|
[124] | 1151 | ENDIF ELSE BEGIN |
---|
[13] | 1152 | tmask = replicate(1b, jpi, jpj, jpk) |
---|
| 1153 | IF keyword_set(fullcgrid) THEN BEGIN |
---|
| 1154 | umaskred = replicate(1b, jpj, jpk) |
---|
| 1155 | vmaskred = replicate(1b, jpi, jpk) |
---|
[124] | 1156 | ENDIF |
---|
[315] | 1157 | fmaskredy = replicate(1b, jpj, jpk) |
---|
| 1158 | fmaskredx = replicate(1b, jpi, jpk) |
---|
[13] | 1159 | ENDELSE |
---|
| 1160 | ; |
---|
[69] | 1161 | IF jpi GT 2 AND jpj GT 2 AND NOT keyword_set(plain) $ |
---|
[435] | 1162 | ;; AND ixminmesh EQ 0l AND ixmaxmesh eq jpiglo-1 $ |
---|
| 1163 | ;; AND iyminmesh EQ 0l AND iymaxmesh eq jpjglo-1 $ |
---|
[69] | 1164 | AND total(tmask[*, 0, *]) EQ 0 AND total(tmask[*, jpj-1, *]) EQ 0 $ |
---|
| 1165 | AND total(tmask[0, *, *]) EQ 0 AND total(tmask[jpi-1, *, *]) EQ 0 THEN BEGIN |
---|
[495] | 1166 | IF NOT keyword_set(key_periodic) THEN BEGIN |
---|
[435] | 1167 | xminmesh = ixminmesh + 1 |
---|
| 1168 | xmaxmesh = ixmaxmesh - 1 |
---|
| 1169 | ENDIF |
---|
| 1170 | yminmesh = iyminmesh + 1 |
---|
| 1171 | ymaxmesh = iymaxmesh - 1 |
---|
| 1172 | ; come back to the original grid before calling computegrid with the |
---|
| 1173 | ; new parameters... |
---|
| 1174 | IF keyword_set(yreverse) THEN BEGIN |
---|
| 1175 | gphit = reverse(gphit, 2) |
---|
| 1176 | glamt = reverse(glamt, 2) |
---|
[495] | 1177 | tmask = reverse(tmask, 2) |
---|
[435] | 1178 | IF jpk EQ 1 THEN tmask = reform(tmask, jpi, jpj, jpk, /over) |
---|
| 1179 | ENDIF |
---|
| 1180 | IF keyword_set(zreverse) THEN BEGIN |
---|
| 1181 | zaxis = reverse(zaxis) |
---|
| 1182 | tmask = reverse(tmask, 3) |
---|
| 1183 | IF jpk EQ 1 THEN tmask = reform(tmask, jpi, jpj, jpk, /over) |
---|
| 1184 | ENDIF |
---|
| 1185 | IF key_shift NE 0 THEN BEGIN |
---|
| 1186 | glamt = shift(glamt, -key_shift, 0) |
---|
| 1187 | gphit = shift(gphit, -key_shift, 0) |
---|
| 1188 | tmask = shift(tmask, -key_shift, 0, 0) |
---|
| 1189 | IF jpk EQ 1 THEN tmask = reform(tmask, jpi, jpj, jpk, /over) |
---|
| 1190 | ENDIF |
---|
| 1191 | IF ixminmesh NE 0 THEN BEGIN |
---|
| 1192 | glamt = [fltarr(ixminmesh, jpj), glamt] |
---|
| 1193 | gphit = [fltarr(ixminmesh, jpj), gphit] |
---|
| 1194 | tmask = [fltarr(ixminmesh, jpj), tmask] |
---|
| 1195 | jpi = jpi+ixminmesh |
---|
| 1196 | IF jpk EQ 1 THEN tmask = reform(tmask, jpi, jpj, jpk, /over) |
---|
| 1197 | ENDIF |
---|
| 1198 | IF ixmaxmesh NE jpiglo-1 THEN BEGIN |
---|
| 1199 | glamt = [glamt, fltarr(jpiglo-1-ixmaxmesh, jpj)] |
---|
| 1200 | gphit = [gphit, fltarr(jpiglo-1-ixmaxmesh, jpj)] |
---|
| 1201 | tmask = [tmask, fltarr(jpiglo-1-ixmaxmesh, jpj)] |
---|
| 1202 | jpi = jpi+jpiglo-1-ixmaxmesh |
---|
| 1203 | IF jpk EQ 1 THEN tmask = reform(tmask, jpi, jpj, jpk, /over) |
---|
| 1204 | ENDIF |
---|
| 1205 | IF iyminmesh NE 0 THEN BEGIN |
---|
| 1206 | glamt = [[fltarr(jpi, iyminmesh)], [glamt]] |
---|
| 1207 | gphit = [[fltarr(jpi, iyminmesh)], [gphit]] |
---|
| 1208 | tmask = [[fltarr(jpi, iyminmesh)], [tmask]] |
---|
| 1209 | jpj = jpj+iyminmesh |
---|
| 1210 | IF jpk EQ 1 THEN tmask = reform(tmask, jpi, jpj, jpk, /over) |
---|
| 1211 | ENDIF |
---|
| 1212 | IF iymaxmesh NE jpjglo-1 THEN BEGIN |
---|
| 1213 | glamt = [[glamt], [fltarr(jpi, jpjglo-1-iymaxmesh)]] |
---|
| 1214 | gphit = [[gphit], [fltarr(jpi, jpjglo-1-iymaxmesh)]] |
---|
| 1215 | tmask = [[tmask], [fltarr(jpi, jpjglo-1-iymaxmesh)]] |
---|
| 1216 | jpj = jpj+jpjglo-1-iymaxmesh |
---|
| 1217 | IF jpk EQ 1 THEN tmask = reform(tmask, jpi, jpj, jpk, /over) |
---|
| 1218 | ENDIF |
---|
| 1219 | computegrid, XAXIS = glamt, YAXIS = gphit, ZAXIS = zaxis $ |
---|
| 1220 | , MASK = tmask, GLAMBOUNDARY = glamboundary $ |
---|
| 1221 | , XMINMESH = xminmesh, XMAXMESH = xmaxmesh $ |
---|
| 1222 | , YMINMESH = yminmesh, YMAXMESH = ymaxmesh $ |
---|
| 1223 | , ZMINMESH = zminmesh, ZMAXMESH = zmaxmesh $ |
---|
| 1224 | , ONEARTH = onearth, PERIODIC = key_periodic $ |
---|
| 1225 | , PLAIN = plain, SHIFT = key_shift, STRIDE = key_stride $ |
---|
| 1226 | , FULLCGRID = fullcgrid, XYINDEX = xyindex $ |
---|
| 1227 | , UBASE2TBASE = ubase2tbase, VBASE2TBASE = vbase2tbase $ |
---|
[491] | 1228 | , FBASE2TBASE = fbase2tbase, FFBASE2TBASE = ffbase2tbase, STRCALLING = strcalling $ |
---|
[435] | 1229 | , ROMSH = romsh, _extra = ex |
---|
| 1230 | return |
---|
[124] | 1231 | ENDIF |
---|
[69] | 1232 | ; |
---|
[13] | 1233 | IF NOT keyword_set(fullcgrid) THEN BEGIN |
---|
| 1234 | umaskred = !values.f_nan |
---|
| 1235 | vmaskred = !values.f_nan |
---|
| 1236 | ENDIF |
---|
| 1237 | ; |
---|
| 1238 | ;==================================================== |
---|
| 1239 | ; stride... |
---|
| 1240 | ;==================================================== |
---|
| 1241 | ; |
---|
| 1242 | IF total(key_stride) GT 3 THEN BEGIN |
---|
| 1243 | IF key_shift NE 0 THEN BEGIN |
---|
| 1244 | ; for explanation, see header of read_ncdf_varget.pro |
---|
| 1245 | jpiright = key_shift |
---|
| 1246 | jpileft = jpi - key_shift - ( (key_stride[0]-1)-((key_shift-1) MOD key_stride[0]) ) |
---|
| 1247 | jpi = ((jpiright-1)/key_stride[0]+1) + ((jpileft-1)/key_stride[0]+1) |
---|
| 1248 | ENDIF ELSE jpi = (jpi-1)/key_stride[0]+1 |
---|
| 1249 | jpj = (jpj-1)/key_stride[1]+1 |
---|
| 1250 | jpk = (jpk-1)/key_stride[2]+1 |
---|
| 1251 | ; |
---|
| 1252 | glamt = (temporary(glamt))[0:*:stride[0], 0:*:stride[1]] |
---|
| 1253 | gphit = (temporary(gphit))[0:*:stride[0], 0:*:stride[1]] |
---|
| 1254 | e1t = (temporary(e1t))[0:*:stride[0], 0:*:stride[1]] |
---|
| 1255 | e2t = (temporary(e2t))[0:*:stride[0], 0:*:stride[1]] |
---|
| 1256 | tmask = (temporary(tmask))[0:*:stride[0], 0:*:stride[1], 0:*:stride[2]] |
---|
| 1257 | gdept = gdept[0:*:stride[2]] |
---|
| 1258 | gdepw = gdepw[0:*:stride[2]] |
---|
| 1259 | e3t = e3t[0:*:stride[2]] |
---|
| 1260 | e3w = e3w[0:*:stride[2]] |
---|
| 1261 | ; we must recompute glamf and gphif... |
---|
[124] | 1262 | IF jpi GT 1 THEN BEGIN |
---|
[74] | 1263 | if (keyword_set(key_onearth) AND keyword_set(xnotsorted)) $ |
---|
[501] | 1264 | OR (keyword_set(key_periodic) AND key_irregular) then BEGIN |
---|
[74] | 1265 | stepxf = (glamt + 720) MOD 360 |
---|
| 1266 | stepxf = shift(stepxf, -1, -1) - stepxf |
---|
| 1267 | stepxf = [ [[stepxf]], [[stepxf + 360]], [[stepxf - 360]] ] |
---|
| 1268 | stepxf = min(abs(stepxf), dimension = 3) |
---|
| 1269 | IF NOT keyword_set(key_periodic) THEN $ |
---|
[501] | 1270 | stepxf[jpi-1, *] = stepxf[jpi-2, *] |
---|
[74] | 1271 | ENDIF ELSE BEGIN |
---|
| 1272 | stepxf = shift(glamt, -1, -1) - glamt |
---|
| 1273 | IF keyword_set(key_periodic) THEN $ |
---|
[501] | 1274 | stepxf[jpi-1, *] = 360 + stepxf[jpi-1, *] $ |
---|
| 1275 | ELSE stepxf[jpi-1, *] = stepxf[jpi-2, *] |
---|
[74] | 1276 | ENDELSE |
---|
[124] | 1277 | IF jpj GT 1 THEN BEGIN |
---|
[74] | 1278 | stepxf[*, jpj-1] = stepxf[*, jpj-2] |
---|
| 1279 | stepxf[jpi-1, jpj-1] = stepxf[jpi-2, jpj-2] |
---|
| 1280 | ENDIF |
---|
| 1281 | glamf = glamt + 0.5 * stepxf |
---|
| 1282 | ENDIF ELSE glamf = glamt + 0.5 |
---|
[252] | 1283 | IF keyword_set(key_periodic) AND (max(glamf)-min(glamt)) GE 360 THEN BEGIN |
---|
| 1284 | IF NOT keyword_set(glamboundary) THEN BEGIN |
---|
| 1285 | bigger = where(glamf GE min(glamt)+360) |
---|
| 1286 | glamf[bigger] = glamf[bigger]-360. |
---|
| 1287 | ENDIF ELSE glamf = glamboundary[0] > temporary(glamf) < glamboundary[1] |
---|
| 1288 | ENDIF |
---|
[124] | 1289 | IF jpj GT 1 THEN BEGIN |
---|
[13] | 1290 | ; we must compute stepyf: y distance between T(i,j) T(i+1,j+1) |
---|
[74] | 1291 | stepyf = shift(gphit, -1, -1) - gphit |
---|
| 1292 | stepyf[*, jpj-1] = stepyf[*, jpj-2] |
---|
| 1293 | IF jpi GT 1 THEN BEGIN |
---|
| 1294 | if NOT keyword_set(key_periodic) THEN $ |
---|
[501] | 1295 | stepyf[jpi-1, *] = stepyf[jpi-2, *] |
---|
[74] | 1296 | stepyf[jpi-1, jpj-1] = stepyf[jpi-2, jpj-2] |
---|
[124] | 1297 | ENDIF |
---|
[74] | 1298 | gphif = gphit + 0.5 * stepyf |
---|
| 1299 | ENDIF ELSE gphif = gphit + 0.5 |
---|
[13] | 1300 | ; |
---|
[124] | 1301 | IF jpj EQ 1 THEN BEGIN |
---|
[74] | 1302 | glamt = reform(glamt, jpi, jpj, /over) |
---|
| 1303 | gphit = reform(gphit, jpi, jpj, /over) |
---|
| 1304 | glamf = reform(glamf, jpi, jpj, /over) |
---|
| 1305 | gphif = reform(gphif, jpi, jpj, /over) |
---|
| 1306 | e1t = reform(e1t, jpi, jpj, /over) |
---|
| 1307 | e2t = reform(e2t, jpi, jpj, /over) |
---|
| 1308 | ENDIF |
---|
| 1309 | ; |
---|
[13] | 1310 | IF keyword_set(fullcgrid) THEN BEGIN |
---|
| 1311 | glamu = (temporary(glamu))[0:*:stride[0], 0:*:stride[1]] |
---|
| 1312 | gphiu = (temporary(gphiu))[0:*:stride[0], 0:*:stride[1]] |
---|
| 1313 | e1u = (temporary(e1u))[0:*:stride[0], 0:*:stride[1]] |
---|
| 1314 | e2u = (temporary(e2u))[0:*:stride[0], 0:*:stride[1]] |
---|
| 1315 | glamv = (temporary(glamv))[0:*:stride[0], 0:*:stride[1]] |
---|
| 1316 | gphiv = (temporary(gphiv))[0:*:stride[0], 0:*:stride[1]] |
---|
| 1317 | e1v = (temporary(e1v))[0:*:stride[0], 0:*:stride[1]] |
---|
| 1318 | e2v = (temporary(e2v))[0:*:stride[0], 0:*:stride[1]] |
---|
| 1319 | e1f = (temporary(e1f))[0:*:stride[0], 0:*:stride[1]] |
---|
| 1320 | e2f = (temporary(e2f))[0:*:stride[0], 0:*:stride[1]] |
---|
| 1321 | umaskred = (temporary(umaskred))[0, 0:*:stride[1], 0:*:stride[2]] |
---|
| 1322 | vmaskred = (temporary(vmaskred))[0:*:stride[0], 0, 0:*:stride[2]] |
---|
| 1323 | fmaskredy = (temporary(fmaskredy))[0, 0:*:stride[1], 0:*:stride[2]] |
---|
| 1324 | fmaskredx = (temporary(fmaskredx))[0:*:stride[0], 0, 0:*:stride[2]] |
---|
[124] | 1325 | IF jpj EQ 1 THEN BEGIN |
---|
[74] | 1326 | glamu = reform(glamu, jpi, jpj, /over) |
---|
| 1327 | gphiu = reform(gphiu, jpi, jpj, /over) |
---|
| 1328 | e1u = reform(e1u, jpi, jpj, /over) |
---|
| 1329 | e2u = reform(e2u, jpi, jpj, /over) |
---|
| 1330 | glamv = reform(glamv, jpi, jpj, /over) |
---|
| 1331 | gphiv = reform(gphiv, jpi, jpj, /over) |
---|
| 1332 | e1v = reform(e1v, jpi, jpj, /over) |
---|
| 1333 | e2v = reform(e2v, jpi, jpj, /over) |
---|
| 1334 | e1f = reform(e1f, jpi, jpj, /over) |
---|
| 1335 | e2f = reform(e2f, jpi, jpj, /over) |
---|
| 1336 | ENDIF |
---|
[13] | 1337 | ENDIF |
---|
| 1338 | ENDIF |
---|
| 1339 | ; |
---|
| 1340 | ;==================================================== |
---|
| 1341 | ; apply all the grid parameters |
---|
| 1342 | ;==================================================== |
---|
| 1343 | ; |
---|
| 1344 | @updateold |
---|
| 1345 | domdef |
---|
| 1346 | ; |
---|
| 1347 | ;==================================================== |
---|
| 1348 | ; Triangulation |
---|
| 1349 | ;==================================================== |
---|
| 1350 | ; |
---|
[69] | 1351 | IF total(tmask) EQ jpi*jpj*jpk $ |
---|
[501] | 1352 | AND NOT keyword_set(key_irregular) THEN triangles_list = -1 $ |
---|
[69] | 1353 | ELSE BEGIN |
---|
| 1354 | ; are we using ORCA2 ? |
---|
[287] | 1355 | IF jpiglo EQ 182 AND jpi EQ 180 AND jpjglo EQ 149 AND jpj EQ 148 THEN $ |
---|
[69] | 1356 | triangles_list = triangule() ELSE triangles_list = triangule(/keep_cont) |
---|
[124] | 1357 | ENDELSE |
---|
[13] | 1358 | ; |
---|
| 1359 | ;==================================================== |
---|
| 1360 | ; time axis (default definition) |
---|
| 1361 | ;==================================================== |
---|
| 1362 | ; |
---|
[124] | 1363 | IF n_elements(time) EQ 0 OR n_elements(jpt) EQ 0 THEN BEGIN |
---|
[69] | 1364 | jpt = 1 |
---|
| 1365 | time = 0 |
---|
| 1366 | ENDIF |
---|
[13] | 1367 | ; |
---|
| 1368 | IF NOT keyword_set(key_forgetold) THEN BEGIN |
---|
| 1369 | @updateold |
---|
| 1370 | ENDIF |
---|
[69] | 1371 | ;==================================================== |
---|
| 1372 | ; grid parameters used by xxx |
---|
| 1373 | ;==================================================== |
---|
| 1374 | ; |
---|
[124] | 1375 | IF NOT keyword_set(strcalling) THEN BEGIN |
---|
[69] | 1376 | IF n_elements(ccmeshparameters) EQ 0 THEN strcalling = 'computegrid' $ |
---|
| 1377 | ELSE strcalling = ccmeshparameters.filename |
---|
| 1378 | ENDIF |
---|
[74] | 1379 | IF n_elements(glamt) GE 2 THEN BEGIN |
---|
| 1380 | glaminfo = moment(glamt) |
---|
| 1381 | IF finite(glaminfo[2]) EQ 0 THEN glaminfo = glaminfo[0:1] |
---|
| 1382 | gphiinfo = moment(gphit) |
---|
| 1383 | IF finite(gphiinfo[2]) EQ 0 THEN gphiinfo = gphiinfo[0:1] |
---|
[124] | 1384 | ENDIF ELSE BEGIN |
---|
[74] | 1385 | glaminfo = glamt |
---|
| 1386 | gphiinfo = gphit |
---|
[124] | 1387 | ENDELSE |
---|
[172] | 1388 | IF keyword_set(romsh) THEN $ |
---|
[399] | 1389 | romszinfos = {h:romsh[xmin:xmax, ymin:ymax], zeta:replicate(0., jpi, jpj), theta_s:-1, theta_b:-1, hc:-1} $ |
---|
[172] | 1390 | ELSE romszinfos = {h:-1, zeta:-1, theta_s:-1, theta_b:-1, hc:-1} |
---|
| 1391 | |
---|
[69] | 1392 | ccmeshparameters = {filename:strcalling $ |
---|
[501] | 1393 | , glaminfo:float(string(glaminfo, format = '(E11.4)')) $ |
---|
| 1394 | , gphiinfo:float(string(gphiinfo, format = '(E11.4)')) $ |
---|
| 1395 | , jpiglo:jpiglo, jpjglo:jpjglo, jpkglo:jpkglo $ |
---|
| 1396 | , jpi:jpi, jpj:jpj, jpk:jpk $ |
---|
| 1397 | , ixminmesh:ixminmesh, ixmaxmesh:ixmaxmesh $ |
---|
| 1398 | , iyminmesh:iyminmesh, iymaxmesh:iymaxmesh $ |
---|
| 1399 | , izminmesh:izminmesh, izmaxmesh:izmaxmesh $ |
---|
| 1400 | , key_shift:key_shift, key_periodic:key_periodic $ |
---|
| 1401 | , key_stride:key_stride, key_gridtype:key_gridtype $ |
---|
| 1402 | , key_yreverse:key_yreverse, key_zreverse:key_zreverse $ |
---|
| 1403 | , key_partialstep:key_partialstep, key_onearth:key_onearth} |
---|
[69] | 1404 | |
---|
| 1405 | ccreadparameters = {funclec_name:'read_ncdf' $ |
---|
[501] | 1406 | , jpidta:jpidta, jpjdta:jpjdta, jpkdta:jpkdta $ |
---|
| 1407 | , ixmindta:ixmindta, ixmaxdta:ixmaxdta $ |
---|
| 1408 | , iymindta:iymindta, iymaxdta:iymaxdta $ |
---|
| 1409 | , izmindta:izmindta, izmaxdta:izmaxdta} |
---|
[2] | 1410 | ;------------------------------------------------------------ |
---|
[13] | 1411 | IF keyword_set(key_performance) EQ 1 THEN $ |
---|
[501] | 1412 | print, 'time computegrid', systime(1)-time1 |
---|
[13] | 1413 | ;------------------------------------------------------------ |
---|
| 1414 | return |
---|
[2] | 1415 | end |
---|