Ignore:
Timestamp:
07/06/06 16:10:25 (18 years ago)
Author:
pinsard
Message:

improvements of Interpolation/*.pro header

File:
1 edited

Legend:

Unmodified
Added
Removed

  • trunk/SRC/Interpolation/compute_fromreg_imoms3_weigaddr.pro

    r118 r125  
    11;+ 
    22; 
    3 ; @file_comments compute the weight and address neede to interpolate data from a 
    4 ;          "regular grid" to any grid using the imoms3 method 
    5 ;    
     3; @file_comments  
     4; compute the weight and address neede to interpolate data from a 
     5; "regular grid" to any grid using the imoms3 method 
     6; 
    67; @categories interpolation 
    78; 
    8 ; @param alonin {in}{required} longitude of the input data  
    9 ; @param alatin  {in}{required} latitude of the input data  
    10 ; @param olonin {in}{required} longitude of the output data  
    11 ; @param olat {in}{required} latitude of the output data  
    12 ; 
    13 ; @keyword /NONORTHERNLINE  
    14 ; @keyword /NOSOUTHERNLINE  
    15 ; activate if you don't whant to take into account the northen/southern line  
     9; @param alonin {in}{required} longitude of the input data 
     10; @param alatin {in}{required} latitude of the input data 
     11; @param olonin {in}{required} longitude of the output data 
     12; @param olat {in}{required} latitude of the output data 
     13; 
     14; @keyword NONORTHERNLINE 
     15; @keyword NOSOUTHERNLINE 
     16; activate if you don't want to take into account the northen/southern line 
    1617; of the input data when perfoming the interpolation. 
    1718; 
    … …  
    2425; 
    2526; @restrictions 
    26 ;  -  the input grid must be a "regular/rectangular grid", defined as a grid for 
    27 ;     which each lontitudes lines have the same latitude and each latitudes columns 
    28 ;     have the same longitude. 
     27;  - the input grid must be a "regular/rectangular grid", defined as a grid for 
     28;  which each lontitudes lines have the same latitude and each latitudes columns 
     29;  have the same longitude. 
    2930;  -  We supposed the data are located on a sphere, with a periodicity along 
    30 ;     the longitude. 
     31;  the longitude. 
    3132;  -  points located between the first/last 2 lines are interpolated 
    32 ;     using a imoms3 interpolation along the longitudinal direction and linear 
    33 ;     interpolation along the latitudinal direction 
    34 ;  -  points located out of the southern and northern boundaries are interpolated 
    35 ;     using a imoms3 interpolation only along the longitudinal direction. 
    36 ;  
     33;  using a imoms3 interpolation along the longitudinal direction and linear 
     34;  interpolation along the latitudinal direction 
     35;  - points located out of the southern and northern boundaries are interpolated 
     36;  using a imoms3 interpolation only along the longitudinal direction. 
     37; 
    3738; @history 
    38 ;  November 2005: Sebastien Masson (smasson\@lodyc.jussieu.fr)  
     39;  November 2005: Sebastien Masson (smasson\@lodyc.jussieu.fr) 
    3940;  March 2006: works for rectangular grids 
    4041; 
    … …  
    4950                                     , NONORTHERNLINE = nonorthernline, NOSOUTHERNLINE = nosouthernline 
    5051; 
    51   compile_opt idl2, strictarrsubs  
     52  compile_opt idl2, strictarrsubs 
    5253; 
    5354  alon = alonin 
    … …  
    6566  IF maxalon-minalon GE 360. THEN stop 
    6667; alon must be monotonically increasing 
    67   IF array_equal(sort(alon), lindgen(jpia)) NE 1 THEN BEGIN  
     68  IF array_equal(sort(alon), lindgen(jpia)) NE 1 THEN BEGIN 
    6869    shiftx = -(where(alon EQ min(alon)))[0] 
    6970    alon = shift(alon, shiftx) 
    … …  
    8990  IF total((step-step[0]) GE 1.e-6) NE 0 THEN noregy = 1 
    9091; 
    91   if keyword_set(nonorthernline) then BEGIN  
     92  if keyword_set(nonorthernline) then BEGIN 
    9293    jpja = jpja - 1L 
    9394    alat = alat[0: jpja-1L] 
    9495  ENDIF 
    95   if keyword_set(nosouthernline) then BEGIN  
     96  if keyword_set(nosouthernline) then BEGIN 
    9697    alat = alat[1: jpja-1L] 
    9798    jpja = jpja - 1L 
    … …  
    123124  indexlat = value_locate(alat, olat) 
    124125; 
    125 ; for the ocean points located below the atm line  
     126; for the ocean points located below the atm line 
    126127; jpja-2 and above the line 1 
    127128; for those points we can always find 16 neighbors 
    … …  
    131132  ilon = indexlon[short] 
    132133  ilat = indexlat[short] 
    133 ;  
    134   IF NOT keyword_set(noregy) THEN BEGIN  
     134; 
     135  IF NOT keyword_set(noregy) THEN BEGIN 
    135136    delta = alat[ilat+1L]-alat[ilat] 
    136137    IF max(abs(delta-delta[0])) GE 1.e-6 THEN stop 
    … …  
    147148    d2 = (alat[ilat+1L]-olat[short])/delta 
    148149    IF min(d2, max = ma) LE 0 THEN stop 
    149     IF ma GT 1 THEN stop   
     150    IF ma GT 1 THEN stop 
    150151    wy2 = imoms3(temporary(d2)) 
    151152    d3 = (alat[ilat+2L]-olat[short])/delta 
    152153    IF min(d3, max = ma) LE 1 THEN stop 
    153     IF ma GT 2 THEN stop   
     154    IF ma GT 2 THEN stop 
    154155    wy3 = imoms3(temporary(d3)) 
    155   ENDIF ELSE BEGIN  
     156  ENDIF ELSE BEGIN 
    156157    nele = n_elements(short) 
    157158    wy0 = fltarr(nele) 
    … …  
    169170      wy3[i] = imoms3(2-newlat) 
    170171    ENDFOR 
    171   ENDELSE  
     172  ENDELSE 
    172173; 
    173174  mi = min(wy0+wy1+wy2+wy3, max = ma) 
    … …  
    175176  IF abs(ma-1) GE 1.e-6 THEN stop 
    176177; 
    177   IF NOT keyword_set(noregx) THEN BEGIN  
     178  IF NOT keyword_set(noregx) THEN BEGIN 
    178179    delta = alon[ilon]-alon[ilon-1L] 
    179180    IF max(abs(delta-delta[0])) GE 1.e-6 THEN stop 
    … …  
    194195    d3 = (alon[ilon+2L]-olon[short])/delta 
    195196    IF min(d3, max = ma) LE 1 THEN stop 
    196     IF ma GT 2 THEN stop   
     197    IF ma GT 2 THEN stop 
    197198    wx3 = imoms3(temporary(d3)) 
    198   ENDIF ELSE BEGIN  
     199  ENDIF ELSE BEGIN 
    199200    nele = n_elements(short) 
    200201    wx0 = fltarr(nele) 
    … …  
    212213      wx3[i] = imoms3(2-newlon) 
    213214    ENDFOR 
    214   ENDELSE  
     215  ENDELSE 
    215216; 
    216217  mi = min(wx0+wx1+wx2+wx3, max = ma) 
    … …  
    220221; line 0 
    221222  xaddr[0, short] = ilon - 1L 
    222   xaddr[1, short] = ilon   
     223  xaddr[1, short] = ilon 
    223224  xaddr[2, short] = ilon + 1L 
    224225  xaddr[3, short] = ilon + 2L 
    … …  
    233234; line 1 
    234235  xaddr[4, short] = ilon - 1L 
    235   xaddr[5, short] = ilon   
     236  xaddr[5, short] = ilon 
    236237  xaddr[6, short] = ilon + 1L 
    237238  xaddr[7, short] = ilon + 2L 
    238   yaddr[4, short] = ilat   
    239   yaddr[5, short] = ilat   
    240   yaddr[6, short] = ilat   
    241   yaddr[7, short] = ilat   
     239  yaddr[4, short] = ilat 
     240  yaddr[5, short] = ilat 
     241  yaddr[6, short] = ilat 
     242  yaddr[7, short] = ilat 
    242243  weig[4, short] = wx0 * wy1 
    243244  weig[5, short] = wx1 * wy1 
    … …  
    246247; line 2 
    247248  xaddr[8, short] = ilon - 1L 
    248   xaddr[9, short] = ilon   
     249  xaddr[9, short] = ilon 
    249250  xaddr[10, short] = ilon + 1L 
    250251  xaddr[11, short] = ilon + 2L 
    … …  
    259260; line 3 
    260261  xaddr[12, short] = ilon - 1L 
    261   xaddr[13, short] = ilon   
     262  xaddr[13, short] = ilon 
    262263  xaddr[14, short] = ilon + 1L 
    263264  xaddr[15, short] = ilon + 2L 
    … …  
    275276  IF abs(ma-1) GE 1.e-6 THEN stop 
    276277; 
    277 ; for the ocean points located between the atm lines  
     278; for the ocean points located between the atm lines 
    278279; jpja-2 and jpja-1 or between the atm lines 0 and 1 
    279280; linear interpolation between line 1 and line 2 
    … …  
    285286; 
    286287    delta = alat[ilat+1L]-alat[ilat] 
    287     IF NOT keyword_set(noregy) THEN BEGIN  
     288    IF NOT keyword_set(noregy) THEN BEGIN 
    288289      IF max(abs(delta-delta[0])) GE 1.e-6 THEN stop 
    289       delta = delta[0]  
    290     ENDIF  
     290      delta = delta[0] 
     291    ENDIF 
    291292; 
    292293    d1 = (alat[ilat   ]-olat[short])/delta 
    … …  
    296297    d2 = (alat[ilat+1L]-olat[short])/delta 
    297298    IF min(d2, max = ma) LE 0 THEN stop 
    298     IF ma GT 1 THEN stop   
     299    IF ma GT 1 THEN stop 
    299300    wy2 = 1.- temporary(d2) 
    300301; 
    … …  
    303304    IF abs(ma-1) GE 1.e-6 THEN stop 
    304305; but imoms3 along the longitude 
    305     IF NOT keyword_set(noregx) THEN BEGIN  
     306    IF NOT keyword_set(noregx) THEN BEGIN 
    306307      delta = alon[ilon]-alon[ilon-1L] 
    307308      IF max(abs(delta-delta[0])) GE 1.e-6 THEN stop 
    … …  
    322323      d3 = (alon[ilon+2L]-olon[short])/delta 
    323324      IF min(d3, max = ma) LE 1 THEN stop 
    324       IF ma GT 2 THEN stop   
     325      IF ma GT 2 THEN stop 
    325326      wx3 = imoms3(temporary(d3)) 
    326     ENDIF ELSE BEGIN  
     327    ENDIF ELSE BEGIN 
    327328      nele = n_elements(short) 
    328329      wx0 = fltarr(nele) 
    … …  
    340341        wx3[i] = imoms3(2-newlon) 
    341342      ENDFOR 
    342     ENDELSE  
     343    ENDELSE 
    343344; 
    344345    mi = min(wx0+wx1+wx2+wx3, max = ma) 
    … …  
    347348; line 1 
    348349    xaddr[0, short] = ilon - 1L 
    349     xaddr[1, short] = ilon   
     350    xaddr[1, short] = ilon 
    350351    xaddr[2, short] = ilon + 1L 
    351352    xaddr[3, short] = ilon + 2L 
    352     yaddr[0, short] = ilat   
    353     yaddr[1, short] = ilat   
    354     yaddr[2, short] = ilat   
    355     yaddr[3, short] = ilat   
     353    yaddr[0, short] = ilat 
     354    yaddr[1, short] = ilat 
     355    yaddr[2, short] = ilat 
     356    yaddr[3, short] = ilat 
    356357    weig[0, short] = wx0 * wy1 
    357358    weig[1, short] = wx1 * wy1 
    … …  
    360361; line 2 
    361362    xaddr[4, short] = ilon - 1L 
    362     xaddr[5, short] = ilon   
     363    xaddr[5, short] = ilon 
    363364    xaddr[6, short] = ilon + 1L 
    364365    xaddr[7, short] = ilon + 2L 
    … …  
    379380; 
    380381; for the ocean points located below the line 0 
    381 ; Interpolation only along the longitude  
     382; Interpolation only along the longitude 
    382383; 
    383384  short = where(indexlat EQ -1) 
    … …  
    385386    ilon = indexlon[short] 
    386387; 
    387     IF NOT keyword_set(noregx) THEN BEGIN  
     388    IF NOT keyword_set(noregx) THEN BEGIN 
    388389      delta = alon[ilon]-alon[ilon-1L] 
    389390      IF max(abs(delta-delta[0])) GE 1.e-6 THEN stop 
    … …  
    404405      d3 = (alon[ilon+2L]-olon[short])/delta 
    405406      IF min(d3, max = ma) LE 1 THEN stop 
    406       IF ma GT 2 THEN stop   
     407      IF ma GT 2 THEN stop 
    407408      wx3 = imoms3(temporary(d3)) 
    408     ENDIF ELSE BEGIN  
     409    ENDIF ELSE BEGIN 
    409410      nele = n_elements(short) 
    410411      wx0 = fltarr(nele) 
    … …  
    422423        wx3[i] = imoms3(2-newlon) 
    423424      ENDFOR 
    424     ENDELSE  
     425    ENDELSE 
    425426; 
    426427    mi = min(wx0+wx1+wx2+wx3, max = ma) 
    … …  
    429430; line 1 
    430431    xaddr[0, short] = ilon - 1L 
    431     xaddr[1, short] = ilon   
     432    xaddr[1, short] = ilon 
    432433    xaddr[2, short] = ilon + 1L 
    433434    xaddr[3, short] = ilon + 2L 
    434     yaddr[0:3, short] = 0  
     435    yaddr[0:3, short] = 0 
    435436    weig[0, short] = wx0 
    436437    weig[1, short] = wx1 
    … …  
    444445  ENDIF 
    445446; 
    446 ; for the ocean points located above jpia-1  
    447 ; Interpolation only along the longitude  
     447; for the ocean points located above jpia-1 
     448; Interpolation only along the longitude 
    448449; 
    449450  short = where(indexlat EQ jpja-1L) 
    … …  
    451452    ilon = indexlon[short] 
    452453; 
    453     IF NOT keyword_set(noregx) THEN BEGIN  
     454    IF NOT keyword_set(noregx) THEN BEGIN 
    454455      delta = alon[ilon]-alon[ilon-1L] 
    455456      IF max(abs(delta-delta[0])) GE 1.e-6 THEN stop 
    … …  
    470471      d3 = (alon[ilon+2L]-olon[short])/delta 
    471472      IF min(d3, max = ma) LE 1 THEN stop 
    472       IF ma GT 2 THEN stop   
     473      IF ma GT 2 THEN stop 
    473474      wx3 = imoms3(temporary(d3)) 
    474     ENDIF ELSE BEGIN  
     475    ENDIF ELSE BEGIN 
    475476      nele = n_elements(short) 
    476477      wx0 = fltarr(nele) 
    … …  
    488489        wx3[i] = imoms3(2-newlon) 
    489490      ENDFOR 
    490     ENDELSE  
     491    ENDELSE 
    491492; 
    492493    mi = min(wx0+wx1+wx2+wx3, max = ma) 
    … …  
    495496; line 1 
    496497    xaddr[0, short] = ilon-1L 
    497     xaddr[1, short] = ilon   
     498    xaddr[1, short] = ilon 
    498499    xaddr[2, short] = ilon+1L 
    499500    xaddr[3, short] = ilon+2L 
    … …  
    531532  IF revy EQ 1 THEN yaddr = jpja - 1L - temporary(yaddr) 
    532533;                        ; 
    533   addr = temporary(yaddr)*jpia+temporary(xaddr)  
     534  addr = temporary(yaddr)*jpia+temporary(xaddr) 
    534535; 
    535536  RETURN 
Note: See TracChangeset for help on using the changeset viewer.