Changeset 473
- Timestamp:
- 02/14/12 16:25:51 (12 years ago)
- Location:
- trunk/SRC/Interpolation
- Files:
-
- 3 edited
-
compute_fromreg_bilinear_weigaddr.pro (modified) (1 diff)
-
extrapolate.pro (modified) (4 diffs)
-
file_interp.pro (modified) (7 diffs)
Legend:
- Unmodified
- Added
- Removed
-
trunk/SRC/Interpolation/compute_fromreg_bilinear_weigaddr.pro
r372 r473 8 8 ; Interpolation 9 9 ; 10 ; @param alonin{in}{required}{type= 2d array}10 ; @param alonin{in}{required}{type=1d array} 11 11 ; longitude of the input data 12 12 ; 13 ; @param alatin {in}{required}{type= 2d array}13 ; @param alatin {in}{required}{type=1d array} 14 14 ; latitude of the input data 15 15 ; -
trunk/SRC/Interpolation/extrapolate.pro
r384 r473 64 64 ;- 65 65 FUNCTION extrapolate, zinput, maskinput, nb_iteration, FILLXDIR = fillxdir, FILLYDIR = fillydir $ 66 , X_PERIODIC = x_periodic, MINVAL=minval, MAXVAL=maxval, GE0=ge0, _EXTRA=ex 66 , X_PERIODIC = x_periodic, MINVAL = minval, MAXVAL = maxval, GE0 = ge0 $ 67 , SMWIN = smwin, NSMOOTH = nsmooth, _EXTRA = ex 67 68 ; 68 69 compile_opt idl2, strictarrsubs … … 104 105 ztmp[temporary(nan)] = 1.e20 105 106 msk = temporary(msk) * temporary(finztmp) 106 ENDIF ELSE finztmp = -1 ; free memory107 ENDIF ELSE finztmp = -1 ; free memory 107 108 z = temporary(ztmp) 108 109 nx2 = nx+2 … … 186 187 keyword_set(fillydir):zcoast = z[nx2+coast]+z[-nx2+coast] 187 188 ELSE:zcoast = z[1+coast]+z[-1+coast]+z[nx2+coast]+z[-nx2+coast] $ 188 +1./sqrt(2)*(z[nx2+1+coast]+z[nx2-1+coast] $189 +z[-nx2+1+coast]+z[-nx2-1+coast])189 +1./sqrt(2)*(z[nx2+1+coast]+z[nx2-1+coast] $ 190 +z[-nx2+1+coast]+z[-nx2-1+coast]) 190 191 ENDCASE 191 192 ; … … 224 225 ; 225 226 ENDWHILE 227 ; 228 z = z[1:nx, 1:ny] 229 ; 230 ;--------------------------------------------------------------- 231 ; smooth the filled values 232 ;--------------------------------------------------------------- 233 IF n_elements(smwin) EQ 0 THEN smwin = 0 234 IF smwin EQ 0 OR smwin EQ 1 THEN return, z 235 IF n_elements(nsmooth) EQ 0 THEN nsmooth = 25 236 IF nsmooth EQ 0 THEN return, z 237 238 ; add extra bands to avoid edge errors... 239 240 one = replicate(1, smwin) 241 IF keyword_set(x_periodic) THEN BEGIN 242 ; nx-smwin nx-1 0 smwin-1 nx-smwin nx-1 0 smwin-1 243 ; |----------| |----------|------------------------------|----------| |----------| 244 ; 0 smwin-1 smwin 2*smwin-1 nx nx+smwin-1 nx+smwin nx+2*smwin-1 245 z = [ z[nx-smwin:nx-1, *], z, z[0:smwin-1, *] ] 246 msk = [ maskinput[nx-smwin:nx-1, *], maskinput, maskinput[0:smwin-1, *] ] 247 ;; IF array_equal(z[0:smwin-1, *], z[nx:nx+smwin-1, *]) NE 1 THEN stop 248 ;; IF array_equal(z[nx+smwin:nx+2*smwin-1, *], z[smwin:2*smwin-1, *]) NE 1 THEN stop 249 ENDIF ELSE BEGIN 250 z = [ one#(z[0, *])[*], z, one#(z[nx-1, *])[*] ] 251 msk = [ one#(maskinput[0, *])[*], maskinput, one#(maskinput[nx-1, *])[*] ] 252 ENDELSE 253 z = [ [ z[*, 0]#one], [ z], [ z[*, ny-1]#one] ] 254 msk = [ [msk[*, 0]#one], [msk], [msk[*, ny-1]#one] ] 255 256 zorg = z 257 mskm1 = 1-msk 258 FOR i = 1, nsmooth DO BEGIN 259 z = smooth( temporary(z)*mskm1 + zorg*msk, smwin ) 260 IF keyword_set(x_periodic) THEN BEGIN 261 z[0:smwin-1, *] = z[nx:nx+smwin-1, *] 262 z[nx+smwin:nx+2*smwin-1, *] = z[smwin:2*smwin-1, *] 263 ENDIF ELSE BEGIN 264 z[0:smwin-1, *] = one#(z[smwin, *])[*] 265 z[nx+smwin:nx+2*smwin-1, *] = one#(z[nx+smwin-1, *])[*] 266 ENDELSE 267 z[*, 0:smwin-1] = z[*, smwin]#one 268 z[*, ny+smwin:ny+2*smwin-1] = z[*, ny+smwin-1]#one 269 ENDFOR 270 z = temporary(z)*mskm1 + zorg*msk 271 z = (temporary(z))[smwin:nx+smwin-1, smwin:ny+smwin-1] 272 ; 226 273 ;--------------------------------------------------------------- 227 274 ; we return the original size of the array 228 275 ;--------------------------------------------------------------- 229 276 ; 230 return, z [1:nx, 1:ny]277 return, z 231 278 END -
trunk/SRC/Interpolation/file_interp.pro
r464 r473 29 29 , DIVX = divx, DIVY = divy, OUTMASK_IND = outmask_ind $ 30 30 , GETHAN = gethan, LETHAN = lethan, SET_OUTMSKVAL = set_outmskval $ 31 , DATA_SMWIN = data_smwin, DATA_NSMOOTH = data_nsmooth $ 32 , EXTRAP_SMWIN = extrap_smwin, EXTRAP_NSMOOTH = extrap_nsmooth $ 31 33 , _EXTRA = ex 32 34 ; … … 34 36 ; 35 37 ; for byte, short and long, convert to double before extrapolation and interpolation 36 intype = size(data, /type) 38 sz = size(data) 39 nx = sz[1] 40 ny = sz[2] 41 intype = sz[3] 37 42 if intype LE 3 THEN data = double(temporary(data)) 38 43 ; 39 44 ; take care of NaN values 40 45 nanmask = finite(data) 41 totnanmask = total(nanmask )46 totnanmask = total(nanmask, /integer) 42 47 IF totnanmask EQ 0 THEN return, !values.f_nan 43 48 IF totnanmask NE n_elements(nanmask) THEN BEGIN … … 57 62 IF mask[0] NE -1 THEN mask = temporary(missmask) * mask ELSE mask = temporary(missmask) 58 63 ENDIF 59 60 64 ; extrapolation 61 65 IF keyword_set(smooth) THEN data = extrapsmooth(temporary(data), mask, /x_periodic, _extra = ex) $ 62 ELSE data = extrapolate(temporary(data), mask, /x_periodic, _extra = ex) 66 ELSE data = extrapolate(temporary(data), mask, /x_periodic, SMWIN = extrap_smwin, NSMOOTH = extrap_nsmooth, _extra = ex) 67 ; smooth 68 IF keyword_set(data_smwin) THEN BEGIN 69 ; add extra bands to avoid edge errors... 70 one = replicate(1, data_smwin) 71 ; 72 ; nx-smwin nx-1 0 smwin-1 nx-smwin nx-1 0 smwin-1 73 ; |----------| |----------|------------------------------|----------| |----------| 74 ; 0 smwin-1 smwin 2*smwin-1 nx nx+smwin-1 nx+smwin nx+2*smwin-1 75 ; define boundaries (extra colomns and rows) 76 data = [ data[nx-data_smwin:nx-1, *], data, data[0:data_smwin-1, *] ] 77 data = [ [data[*, 0]#one], [data], [data[*, ny-1]#one] ] 78 FOR i = 1, data_nsmooth DO BEGIN 79 data = smooth( temporary(data), data_smwin ) 80 ; update boundaries 81 data[0:data_smwin-1, *] = data[nx:nx+data_smwin-1, *] 82 data[nx+data_smwin:nx+2*data_smwin-1, *] = data[data_smwin:2*data_smwin-1, *] 83 data[*, 0:data_smwin-1] = data[*, data_smwin]#one 84 data[*, ny+data_smwin:ny+2*data_smwin-1] = data[*, ny+data_smwin-1]#one 85 ENDFOR 86 ; keep only the interior domain 87 data = (temporary(data))[data_smwin:nx+data_smwin-1, data_smwin:ny+data_smwin-1] 88 ENDIF 63 89 ; interpolation 64 90 IF method EQ 'boxmean' THEN BEGIN … … 314 340 , OUTXAXISNAME = outxaxisname, OUTYAXISNAME = outyaxisname $ 315 341 , GETHAN = gethan, LETHAN = lethan, SET_OUTMSKVAL = set_outmskval $ 342 , DATA_SMWIN = data_smwin, DATA_NSMOOTH = data_nsmooth $ 343 , EXTRAP_SMWIN = extrap_smwin, EXTRAP_NSMOOTH = extrap_nsmooth $ 316 344 , _EXTRA = ex 317 345 ; … … 555 583 , SET_OUTMSKVAL = outmiss[i] $ 556 584 , DIVX = divx, DIVY = divy, OUTMASK_IND = outmask_ind $ 557 , GETHAN = gethan, LETHAN = lethan, _extra = ex) 585 , GETHAN = gethan, LETHAN = lethan $ 586 , DATA_SMWIN = data_smwin, DATA_NSMOOTH = data_nsmooth $ 587 , EXTRAP_SMWIN = extrap_smwin, EXTRAP_NSMOOTH = extrap_nsmooth $ 588 , _extra = ex) 558 589 IF interp AND n_elements(data) EQ 1 THEN data = replicate(data, jpiout, jpjout) 559 590 ncdf_varput, outid, outvarid[i], temporary(data) … … 575 606 , SET_OUTMSKVAL = outmiss[i] $ 576 607 , DIVX = divx, DIVY = divy, OUTMASK_IND = outmask_ind $ 577 , GETHAN = gethan, LETHAN = lethan, _extra = ex) 608 , GETHAN = gethan, LETHAN = lethan $ 609 , DATA_SMWIN = data_smwin, DATA_NSMOOTH = data_nsmooth $ 610 , EXTRAP_SMWIN = extrap_smwin, EXTRAP_NSMOOTH = extrap_nsmooth $ 611 , _extra = ex) 578 612 IF interp AND n_elements(data) EQ 1 THEN data = replicate(data, jpiout, jpjout) 579 613 ncdf_varput, outid, outvarid[i], temporary(data), offset = off, count = outcnt … … 596 630 , SET_OUTMSKVAL = outmiss[i] $ 597 631 , DIVX = divx, DIVY = divy, OUTMASK_IND = outmask_ind $ 598 , GETHAN = gethan, LETHAN = lethan, _extra = ex) 632 , GETHAN = gethan, LETHAN = lethan $ 633 , DATA_SMWIN = data_smwin, DATA_NSMOOTH = data_nsmooth $ 634 , EXTRAP_SMWIN = extrap_smwin, EXTRAP_NSMOOTH = extrap_nsmooth $ 635 , _extra = ex) 599 636 IF interp AND n_elements(data) EQ 1 THEN data = replicate(data, jpiout, jpjout) 600 637 ncdf_varput, outid, outvarid[i], temporary(data), offset = off, count = outcnt
Note: See TracChangeset
for help on using the changeset viewer.
