Changeset 69 for trunk/Interpolation


Ignore:
Timestamp:
05/11/06 12:35:53 (18 years ago)
Author:
smasson
Message:

debug + new xxx

Location:
trunk/Interpolation
Files:
3 added
7 edited

Legend:

Unmodified
Added
Removed
  • trunk/Interpolation/.idlwave_catalog

    r68 r69  
    22;; IDLWAVE catalog for library saxo 
    33;; Automatically Generated -- do not edit. 
    4 ;; Created by idlwave_catalog on Tue Jan 24 17:06:19 2006 
     4;; Created by idlwave_catalog on Thu May  4 08:42:55 2006 
    55;; 
    66(setq idlwave-library-catalog-libname "saxo") 
     
    1010   ("clickincell" fun nil (lib "clickincell.pro" nil "saxo") "Result = %s" (nil ("_EXTRA") ("CELLTYPE") ("COLOR") ("DRAWCELL") ("IJ") ("ORIGINAL"))) 
    1111   ("compute_fromreg_bilinear_weigaddr" pro nil (lib "compute_fromreg_bilinear_weigaddr.pro" nil "saxo") "%s, alonin, alatin, olonin, olat, weig, addr" (nil ("NONORTHERNLINE") ("NOSOUTHERNLINE"))) 
    12    ("compute_fromreg_imoms3_weigaddr" pro nil (lib "compute_fromreg_imoms3_weigaddr.pro" nil "saxo") "%s, alonin, alatin, olonin, olat, weig, addr" (nil ("FORCE") ("NONORTHERNLINE") ("NOSOUTHERNLINE"))) 
     12   ("compute_fromreg_imoms3_weigaddr" pro nil (lib "compute_fromreg_imoms3_weigaddr.pro" nil "saxo") "%s, alonin, alatin, olonin, olat, weig, addr" (nil ("NONORTHERNLINE") ("NOSOUTHERNLINE"))) 
    1313   ("cutpar" fun nil (lib "cutpar.pro" nil "saxo") "Result = %s(x0, y0, x1, y1, x2, y2, x3, y3, n)" (nil ("endpoints") ("onsphere"))) 
    1414   ("cutsegment" fun nil (lib "cutsegment.pro" nil "saxo") "Result = %s(x0, y0, x1, y1, n)" (nil ("endpoints") ("onsphere"))) 
     
    2323   ("neighbor" fun nil (lib "neighbor.pro" nil "saxo") "Result = %s(p0lon, p0lat, neighlon, neighlat)" (nil ("distance") ("radians") ("sphere"))) 
    2424   ("quadrilateral2square" fun nil (lib "quadrilateral2square.pro" nil "saxo") "Result = %s(x0in, y0in, x1in, y1in, x2in, y2in, x3in, y3in, xxin, yyin)" (nil ("PERF"))) 
     25   ("spl_fstdrv" fun nil (lib "spl_fstdrv.pro" nil "saxo") "Result = %s(x, y, yscd, x2)" (nil)) 
     26   ("pure_concave" fun nil (lib "spl_incr.pro" nil "saxo") "Result = %s(x1, x2, y1, y2, der2, x)" (nil)) 
     27   ("pure_convex" fun nil (lib "spl_incr.pro" nil "saxo") "Result = %s(x1, x2, y1, y2, der2, x)" (nil)) 
     28   ("spl_incr" fun nil (lib "spl_incr.pro" nil "saxo") "Result = %s(x, y, x2)" (nil ("YP0") ("YPN_1"))) 
     29   ("spl_keep_mean" fun nil (lib "spl_keep_mean.pro" nil "saxo") "Result = %s(x, yin, x2)" (nil ("GE0") ("YP0") ("YPN_1"))) 
    2530   ("square2quadrilateral" fun nil (lib "square2quadrilateral.pro" nil "saxo") "Result = %s(x0in, y0in, x1in, y1in, x2in, y2in, x3in, y3in, xxin, yyin)" (nil)) 
    2631   ("testinterp" pro nil (lib "testinterp.pro" nil "saxo") "%s" (nil)))) 
  • trunk/Interpolation/clickincell.pro

    r59 r69  
    1919;     CELLTYPE = 'T', 'W', 'U', 'V' or 'F': This this the type of point 
    2020;     that is located in the center of the cell which the click is 
    21 ;     located. defalut is T type of cell (with corner defined by F 
     21;     located. default is T type of cell (with corner defined by F 
    2222;     points). 
    2323; 
  • trunk/Interpolation/compute_fromreg_imoms3_weigaddr.pro

    r59 r69  
    3131; 
    3232; RESTRICTIONS: 
    33 ;  -  the input grid must be a "regular grid", defined as a grid for which 
    34 ;     lontitudes and latitudes are regularly spaced. 
     33;  -  the input grid must be a "regular/rectangular grid", defined as a grid for 
     34;     which each lontitudes lines have the same latitude and each latitudes columns 
     35;     have the same longitude. 
    3536;  -  We supposed the data are located on a sphere, with a periodicity along 
    3637;     the longitude. 
     
    4546; MODIFICATION HISTORY: 
    4647;  November 2005: Sebastien Masson (smasson@lodyc.jussieu.fr)  
    47 ;  
     48;  March 2006: works for rectangular grids 
    4849;- 
    4950; 
     
    5253; 
    5354PRO compute_fromreg_imoms3_weigaddr, alonin, alatin, olonin, olat, weig, addr $ 
    54   , NONORTHERNLINE = nonorthernline, NOSOUTHERNLINE = nosouthernline, FORCE = force 
     55                                     , NONORTHERNLINE = nonorthernline, NOSOUTHERNLINE = nosouthernline 
    5556; 
    5657  compile_opt strictarr, strictarrsubs  
     
    7576    IF array_equal(sort(alon), lindgen(jpia)) NE 1 THEN stop 
    7677  ENDIF ELSE shiftx = 0 
    77 ; alon must be regularly spaced 
     78; alon is it regularly spaced? 
    7879  step = alon-shift(alon, 1) 
    7980  step[0] = step[0] + 360. 
    80   IF total((step-step[0]) GE 1.e-6) NE 0 THEN BEGIN  
    81     print,  'input longitude must be regularly spaced, we stop...' 
    82     stop 
    83   ENDIF 
     81  IF total((step-step[0]) GE 1.e-6) NE 0 THEN noregx = 1 
    8482; we extend the longitude range of alon (-> easy interpolation even 
    8583; near minalon et maxalon) 
     
    9290; alat must be monotonically increasing 
    9391  IF array_equal(sort(alat), lindgen(jpja)) NE 1 THEN stop 
    94 ; alat must be regularly spaced 
     92; alat is it regularly spaced? 
    9593  step = alat-shift(alat, 1) 
    9694  step = step[1:jpja - 1L] 
    97   IF total((step-step[0]) GE 1.e-6) NE 0 AND NOT keyword_set(force) THEN BEGIN  
    98     print,  'input latitude must be regularly spaced, we stop...' 
    99     stop 
    100   ENDIF 
     95  IF total((step-step[0]) GE 1.e-6) NE 0 THEN noregy = 1 
    10196; 
    10297  if keyword_set(nonorthernline) then BEGIN  
     
    143138  ilat = indexlat[short] 
    144139;  
    145   delta = alat[ilat+1L]-alat[ilat] 
    146   IF max(abs(delta-delta[0])) GE 1.e-6 THEN stop 
    147   delta = delta[0] 
    148 ; 
    149   d0 = (alat[ilat-1]-olat[short])/delta 
    150   IF min(d0, max = ma) LE -2 THEN stop 
    151   IF ma GT -1 THEN stop 
    152   wy0 = imoms3(temporary(d0)) 
    153   d1 = (alat[ilat  ]-olat[short])/delta 
    154   IF min(d1, max = ma) LE -1 THEN stop 
    155   IF ma GT 0 THEN stop 
    156   wy1 = imoms3(temporary(d1)) 
    157   d2 = (alat[ilat+1]-olat[short])/delta 
    158   IF min(d2, max = ma) LE 0 THEN stop 
    159   IF ma GT 1 THEN stop   
    160   wy2 = imoms3(temporary(d2)) 
    161   d3 = (alat[ilat+2]-olat[short])/delta 
    162   IF min(d3, max = ma) LE 1 THEN stop 
    163   IF ma GT 2 THEN stop   
    164   wy3 = imoms3(temporary(d3)) 
     140  IF NOT keyword_set(noregy) THEN BEGIN  
     141    delta = alat[ilat+1L]-alat[ilat] 
     142    IF max(abs(delta-delta[0])) GE 1.e-6 THEN stop 
     143    delta = delta[0] 
     144; 
     145    d0 = (alat[ilat-1L]-olat[short])/delta 
     146    IF min(d0, max = ma) LE -2 THEN stop 
     147    IF ma GT -1 THEN stop 
     148    wy0 = imoms3(temporary(d0)) 
     149    d1 = (alat[ilat   ]-olat[short])/delta 
     150    IF min(d1, max = ma) LE -1 THEN stop 
     151    IF ma GT 0 THEN stop 
     152    wy1 = imoms3(temporary(d1)) 
     153    d2 = (alat[ilat+1L]-olat[short])/delta 
     154    IF min(d2, max = ma) LE 0 THEN stop 
     155    IF ma GT 1 THEN stop   
     156    wy2 = imoms3(temporary(d2)) 
     157    d3 = (alat[ilat+2L]-olat[short])/delta 
     158    IF min(d3, max = ma) LE 1 THEN stop 
     159    IF ma GT 2 THEN stop   
     160    wy3 = imoms3(temporary(d3)) 
     161  ENDIF ELSE BEGIN  
     162    nele = n_elements(short) 
     163    wy0 = fltarr(nele) 
     164    wy1 = fltarr(nele) 
     165    wy2 = fltarr(nele) 
     166    wy3 = fltarr(nele) 
     167    FOR i = 0L, nele-1 DO BEGIN 
     168      IF i MOD 10000 EQ 0 THEN print, i 
     169      newlat = spl_incr(alat[ilat[i]-1L:ilat[i]+2L], [-1., 0., 1., 2.], olat[short[i]]) 
     170      IF newlat LE 0 THEN stop 
     171      IF newlat GT 1 THEN stop 
     172      wy0[i] = imoms3(newlat+1) 
     173      wy1[i] = imoms3(newlat) 
     174      wy2[i] = imoms3(1-newlat) 
     175      wy3[i] = imoms3(2-newlat) 
     176    ENDFOR 
     177  ENDELSE  
    165178; 
    166179  mi = min(wy0+wy1+wy2+wy3, max = ma) 
     
    168181  IF abs(ma-1) GE 1.e-6 THEN stop 
    169182; 
    170   delta = alon[ilon]-alon[ilon-1] 
    171   IF max(abs(delta-delta[0])) GE 1.e-6 THEN stop 
    172   delta = delta[0] 
    173 ; 
    174   d0 = (alon[ilon-1]-olon[short])/delta 
    175   IF min(d0, max = ma) LE -2 THEN stop 
    176   IF ma GT -1 THEN stop 
    177   wx0 = imoms3(temporary(d0)) 
    178   d1 = (alon[ilon  ]-olon[short])/delta 
    179   IF min(d1, max = ma) LE -1 THEN stop 
    180   IF ma GT 0 THEN stop 
    181   wx1 = imoms3(temporary(d1)) 
    182   d2 = (alon[ilon+1]-olon[short])/delta 
    183   IF min(d2, max = ma) LE 0 THEN stop 
    184   IF ma GT 1 THEN stop 
    185   wx2 = imoms3(temporary(d2)) 
    186   d3 = (alon[ilon+2]-olon[short])/delta 
    187   IF min(d3, max = ma) LE 1 THEN stop 
    188   IF ma GT 2 THEN stop   
    189   wx3 = imoms3(temporary(d3)) 
     183  IF NOT keyword_set(noregx) THEN BEGIN  
     184    delta = alon[ilon]-alon[ilon-1L] 
     185    IF max(abs(delta-delta[0])) GE 1.e-6 THEN stop 
     186    delta = delta[0] 
     187; 
     188    d0 = (alon[ilon-1L]-olon[short])/delta 
     189    IF min(d0, max = ma) LE -2 THEN stop 
     190    IF ma GT -1 THEN stop 
     191    wx0 = imoms3(temporary(d0)) 
     192    d1 = (alon[ilon   ]-olon[short])/delta 
     193    IF min(d1, max = ma) LE -1 THEN stop 
     194    IF ma GT 0 THEN stop 
     195    wx1 = imoms3(temporary(d1)) 
     196    d2 = (alon[ilon+1L]-olon[short])/delta 
     197    IF min(d2, max = ma) LE 0 THEN stop 
     198    IF ma GT 1 THEN stop 
     199    wx2 = imoms3(temporary(d2)) 
     200    d3 = (alon[ilon+2L]-olon[short])/delta 
     201    IF min(d3, max = ma) LE 1 THEN stop 
     202    IF ma GT 2 THEN stop   
     203    wx3 = imoms3(temporary(d3)) 
     204  ENDIF ELSE BEGIN  
     205    nele = n_elements(short) 
     206    wx0 = fltarr(nele) 
     207    wx1 = fltarr(nele) 
     208    wx2 = fltarr(nele) 
     209    wx3 = fltarr(nele) 
     210    FOR i = 0L, nele-1 DO BEGIN 
     211      IF i MOD 10000 EQ 0 THEN print, i 
     212      newlon = spl_incr(alon[ilon[i]-1L:ilon[i]+2L], [-1., 0., 1., 2.], olon[short[i]]) 
     213      IF newlon LE 0 THEN stop 
     214      IF newlon GT 1 THEN stop 
     215      wx0[i] = imoms3(newlon+1) 
     216      wx1[i] = imoms3(newlon) 
     217      wx2[i] = imoms3(1-newlon) 
     218      wx3[i] = imoms3(2-newlon) 
     219    ENDFOR 
     220  ENDELSE  
    190221; 
    191222  mi = min(wx0+wx1+wx2+wx3, max = ma) 
     
    224255  xaddr[10, short] = ilon + 1L 
    225256  xaddr[11, short] = ilon + 2L 
    226   yaddr[8 , short] = ilat + 1L 
    227   yaddr[9 , short] = yaddr[8 , short] 
    228   yaddr[10, short] = yaddr[8 , short] 
    229   yaddr[11, short] = yaddr[8 , short] 
    230   weig[8 , short] = wx0 * wy2 
    231   weig[9 , short] = wx1 * wy2 
     257  yaddr[8, short] = ilat + 1L 
     258  yaddr[9, short] = yaddr[8, short] 
     259  yaddr[10, short] = yaddr[8, short] 
     260  yaddr[11, short] = yaddr[8, short] 
     261  weig[8, short] = wx0 * wy2 
     262  weig[9, short] = wx1 * wy2 
    232263  weig[10, short] = wx2 * wy2 
    233264  weig[11, short] = wx3 * wy2 
     
    254285; linear interpolation between line 1 and line 2 
    255286; 
    256   short = where(indexlat EQ jpja-2 OR indexlat EQ 0) 
     287  short = where(indexlat EQ jpja-2L OR indexlat EQ 0) 
    257288  IF short[0] NE -1 THEN BEGIN 
    258289    ilon = indexlon[short] 
    259290    ilat = indexlat[short] 
    260291; 
    261     delta = alat[ilat+1]-alat[ilat] 
    262     IF max(abs(delta-delta[0])) GE 1.e-6 THEN stop 
    263     delta = delta[0] 
    264 ; 
    265     d1 = (alat[ilat  ]-olat[short])/delta 
     292    delta = alat[ilat+1L]-alat[ilat] 
     293    IF NOT keyword_set(noregy) THEN BEGIN  
     294      IF max(abs(delta-delta[0])) GE 1.e-6 THEN stop 
     295      delta = delta[0]  
     296    ENDIF  
     297; 
     298    d1 = (alat[ilat   ]-olat[short])/delta 
    266299    IF min(d1, max = ma) LE -1 THEN stop 
    267300    IF ma GT 0 THEN stop 
    268301    wy1 = 1.+ temporary(d1) 
    269     d2 = (alat[ilat+1]-olat[short])/delta 
     302    d2 = (alat[ilat+1L]-olat[short])/delta 
    270303    IF min(d2, max = ma) LE 0 THEN stop 
    271304    IF ma GT 1 THEN stop   
     
    276309    IF abs(ma-1) GE 1.e-6 THEN stop 
    277310; but imoms3 along the longitude 
    278     delta = alon[ilon]-alon[ilon-1] 
    279     IF max(abs(delta-delta[0])) GE 1.e-6 THEN stop 
    280     delta = delta[0] 
    281 ; 
    282     d0 = (alon[ilon-1]-olon[short])/delta 
    283     IF min(d0, max = ma) LE -2 THEN stop 
    284     IF ma GT -1 THEN stop 
    285     wx0 = imoms3(temporary(d0)) 
    286     d1 = (alon[ilon  ]-olon[short])/delta 
    287     IF min(d1, max = ma) LE -1 THEN stop 
    288     IF ma GT 0 THEN stop 
    289     wx1 = imoms3(temporary(d1)) 
    290     d2 = (alon[ilon+1]-olon[short])/delta 
    291     IF min(d2, max = ma) LE 0 THEN stop 
    292     IF ma GT 1 THEN stop 
    293     wx2 = imoms3(temporary(d2)) 
    294     d3 = (alon[ilon+2]-olon[short])/delta 
    295     IF min(d3, max = ma) LE 1 THEN stop 
    296     IF ma GT 2 THEN stop   
    297     wx3 = imoms3(temporary(d3)) 
     311    IF NOT keyword_set(noregx) THEN BEGIN  
     312      delta = alon[ilon]-alon[ilon-1L] 
     313      IF max(abs(delta-delta[0])) GE 1.e-6 THEN stop 
     314      delta = delta[0] 
     315; 
     316      d0 = (alon[ilon-1L]-olon[short])/delta 
     317      IF min(d0, max = ma) LE -2 THEN stop 
     318      IF ma GT -1 THEN stop 
     319      wx0 = imoms3(temporary(d0)) 
     320      d1 = (alon[ilon   ]-olon[short])/delta 
     321      IF min(d1, max = ma) LE -1 THEN stop 
     322      IF ma GT 0 THEN stop 
     323      wx1 = imoms3(temporary(d1)) 
     324      d2 = (alon[ilon+1L]-olon[short])/delta 
     325      IF min(d2, max = ma) LE 0 THEN stop 
     326      IF ma GT 1 THEN stop 
     327      wx2 = imoms3(temporary(d2)) 
     328      d3 = (alon[ilon+2L]-olon[short])/delta 
     329      IF min(d3, max = ma) LE 1 THEN stop 
     330      IF ma GT 2 THEN stop   
     331      wx3 = imoms3(temporary(d3)) 
     332    ENDIF ELSE BEGIN  
     333      nele = n_elements(short) 
     334      wx0 = fltarr(nele) 
     335      wx1 = fltarr(nele) 
     336      wx2 = fltarr(nele) 
     337      wx3 = fltarr(nele) 
     338      FOR i = 0L, nele-1 DO BEGIN 
     339        IF i MOD 10000 EQ 0 THEN print, i 
     340        newlon = spl_incr(alon[ilon[i]-1L:ilon[i]+2L], [-1., 0., 1., 2.], olon[short[i]]) 
     341        IF newlon LE 0 THEN stop 
     342        IF newlon GT 1 THEN stop 
     343        wx0[i] = imoms3(newlon+1) 
     344        wx1[i] = imoms3(newlon) 
     345        wx2[i] = imoms3(1-newlon) 
     346        wx3[i] = imoms3(2-newlon) 
     347      ENDFOR 
     348    ENDELSE  
    298349; 
    299350    mi = min(wx0+wx1+wx2+wx3, max = ma) 
     
    340391    ilon = indexlon[short] 
    341392; 
    342     delta = alon[ilon]-alon[ilon-1] 
    343     IF max(abs(delta-delta[0])) GE 1.e-6 THEN stop 
    344     delta = delta[0] 
    345 ; 
    346     d0 = (alon[ilon-1]-olon[short])/delta 
    347     IF min(d0, max = ma) LE -2 THEN stop 
    348     IF ma GT -1 THEN stop 
    349     wx0 = imoms3(temporary(d0)) 
    350     d1 = (alon[ilon  ]-olon[short])/delta 
    351     IF min(d1, max = ma) LE -1 THEN stop 
    352     IF ma GT 0 THEN stop 
    353     wx1 = imoms3(temporary(d1)) 
    354     d2 = (alon[ilon+1]-olon[short])/delta 
    355     IF min(d2, max = ma) LE 0 THEN stop 
    356     IF ma GT 1 THEN stop 
    357     wx2 = imoms3(temporary(d2)) 
    358     d3 = (alon[ilon+2]-olon[short])/delta 
    359     IF min(d3, max = ma) LE 1 THEN stop 
    360     IF ma GT 2 THEN stop   
    361     wx3 = imoms3(temporary(d3)) 
     393    IF NOT keyword_set(noregx) THEN BEGIN  
     394      delta = alon[ilon]-alon[ilon-1L] 
     395      IF max(abs(delta-delta[0])) GE 1.e-6 THEN stop 
     396      delta = delta[0] 
     397; 
     398      d0 = (alon[ilon-1L]-olon[short])/delta 
     399      IF min(d0, max = ma) LE -2 THEN stop 
     400      IF ma GT -1 THEN stop 
     401      wx0 = imoms3(temporary(d0)) 
     402      d1 = (alon[ilon   ]-olon[short])/delta 
     403      IF min(d1, max = ma) LE -1 THEN stop 
     404      IF ma GT 0 THEN stop 
     405      wx1 = imoms3(temporary(d1)) 
     406      d2 = (alon[ilon+1L]-olon[short])/delta 
     407      IF min(d2, max = ma) LE 0 THEN stop 
     408      IF ma GT 1 THEN stop 
     409      wx2 = imoms3(temporary(d2)) 
     410      d3 = (alon[ilon+2L]-olon[short])/delta 
     411      IF min(d3, max = ma) LE 1 THEN stop 
     412      IF ma GT 2 THEN stop   
     413      wx3 = imoms3(temporary(d3)) 
     414    ENDIF ELSE BEGIN  
     415      nele = n_elements(short) 
     416      wx0 = fltarr(nele) 
     417      wx1 = fltarr(nele) 
     418      wx2 = fltarr(nele) 
     419      wx3 = fltarr(nele) 
     420      FOR i = 0L, nele-1 DO BEGIN 
     421        IF i MOD 10000 EQ 0 THEN print, i 
     422        newlon = spl_incr(alon[ilon[i]-1L:ilon[i]+2L], [-1., 0., 1., 2.], olon[short[i]]) 
     423        IF newlon LE 0 THEN stop 
     424        IF newlon GT 1 THEN stop 
     425        wx0[i] = imoms3(newlon+1) 
     426        wx1[i] = imoms3(newlon) 
     427        wx2[i] = imoms3(1-newlon) 
     428        wx3[i] = imoms3(2-newlon) 
     429      ENDFOR 
     430    ENDELSE  
    362431; 
    363432    mi = min(wx0+wx1+wx2+wx3, max = ma) 
     
    384453; Interpolation only along the longitude  
    385454; 
    386   short = where(indexlat EQ jpja-1) 
     455  short = where(indexlat EQ jpja-1L) 
    387456  IF short[0] NE -1 THEN BEGIN 
    388457    ilon = indexlon[short] 
    389458; 
    390     delta = alon[ilon]-alon[ilon-1] 
    391     IF max(abs(delta-delta[0])) GE 1.e-6 THEN stop 
    392     delta = delta[0] 
    393 ; 
    394     d0 = (alon[ilon-1]-olon[short])/delta 
    395     IF min(d0, max = ma) LE -2 THEN stop 
    396     IF ma GT -1 THEN stop 
    397     wx0 = imoms3(temporary(d0)) 
    398     d1 = (alon[ilon  ]-olon[short])/delta 
    399     IF min(d1, max = ma) LE -1 THEN stop 
    400     IF ma GT 0 THEN stop 
    401     wx1 = imoms3(temporary(d1)) 
    402     d2 = (alon[ilon+1]-olon[short])/delta 
    403     IF min(d2, max = ma) LE 0 THEN stop 
    404     IF ma GT 1 THEN stop 
    405     wx2 = imoms3(temporary(d2)) 
    406     d3 = (alon[ilon+2]-olon[short])/delta 
    407     IF min(d3, max = ma) LE 1 THEN stop 
    408     IF ma GT 2 THEN stop   
    409     wx3 = imoms3(temporary(d3)) 
     459    IF NOT keyword_set(noregx) THEN BEGIN  
     460      delta = alon[ilon]-alon[ilon-1L] 
     461      IF max(abs(delta-delta[0])) GE 1.e-6 THEN stop 
     462      delta = delta[0] 
     463; 
     464      d0 = (alon[ilon-1L]-olon[short])/delta 
     465      IF min(d0, max = ma) LE -2 THEN stop 
     466      IF ma GT -1 THEN stop 
     467      wx0 = imoms3(temporary(d0)) 
     468      d1 = (alon[ilon   ]-olon[short])/delta 
     469      IF min(d1, max = ma) LE -1 THEN stop 
     470      IF ma GT 0 THEN stop 
     471      wx1 = imoms3(temporary(d1)) 
     472      d2 = (alon[ilon+1L]-olon[short])/delta 
     473      IF min(d2, max = ma) LE 0 THEN stop 
     474      IF ma GT 1 THEN stop 
     475      wx2 = imoms3(temporary(d2)) 
     476      d3 = (alon[ilon+2L]-olon[short])/delta 
     477      IF min(d3, max = ma) LE 1 THEN stop 
     478      IF ma GT 2 THEN stop   
     479      wx3 = imoms3(temporary(d3)) 
     480    ENDIF ELSE BEGIN  
     481      nele = n_elements(short) 
     482      wx0 = fltarr(nele) 
     483      wx1 = fltarr(nele) 
     484      wx2 = fltarr(nele) 
     485      wx3 = fltarr(nele) 
     486      FOR i = 0L, nele-1 DO BEGIN 
     487        IF i MOD 10000 EQ 0 THEN print, i 
     488        newlon = spl_incr(alon[ilon[i]-1L:ilon[i]+2L], [-1., 0., 1., 2.], olon[short[i]]) 
     489        IF newlon LE 0 THEN stop 
     490        IF newlon GT 1 THEN stop 
     491        wx0[i] = imoms3(newlon+1) 
     492        wx1[i] = imoms3(newlon) 
     493        wx2[i] = imoms3(1-newlon) 
     494        wx3[i] = imoms3(2-newlon) 
     495      ENDFOR 
     496    ENDELSE  
    410497; 
    411498    mi = min(wx0+wx1+wx2+wx3, max = ma) 
  • trunk/Interpolation/extrapolate.pro

    r59 r69  
    4242  sqrtinv = 1./sqrt(2) 
    4343; 
    44   cnt = 0 
     44  cnt = 1 
    4545; When we look for the coast line, we don't whant to select the 
    4646; borderlines of the array. -> we force the value of the mask for 
  • trunk/Interpolation/fromreg.pro

    r59 r69  
    22; NAME: fromreg 
    33; 
    4 ; PURPOSE: interpolate data from a "regular grid" to any grid. 
     4; PURPOSE: interpolate data from a "regular/rectangular grid" to any grid. 
    55;   2 metods availables: bilinear and imoms3  
    6 ;  
    7 ;   In the case on 'bilinear interpolation': 
    8 ;   A "regular grid" is defined as a grid for which each lontitudes lines have  
     6;   A "regular/rectangular grid" is defined as a grid for which each lontitudes lines have  
    97;   the same latitude and each latitudes columns have the same longitude. 
    10 ;    
    11 ;  
    12 ;   In the case on 'imoms3 interpolation': 
    13 ;   A "regular grid" is defined as a grid for which lontitudes and latitudes 
    14 ;   are regularly spaced. 
    158;    
    169; CATEGORY:interpolation 
  • trunk/Interpolation/get_gridparams.pro

    r59 r69  
    66;      and make sure it is 1D or 2D arrays  
    77;    
    8 ;   or 2) given longitude and latidude arrays get their dimensions and make sure  
     8;   or 2) given longitude and latitude arrays get their dimensions and make sure  
    99;      they are 1D or 2D arrays  
    1010; 
     
    1313; CALLING SEQUENCE: 
    1414;  
    15 ; 1) ncdf_get_gridparams, file, lonname, latname, lon, lat, jpi, jpj, n_dimensions 
     15; 1) get_gridparams, file, lonname, latname, lon, lat, jpi, jpj, n_dimensions 
    1616; 
    1717; or  
    1818; 
    19 ; 2) ncdf_get_gridparams, lon, lat, jpi, jpj, n_dimensions 
     19; 2) get_gridparams, lon, lat, jpi, jpj, n_dimensions 
    2020; 
    2121; INPUTS: 
     
    2929; 
    3030; 2) lon and lat: 1d or 2D arrays defining longitudes and latitudes.  
    31 ;    Note that these arrays are also outpust and can therefore be modified.  
     31;    Note that these arrays are also outputs and can therefore be modified.  
    3232 
    3333; KEYWORD PARAMETERS: none 
     
    3535; OUTPUTS: 
    3636;  lon the variable that will contain the longitudes 
    37 ;  lat the variable that will contain the latidutes 
     37;  lat the variable that will contain the latitudes 
    3838;  jpi the number of points in the longitudinal direction 
    3939;  jpj the number of points in the latitudinal direction 
  • trunk/Interpolation/imoms3.pro

    r59 r69  
    77  IF test1[0] NE -1 THEN BEGIN 
    88    xtmp = x[test1] 
    9     y[test1] = .5*xtmp^3 -xtmp^2 -.5*xtmp+1 
     9    y[test1] = 0.5*xtmp*xtmp*xtmp -xtmp*xtmp -0.5*xtmp + 1 
    1010  ENDIF 
    1111  test1 = where(x LT 2 AND x GE 1) 
    1212  IF test1[0] NE -1 THEN BEGIN 
    1313    xtmp = x[test1] 
    14     y[test1] = (-1./6.)*xtmp^3 + xtmp^2 - (11./6.)*xtmp+1 
     14    y[test1] = (-1./6.)*xtmp*xtmp*xtmp + xtmp*xtmp - (11./6.)*xtmp + 1 
    1515  ENDIF 
    1616 
Note: See TracChangeset for help on using the changeset viewer.