Changeset 44 for trunk/ToBeReviewed


Ignore:
Timestamp:
05/02/06 18:13:44 (18 years ago)
Author:
pinsard
Message:

upgrade of LECTURE/ReadWrite according to cerbere.lodyc.jussieu.fr: /usr/home/smasson/SAXO_RD/ : files

Location:
trunk/ToBeReviewed/LECTURE
Files:
3 added
7 copied

Legend:

Unmodified
Added
Removed

  • trunk/ToBeReviewed/LECTURE/litchamp.pro

    r42 r44  
    1313;   Au passage litchamp regarde les autres elements de la structure et 
    1414;   met a jour si besoin les variables globales qui se rapportent au 
    15 ;   champ: vargrid, varname, varunit, vardate, varexp , valmask, 
    16 ;   niveau et time  
     15;   champ: vargrid, varname, varunit, vardate, varexp , valmask et time  
    1716; 
    1817; CATEGORY:permet d''appeler plt, pltz, pltt ... avec un tableau ou une 
    … …  
    4342;             n  pour actualiser varname 
    4443;             m  pour actualiser valmask 
    45 ;             l  pour actualiser niveau (l pour level! et non pas n) 
    4644; 
    4745; KEYWORD PARAMETERS: 
    … …  
    8179; 
    8280; SIDE EFFECTS: actualise au besion les variables globales vargrid, 
    83 ; varname, varunit, vardate, varexp, valmask, niveau et time. 
     81; varname, varunit, vardate, varexp, valmask et time. 
    8482; 
    8583; RESTRICTIONS: 
    … …  
    150148            varexp = struct.(i) 
    151149         END 
    152          'l':BEGIN 
    153             if keyword_set(level) then return, struct.(i) 
    154             niveau = struct.(i) 
    155          END 
    156150         'm':BEGIN 
    157151            if keyword_set(mask) then return, struct.(i) 

  • trunk/ToBeReviewed/LECTURE/read_ncdf.pro

    r42 r44  
    2929; KEYWORD PARAMETERS: utilisables hors du contexte des widgets 
    3030; 
    31 ;        BOITE: contient la boite sur laquelle on doit faire la 
     31;        BOXZOOM: contient la boxzoom sur laquelle on doit faire la 
    3232;        lecture 
    3333;        FILENAME: string contennant le nom du fichier 
    34 ;        IODIRECTORY: string contennant le nom du repertoire ds lequel 
    35 ;        se trouve le fichier 
     34;        /INIT; to call automatically initncdf, filename and thus 
     35;        redefine all the grid parameters 
     36;        GRID='[UTVWF]' to specify the type of grid. Defaut is (1) 
     37;        based on the name of the file if the file ends by 
     38;        GRID[._][TUVFW].NC (not case sensible) or (2) T if case (1) 
     39;        is not found. 
     40;        IODIRECTORY;a string giving the name of iodirectory (see 
     41;        isafile.pro for all possibilities). default value is common  
     42;        variable iodir 
    3643;        TIMESTEP:activer pour specifier que debut et fin font 
    3744;        reference a des indices de l''axe du temps et non pas a des 
    3845;        dates.  
    3946;        TOUT: activer si on veut lire le ficher sur l''ensemble du 
    40 ;        domaine sans tenir compte du sous domaine definit par boite 
    41 ;        ou lon1,lon2,lat1,lat2,prof1,prof2. 
    42 ;        NOSTRUCT: activer si on ne veut pas que read_ncdf retourne 
     47;        domaine sans tenir compte du sous domaine definit par boxzoom 
     48;        ou lon1,lon2,lat1,lat2,vert1,vert2. 
     49;        NOSTRUCT: activer si on ne veut pas que read_ncdf reourne 
    4350;        une structure mais uniquement le tableau se rapportant au 
    4451;        champ.  
     52;        TIMEVAR: a string to define the name of the variable that 
     53;        contains the time axis. This keyword can be usefull if there 
     54;        is no unlimited dimension of if the time axis selected by defaut 
     55;        (the first 1D array with unlimited dimension is not the good one) 
     56;  
    4557; 
    4658; OUTPUTS:une stucture lisible par litchamp.pro ou un simple tableau 
    … …  
    6173;------------------------------------------------------------ 
    6274;------------------------------------------------------------ 
    63 FUNCTION read_ncdf,name,debut,fin, pour_etre_compatible, BOITE=boite, FILENAME = filename $ 
    64              , IODIRECTORY = iodirectory, PARENTIN = parentin, TIMESTEP = timestep $ 
    65              , TOUT = tout, NOSTRUCT = nostruct, CONT_NOFILL = CONT_NOFILL, _EXTRA = ex 
    66 @common 
     75FUNCTION read_ncdf, name, debut, fin, pour_etre_compatible, BOXZOOM = boxzoom, FILENAME = filename $ 
     76                    , PARENTIN = parentin, TIMESTEP = timestep, TIMEVAR = timevar $ 
     77                    , TOUT = tout, NOSTRUCT = nostruct, CONT_NOFILL = CONT_NOFILL, INIT = init $ 
     78                    , GRID = grid, FBASE2TBASE = fbase2tbase, _EXTRA = ex 
     79;--------------------------------------------------------- 
     80@cm_4mesh 
     81@cm_4data 
     82@cm_4cal 
     83  IF NOT keyword_set(key_forgetold) THEN BEGIN 
     84@updatenew 
     85@updatekwd 
     86  ENDIF 
    6787;------------------------------------------------------------ 
    6888; we find the filename. 
    6989;------------------------------------------------------------ 
     90;  print,filename  
    7091; is parent a valid widget ? 
    71    if keyword_set(parentin) then BEGIN 
    72       parent = long(parentin) 
    73       parent = parent*widget_info(parent,/managed) 
    74    ENDIF 
    75    if n_elements(iodirectory) EQ 0 AND n_elements(iodir) NE 0 then iodirectory = iodir 
    76    filename = isafile(filename = filename, IODIRECTORY = iodirectory) 
     92  if keyword_set(parentin) then BEGIN 
     93    parent = long(parentin) 
     94    parent = parent*widget_info(parent, /managed) 
     95  ENDIF 
     96  filename = isafile(filename = filename, IODIRECTORY = iodir, _EXTRA = ex) 
    7797;------------------------------------------------------------ 
    7898; ouverture du fichier nom 
    7999;------------------------------------------------------------ 
    80    if size(filename, /type) NE 7 then $ 
     100  if size(filename, /type) NE 7 then $ 
    81101    return, report('read_ncdf cancelled') 
    82    cdfid=ncdf_open(filename) 
    83    contient=ncdf_inquire(cdfid) 
     102  IF !version.OS_FAMILY EQ 'unix' THEN spawn, '\file '+filename+' > /dev/null' 
     103  cdfid = ncdf_open(filename) 
     104  contient = ncdf_inquire(cdfid) 
    84105;------------------------------------------------------------ 
    85106; we check if the variable name exists in the file. 
    86107;------------------------------------------------------------ 
    87    if ncdf_varid(cdfid,name) EQ -1 then $ 
     108  if ncdf_varid(cdfid, name) EQ -1 then BEGIN 
     109    ncdf_close, cdfid 
    88110    return, report('variable '+name+' !C not found in the file '+filename) 
    89 ; 
    90    varcontient=ncdf_varinq(cdfid,name) 
     111  ENDIF 
     112  varcontient = ncdf_varinq(cdfid, name) 
     113;------------------------------------------------------------ 
     114; shall we redefine the grid parameters 
     115;------------------------------------------------------------ 
     116  if keyword_set(init) THEN initncdf, filename, _extra = ex 
    91117;------------------------------------------------------------ 
    92118; check the time axis and the debut and fin dates 
    93119;------------------------------------------------------------ 
    94    if n_elements(debut) EQ 0 then begin 
    95       debut = 0 
    96       timestep = 1 
    97    endif 
    98    if keyword_set(timestep) then begin 
    99       premiertps = debut 
    100       if n_elements(fin) NE 0 then derniertps = fin ELSE derniertps = premiertps 
    101       jpt = derniertps-premiertps+1 
    102       time = lindgen(jpt) 
    103    ENDIF ELSE BEGIN 
    104       date1 = juldate(debut, /vraidate) 
    105       if n_elements(fin) NE 0 then date2 = juldate(fin, /vraidate) ELSE date2 = date1 
    106       if keyword_set(parent) then BEGIN 
    107          widget_control, parent, get_uvalue = top_uvalue 
    108          filelist = extractatt(top_uvalue,'filelist') 
    109          currentfile = (where(filelist EQ filename))[0] 
    110          time = (*(extractatt(top_uvalue, 'fileparameters'))[currentfile]).time_counter 
    111          premiertps = where(time EQ date1) & premiertps = premiertps[0] 
    112          derniertps = where(time EQ date2) & derniertps = derniertps[0] 
    113       ENDIF ELSE BEGIN 
     120  if n_elements(debut) EQ 0 then begin 
     121    debut = 0 
     122    timestep = 1 
     123  endif 
     124  if keyword_set(timestep) then begin 
     125    firsttps = debut[0] 
     126    if n_elements(fin) NE 0 then lasttps = fin[0] ELSE lasttps = firsttps 
     127    jpt = lasttps-firsttps+1 
     128    time = lindgen(jpt) 
     129  ENDIF ELSE BEGIN 
     130    if keyword_set(parent) then BEGIN 
     131      widget_control, parent, get_uvalue = top_uvalue 
     132      filelist = extractatt(top_uvalue, 'filelist') 
     133      IF filelist[0] EQ 'many !' THEN filelist = filename 
     134      currentfile = (where(filelist EQ filename))[0] 
     135      time = (*(extractatt(top_uvalue, 'fileparameters'))[currentfile]).time_counter 
     136      date1 = date2jul(debut[0]) 
     137      if n_elements(fin) NE 0 then date2 = date2jul(fin[0]) ELSE date2 = date1 
     138      firsttps = where(time EQ date1) & firsttps = firsttps[0] 
     139      lasttps = where(time EQ date2) & lasttps = lasttps[0] 
     140    ENDIF ELSE BEGIN 
     141      IF keyword_set(timevar) THEN BEGIN  
     142        timeid = ncdf_varid(cdfid, timevar) 
     143        IF timeid EQ -1 THEN BEGIN 
     144          ncdf_close, cdfid 
     145          return, report('the file '+filename+' as no variable ' + timevar $ 
     146                         + '. !C Use the TIMESTEP keyword') 
     147        endif 
     148        timecontient = ncdf_varinq(cdfid, timeid) 
     149        contient.recdim = timecontient.dim[0] 
     150      ENDIF ELSE BEGIN  
    114151; we find the infinite dimension 
    115          timedim = contient.recdim 
    116          if timedim EQ -1 then $ 
    117           return, report('the file '+filename+' as no infinite dimension. !C Use the TIMESTEP keyword') 
     152        timedim = contient.recdim 
     153        if timedim EQ -1 then BEGIN 
     154          ncdf_close, cdfid 
     155          return, report('the file '+filename+' as no infinite dimension. !C Use TIMESTEP or TIMEVAR keyword') 
     156        endif 
    118157; we find the FIRST time axis       
    119          timeid = 0 
    120          repeat BEGIN           ; tant que l''on a pas trouve une variable qui n''a qu'' 
     158        timeid = 0 
     159        repeat BEGIN       ; tant que l''on a pas trouve une variable qui n''a qu'' 
    121160                                ; une dimension: la dimension infinie 
    122             timecontient=ncdf_varinq(cdfid,timeid) ; que contient la variable 
    123             timeid = timeid+1 
    124          endrep until (n_elements(timecontient.dim) EQ 1 AND timecontient.dim[0] EQ contient.recdim) $ 
     161          timecontient = ncdf_varinq(cdfid, timeid) ; que contient la variable 
     162          timeid = timeid+1 
     163        endrep until (n_elements(timecontient.dim) EQ 1 $ 
     164                      AND timecontient.dim[0] EQ contient.recdim) $ 
    125165          OR timeid EQ contient.nvars+1 
    126166; 
    127          if timeid EQ contient.nvars+1 then $ 
     167        if timeid EQ contient.nvars+1 then BEGIN 
     168          ncdf_close, cdfid 
    128169          return, report('the file '+filename+' as no time axis variable. !C Use the TIMESTEP keyword') 
     170        endif 
     171        timeid = timeid-1 
     172      ENDELSE  
    129173; we must found the time origin of the julian calendar used in the 
    130174; time axis.  
    131 ; does the attribut units exist for the variable time axis? 
    132          timeid = timeid-1 
    133          if timecontient.natts EQ 0 then $ 
    134           return, report('the variable '+timecontient.name+' has no attribut.!C Use the TIMESTEP keyword or add the attribut ''units'' to the variable') 
    135          attnames = strarr(timecontient.natts) 
    136          for attiq=0,timecontient.natts-1 do attnames[attiq]=ncdf_attname(cdfid,timeid,attiq) 
    137          if (where(attnames EQ 'units'))[0] EQ -1 then $ 
    138           return, report('Attribut ''units'' not found for the variable '+timecontient.name+' !C Use the TIMESTEP keyword') 
    139          ncdf_attget,cdfid,timeid,'units',value 
     175; does the attribut units an dcalendar exist for the variable time axis? 
     176      if timecontient.natts EQ 0 then BEGIN 
     177        ncdf_close, cdfid 
     178        return, report('the variable '+timecontient.name+' has no attribut.!C Use the TIMESTEP keyword or add the attribut ''units'' to the variable') 
     179      endif 
     180      attnames = strarr(timecontient.natts) 
     181      for attiq = 0, timecontient.natts-1 do attnames[attiq] = ncdf_attname(cdfid, timeid, attiq) 
     182      if (where(attnames EQ 'units'))[0] EQ -1 then BEGIN 
     183        ncdf_close, cdfid 
     184        return, report('Attribut ''units'' not found for the variable '+timecontient.name+' !C Use the TIMESTEP keyword') 
     185      ENDIF 
     186; 
     187; now we try to find the attribut called calendar... 
     188; the the attribute "calendar" exists? 
     189; If no, we suppose that the calendar is gregorian calendar 
     190; 
     191      if (where(attnames EQ 'calendar'))[0] NE -1 then BEGIN 
     192        ncdf_attget, cdfid, timeid, 'calendar', value 
     193        value = string(value) 
     194        CASE value OF 
     195          '360d':key_caltype = '360d' 
     196          'greg':IF n_elements(key_caltype) EQ 0 THEN key_caltype = 'greg' 
     197          ELSE:BEGIN 
     198;            notused = report('Unknown calendar: '+value+', we use greg calendar.')  
     199            key_caltype = 'greg' 
     200          END 
     201        ENDCASE 
     202      ENDIF ELSE BEGIN 
     203;        notused = report('Unknown calendar, we use '+key_caltype+' calendar.')  
     204        IF n_elements(key_caltype) EQ 0 THEN key_caltype = 'greg' 
     205      ENDELSE 
     206; 
     207; now we take acre of units attribut 
     208      ncdf_attget, cdfid, timeid, 'units', value 
    140209; 
    141210; time_counter:units = "seconds since 0001-01-01 00:00:00" ; 
    … …  
    146215; 
    147216; we decript the "units" attribut to find the time origin 
    148          value = strtrim(strcompress(string(value)), 2) 
    149          mots = str_sep(value, ' ') 
    150          unite = mots[0] 
    151          depart = str_sep(mots[2], '-') 
    152          ncdf_varget, cdfid, timeid, time 
    153          time = long(time) 
    154          case strlowcase(unite) of 
    155             'seconds':time = julday(depart[1], depart[2], depart[0])+time/(long(24)*3600) 
    156             'hours':time = julday(depart[1], depart[2], depart[0])+time/long(24) 
    157             'days':time = julday(depart[1], depart[2], depart[0])+time 
    158             'months':BEGIN  
    159                for t = 0, n_elements(time)-1  do begin 
    160                   time[t] = julday(depart[1]+time[t], depart[2], depart[0]) 
    161                endfor 
    162             END 
    163             'years':BEGIN 
    164                for t = 0, n_elements(time)-1  do begin 
    165                   time[t] = julday(depart[1], depart[2], depart[0]+time[t]) 
    166                endfor 
    167             END 
    168             ELSE:return, report('The "units" attribu of the time axis must be something like: !C "seconds since 0001-01-01 ..." !C "days since 1979-01-01 ..." !C "months since 1979-01-01 ..." !C "years since 1979-01-01 ..."') 
    169          ENDCASE 
    170          premiertps = where(time GE date1) & premiertps = premiertps[0] 
    171          if premiertps EQ -1 THEN $ 
    172           return, report('date 1: '+strtrim(date1, 1)+' is not found in the time axis.') 
    173          derniertps = where(time LE date2) 
    174          if derniertps[0] EQ -1 THEN $ 
    175           return, report('the time axis as no date before date 2: '+strtrim(date2, 1)) 
    176          derniertps = derniertps[n_elements(derniertps)-1] 
    177       ENDELSE 
    178       time = time[premiertps:derniertps] 
    179       jpt = derniertps-premiertps+1 
    180    ENDELSE 
     217      value = strtrim(strcompress(string(value)), 2) 
     218      mots = str_sep(value, ' ') 
     219      unite = mots[0] 
     220      depart = str_sep(mots[2], '-') 
     221      ncdf_varget, cdfid, timeid, time 
     222      time = double(time) 
     223      unite = strlowcase(unite) 
     224      IF strpos(unite, 's', strlen(unite)-1) NE -1 THEN unite = strmid(unite, 0, strlen(unite)-1) 
     225      IF strpos(unite, 'julian_') NE -1 THEN unite = strmid(unite, 7) 
     226      case unite of 
     227        'second':time = julday(depart[1], depart[2], depart[0])+time/(long(24)*3600) 
     228        'hour':time = julday(depart[1], depart[2], depart[0])+time/long(24) 
     229        'day':time = julday(depart[1], depart[2], depart[0])+time 
     230        'month':BEGIN  
     231          if total(fix(time) NE time) NE 0 then $ ; we switch to days with 30d/m 
     232            time = julday(depart[1], depart[2], depart[0])+round(time*30) $ 
     233            ELSE for t = 0, n_elements(time)-1 DO $ 
     234            time[t] = julday(depart[1]+time[t], depart[2], depart[0]) 
     235        END 
     236        'year':BEGIN 
     237          if total(fix(time) NE time) NE 0 then $ ; we switch to days with 365d/y 
     238            time = julday(depart[1], depart[2], depart[0])+round(time*365) $ 
     239            ELSE for t = 0, n_elements(time)-1 do $ 
     240            time[t] = julday(depart[1], depart[2], depart[0]+time[t]) 
     241        END 
     242        ELSE:BEGIN 
     243          ncdf_close, cdfid 
     244          return, report('The "units" attribu of the time axis must be something like: !C "seconds since 0001-01-01 ..." !C "days since 1979-01-01 ..." !C "months since 1979-01-01 ..." !C "years since 1979-01-01 ..."') 
     245        end 
     246      ENDCASE 
     247      date1 = date2jul(debut[0]) 
     248      if n_elements(fin) NE 0 then date2 = date2jul(fin[0]) ELSE date2 = date1 
     249      time = long(time) 
     250      firsttps = where(time GE date1) & firsttps = firsttps[0] 
     251      if firsttps EQ -1 THEN BEGIN 
     252        ncdf_close, cdfid 
     253        return, report('date 1: '+strtrim(jul2date(date1), 1)+' is not found in the time axis.') 
     254      endif 
     255      lasttps = where(time LE date2) 
     256      if lasttps[0] EQ -1 THEN BEGIN 
     257        ncdf_close, cdfid 
     258        return, report('the time axis as no date before date 2: '+strtrim(jul2date(date2), 1)) 
     259      endif 
     260      lasttps = lasttps[n_elements(lasttps)-1] 
     261      if lasttps LT firsttps then BEGIN 
     262        ncdf_close, cdfid 
     263        return, report('the time axis as no dates between date1 and  date 2: '+strtrim(jul2date(date1), 1)+' '+strtrim(jul2date(date2), 1)) 
     264      endif 
     265    ENDELSE 
     266    time = time[firsttps:lasttps] 
     267    jpt = lasttps-firsttps+1 
     268  ENDELSE 
    181269;------------------------------------------------------------ 
    182270; nom de la grille a laquelle se rapporte le champ 
    183271;------------------------------------------------------------ 
    184    vargrid = strupcase(strmid(filename, strlen(filename)-9)) 
    185    if vargrid EQ 'GRID.T.NC' OR vargrid EQ 'GRID.U.NC' $ 
    186     OR vargrid EQ 'GRID.V.NC' OR vargrid EQ 'GRID.F.NC' $ 
    187     OR vargrid eq 'GRID.W.NC' $ 
    188     OR vargrid EQ 'GRID_T.NC' OR vargrid EQ 'GRID_U.NC' $ 
    189     OR vargrid EQ 'GRID_V.NC' OR vargrid EQ 'GRID_F.NC' $ 
    190     OR vargrid eq 'GRID_W.NC' then $ 
    191     vargrid = strupcase(strmid(filename, strlen(filename)-4, 1)) ELSE vargrid='T' 
     272  IF keyword_set(grid) THEN vargrid = strupcase(grid) ELSE BEGIN 
     273    vargrid = 'T'               ; default definition 
     274    pattern = ['GRID.', 'GRID_', 'GRID', 'UPID_', '30ID_'] 
     275    gdtype = ['T', 'U', 'V', 'W', 'F'] 
     276    fnametest = strupcase(filename) 
     277    FOR i = 0, n_elements(pattern)-1 DO BEGIN 
     278      FOR j = 0, n_elements(gdtype)-1 DO BEGIN 
     279        substr = pattern[i]+gdtype[j] 
     280        pos = strpos(fnametest, substr) 
     281        IF pos NE -1 THEN $ 
     282          vargrid = strmid(fnametest, pos+strlen(substr)-1, 1) 
     283      ENDFOR 
     284    ENDFOR 
     285  ENDELSE 
     286;--------------------------------------------------------------- 
     287; call the init function ? 
     288;--------------------------------------------------------------- 
     289 
    192290;--------------------------------------------------------------- 
    193291; redefinition du domaine  
    194292;--------------------------------------------------------------- 
    195    if keyword_set(tout) then begin 
    196       nx=jpi 
    197       ny=jpj 
    198       nz=jpk 
    199       premierx = 0 
    200       premiery = 0 
    201       premierz = 0 
    202       dernierx = jpi-1 
    203       derniery = jpj-1 
    204       dernierz = jpk-1 
    205       case strupcase(vargrid) of 
    206          'T':mask = tmask 
    207          'U':mask = umask() 
    208          'V':mask = vmask() 
    209          'W':mask = tmask 
    210          'F':mask = fmask() 
    211       endcase 
    212    ENDIF ELSE BEGIN 
    213       if keyword_set(boite) then BEGIN  
    214          oldboite = [lon1, lon2, lat1, lat2, prof1, prof2] 
    215          Case 1 Of 
    216             N_Elements(Boite) Eq 1:bte=[lon1, lon2, lat1, lat2, 0.,boite[0]] 
    217             N_Elements(Boite) Eq 2:bte=[lon1, lon2, lat1, lat2, boite[0],boite[1]] 
    218             N_Elements(Boite) Eq 4:bte=[Boite, prof1, prof2] 
    219             N_Elements(Boite) Eq 5:bte=[Boite[0:3], 0, Boite[4]] 
    220             N_Elements(Boite) Eq 6:bte=Boite 
    221             Else: return, report('Mauvaise Definition de Boite') 
    222          ENDCASE 
    223          domdef, bte,GRILLE=['T', vargrid], _extra = ex 
    224       ENDIF 
    225       grille, mask, glam, gphi, gdep, nx, ny,nz,premierx,premiery,premierz,dernierx, derniery, dernierz 
    226       undefine, glam & undefine, gphi & ; on libere un peu de memoire! 
    227    ENDELSE 
     293  if keyword_set(tout) then begin 
     294    nx = jpi 
     295    ny = jpj 
     296    nz = jpk 
     297    firstx = 0 
     298    firsty = 0 
     299    firstz = 0 
     300    lastx = jpi-1 
     301    lasty = jpj-1 
     302    lastz = jpk-1 
     303    case strupcase(vargrid) of 
     304      'T':mask = tmask 
     305      'U':mask = umask() 
     306      'V':mask = vmask() 
     307      'W':mask = tmask 
     308      'F':mask = fmask() 
     309    endcase 
     310  ENDIF ELSE BEGIN 
     311    if keyword_set(boxzoom) then BEGIN  
     312      Case 1 Of 
     313        N_Elements(Boxzoom) Eq 1:bte = [lon1, lon2, lat1, lat2, 0., boxzoom[0]] 
     314        N_Elements(Boxzoom) Eq 2:bte = [lon1, lon2, lat1, lat2, boxzoom[0], boxzoom[1]] 
     315        N_Elements(Boxzoom) Eq 4:bte = [Boxzoom, vert1, vert2] 
     316        N_Elements(Boxzoom) Eq 5:bte = [Boxzoom[0:3], 0, Boxzoom[4]] 
     317        N_Elements(Boxzoom) Eq 6:bte = Boxzoom 
     318        Else: BEGIN 
     319          ncdf_close, cdfid 
     320          return, report('Wrong Definition of Boxzoom') 
     321        end 
     322      ENDCASE 
     323      savedbox = 1b 
     324      saveboxparam, 'boxparam4rdncdf.dat' 
     325      domdef, bte, GRIDTYPE = ['T', vargrid], _extra = ex 
     326    ENDIF 
     327    grille, mask, glam, gphi, gdep, nx, ny, nz, firstx, firsty, firstz, lastx, lasty, lastz 
     328    undefine, glam & undefine, gphi & ; on libere un peu de memoire! 
     329  ENDELSE 
    228330;--------------------------------------------------------------------- 
    229331; on initialise les ixmindta, iymindta au besoin 
    230332;--------------------------------------------------------------------- 
    231    if n_elements(jpidta) EQ 0 THEN jpidta = jpiglo 
    232    if n_elements(jpjdta) EQ 0 THEN jpjdta = jpjglo 
    233    if n_elements(jpkdta) EQ 0 THEN jpkdta = jpkglo 
    234    if n_elements(ixmindta) EQ 0 THEN ixmindta = 0 
    235    if n_elements(ixmaxdta) EQ 0 then ixmaxdta = jpidta-1 
    236    if ixmindta EQ -1 THEN ixmindta = 0 
    237    IF ixmaxdta EQ -1 then ixmaxdta = jpidta-1 
    238    if n_elements(iymindta) EQ 0 THEN iymindta = 0 
    239    IF n_elements(iymaxdta) EQ 0 then iymaxdta = jpjdta-1 
    240    if iymindta EQ -1 THEN iymindta = 0 
    241    IF iymaxdta EQ -1 then iymaxdta = jpjdta-1 
    242    if n_elements(izmindta) EQ 0 THEN izmindta = 0 
    243    IF n_elements(izmaxdta) EQ 0 then izmaxdta = jpkdta-1 
    244    if izmindta EQ -1 THEN izmindta = 0 
    245    IF izmaxdta EQ -1 then izmaxdta = jpkdta-1 
     333  if n_elements(jpidta) EQ 0 THEN jpidta = jpiglo 
     334  if n_elements(jpjdta) EQ 0 THEN jpjdta = jpjglo 
     335  if n_elements(jpkdta) EQ 0 THEN jpkdta = jpkglo 
     336  if n_elements(ixmindta) EQ 0 THEN ixmindta = 0 
     337  if n_elements(ixmaxdta) EQ 0 then ixmaxdta = jpidta-1 
     338  if ixmindta EQ -1 THEN ixmindta = 0 
     339  IF ixmaxdta EQ -1 then ixmaxdta = jpidta-1 
     340  if n_elements(iymindta) EQ 0 THEN iymindta = 0 
     341  IF n_elements(iymaxdta) EQ 0 then iymaxdta = jpjdta-1 
     342  if iymindta EQ -1 THEN iymindta = 0 
     343  IF iymaxdta EQ -1 then iymaxdta = jpjdta-1 
     344  if n_elements(izmindta) EQ 0 THEN izmindta = 0 
     345  IF n_elements(izmaxdta) EQ 0 then izmaxdta = jpkdta-1 
     346  if izmindta EQ -1 THEN izmindta = 0 
     347  IF izmaxdta EQ -1 then izmaxdta = jpkdta-1 
    246348;---------------------------------------------------------------- 
    247349; on va lire le fichier 
    248350;--------------------------------------------------------------- 
    249    if n_elements(key_stride) LE 2 then key_stride = [1, 1, 1] 
    250    key_stride = 1l > long(key_stride) 
    251    jpitotal = long(ixmaxmesh-ixminmesh+1) 
    252    key_shift = long(testvar(var = key_shift)) 
    253 ; 
    254    key = long(key_shift MOD jpitotal) 
    255    if key LT 0 then key = key+jpitotal 
    256    ixmin = ixminmesh-ixmindta 
    257 ; 
     351  if n_elements(key_stride) LE 2 then key_stride = [1, 1, 1] 
     352  key_stride = 1l > long(key_stride) 
     353  key_shift = long(testvar(var = key_shift)) 
     354; 
     355  IF n_elements(key_yreverse) EQ 0 THEN key_yreverse = 0 
     356  IF keyword_set(key_yreverse) THEN BEGIN 
     357    tmp = jpj-1-firsty 
     358    firsty = jpj-1-lasty 
     359    lasty = tmp 
     360  ENDIF 
     361; 
     362  IF keyword_set(fbase2tbase) THEN BEGIN 
     363    case strupcase(vargrid) of 
     364      'U':BEGIN 
     365        IF NOT keyword_set(key_periodic) THEN BEGIN 
     366          firstx = firstx+1 
     367          lastx = lastx+1 
     368        ENDIF 
     369      END 
     370      'V':BEGIN 
     371        firsty = firsty+1 
     372        lasty = lasty+1 
     373      END 
     374      'F':BEGIN 
     375        firsty = firsty+1 
     376        lasty = lasty+1 
     377        IF NOT keyword_set(key_periodic) THEN BEGIN 
     378          firstx = firstx+1 
     379          lastx = lastx+1 
     380        ENDIF 
     381      END 
     382      ELSE: 
     383    endcase 
     384  ENDIF 
     385; 
     386  IF keyword_set(fbase2tbase) AND keyword_set(key_periodic) $ 
     387    AND (strupcase(vargrid) EQ 'U' OR strupcase(vargrid) EQ 'F') THEN key_shift =  key_shift-1 
     388; 
     389;--------------------------------------------------------------------- 
     390;--------------------------------------------------------------------- 
    258391@read_ncdf_varget 
    259392;--------------------------------------------------------------------- 
     393;--------------------------------------------------------------------- 
     394; 
     395  IF keyword_set(fbase2tbase) AND keyword_set(key_periodic) $ 
     396    AND (strupcase(vargrid) EQ 'U' OR strupcase(vargrid) EQ 'F') THEN key_shift =  key_shift+1 
     397;--------------------------------------------------------------------- 
    260398; on definit les variables globales rattachees a la variable 
    261399;--------------------------------------------------------------------- 
    262400; varname 
    263    varname = name 
     401  varname = name 
    264402; varunit 
    265    if varcontient.natts NE 0 then begin 
    266       attnames = strarr(varcontient.natts) 
    267       for attiq=0,varcontient.natts-1 do attnames[attiq]=strlowcase(ncdf_attname(cdfid,name,attiq)) 
    268 ; 
    269       if (where(attnames EQ 'units'))[0] EQ -1 then value = '' $ 
    270       ELSE ncdf_attget,cdfid,name,'units',value 
    271       varunit = strtrim(string(value), 2) 
    272 ; 
    273       if (where(attnames EQ 'add_offset'))[0] EQ -1 then add_offset = 0. $ 
    274       ELSE ncdf_attget,cdfid,name,'add_offset',add_offset 
    275 ; 
    276       if (where(attnames EQ 'scale_factor'))[0] EQ -1 then scale_factor = 1. $ 
    277       ELSE ncdf_attget,cdfid,name,'scale_factor',scale_factor 
    278 ; 
    279       if (where(attnames EQ 'missing_value'))[0] EQ -1 then missing_value = 'no' $ 
    280       ELSE ncdf_attget,cdfid,name,'missing_value',missing_value 
    281 ; 
    282    ENDIF ELSE BEGIN  
    283       varunit = '' 
    284       add_offset = 0. 
    285       scale_factor = 1. 
    286       missing_value = 'no' 
    287    ENDELSE 
     403  if varcontient.natts NE 0 then begin 
     404    attnames = strarr(varcontient.natts) 
     405    for attiq = 0, varcontient.natts-1 do attnames[attiq] = ncdf_attname(cdfid, name, attiq) 
     406    lowattnames = strlowcase(attnames) 
     407; 
     408    found = (where(lowattnames EQ 'units'))[0] 
     409    IF found NE -1 then ncdf_attget, cdfid, name, attnames[found], value ELSE value = '' 
     410    varunit = strtrim(string(value), 2) 
     411; 
     412    found = (where(lowattnames EQ 'add_offset'))[0] 
     413    if found NE -1 then ncdf_attget, cdfid, name, attnames[found], add_offset ELSE add_offset = 0. 
     414; 
     415    found = (where(lowattnames EQ 'scale_factor'))[0] 
     416    if found NE -1 then ncdf_attget, cdfid, name, attnames[found], scale_factor ELSE scale_factor = 1. 
     417; 
     418    missing_value = 'no' 
     419    found = (where(lowattnames EQ '_fillvalue'))[0] 
     420    if found NE -1 then ncdf_attget, cdfid, name, attnames[found], missing_value 
     421    found = (where(lowattnames EQ 'missing_value'))[0] 
     422    if found NE -1 then ncdf_attget, cdfid, name, attnames[found], missing_value 
     423; 
     424  ENDIF ELSE BEGIN  
     425    varunit = '' 
     426    add_offset = 0. 
     427    scale_factor = 1. 
     428    missing_value = 'no' 
     429  ENDELSE 
    288430; vardate 
    289431; on construit une belle date lisible en fonction du langage specifie !!! 
    290    nothing = juldate(debut,/vrai) ; bidouille pour recuperer year, month, day 
    291    if keyword_set(language) then BEGIN  
    292       if language eq 'gb' THEN $  
    293        vardate = string(format='(C(CMoA))',31*(month-1))+' '+strtrim(day, 1)+' '+strtrim(year, 1) $ 
    294       ELSE vardate = strtrim(day, 1)+' '+string(format='(C(CMoA))',31*(month-1))+' '+strtrim(year, 1) 
    295    ENDIF ELSE vardate = strtrim(day, 1)+' '+string(format='(C(CMoA))',31*(month-1))+' '+strtrim(year, 1) 
     432  year = long(debut[0])/10000 
     433  month = (long(debut[0])/100) MOD 100 
     434  day = (long(debut[0]) MOD 100) 
     435  vardate = string(format = '(C(CMoA))', 31*(month-1))+' '+strtrim(day, 1)+', '+strtrim(year, 1) 
    296436; varexp ??? 
     437 
     438; we apply reverse 
     439  if keyword_set(key_yreverse) then res = reverse(temporary(res),  2) 
     440  if keyword_set(key_zreverse) AND (size(res))[0] EQ 3 AND jpt EQ 1 then res = reverse(temporary(res), 3) 
     441  if keyword_set(key_zreverse) AND (size(res))[0] EQ 4 THEN res = reverse(temporary(res), 3) 
     442 
    297443; on applique la valeur valmask sur les points terre 
    298    if NOT keyword_set(cont_nofill) then begin 
    299       valmask = 1e20 
    300       case 1 of 
    301          varcontient.ndims eq 2:BEGIN ;xy array 
    302             mask = mask[*, *, 0] 
    303             earth = where(mask EQ 0) 
    304             if earth[0] NE -1 then res[earth] = 1e20 
    305          END 
    306          varcontient.ndims eq 3 AND (where(varcontient.dim EQ contient.recdim))[0] EQ -1:BEGIN ;xyz array 
    307             earth = where(mask EQ 0) 
    308             if earth[0] NE -1 then res[earth] = 1e20 
    309          END 
    310          varcontient.ndims eq 3 AND (where(varcontient.dim EQ contient.recdim))[0] NE -1:BEGIN ;xyt array 
    311             mask = mask[*, *, 0] 
    312             earth = where(mask EQ 0) 
    313             if earth[0] NE -1 then BEGIN 
    314                earth = earth#replicate(1, jpt)+replicate(nx*ny, n_elements(earth))#lindgen(jpt) 
    315                res[earth] = 1e20 
    316             END 
    317          END 
    318          varcontient.ndims eq 4:BEGIN ;xyzt array 
    319             earth = where(mask EQ 0) 
    320             if earth[0] NE -1 then BEGIN 
    321                earth = earth#replicate(1, jpt)+replicate(nx*ny*nz, n_elements(earth))#lindgen(jpt) 
    322                res[earth] = 1e20 
    323             END 
    324          END 
    325       endcase 
    326    endif 
     444  if NOT keyword_set(cont_nofill) then begin 
     445    valmask = 1e20 
     446    case 1 of 
     447      varcontient.ndims eq 2:BEGIN ;xy array 
     448        mask = mask[*, *, 0] 
     449        earth = where(mask EQ 0) 
     450      END 
     451      varcontient.ndims eq 3 AND (where(varcontient.dim EQ contient.recdim))[0] EQ -1:BEGIN ;xyz array 
     452        earth = where(mask EQ 0) 
     453      END 
     454      varcontient.ndims eq 3 AND (where(varcontient.dim EQ contient.recdim))[0] NE -1:BEGIN ;xyt array 
     455        mask = mask[*, *, 0] 
     456        earth = where(mask EQ 0) 
     457        if earth[0] NE -1 then BEGIN 
     458          earth = earth#replicate(1, jpt)+replicate(nx*ny, n_elements(earth))#lindgen(jpt) 
     459        END 
     460      END 
     461      varcontient.ndims eq 4:BEGIN ;xyzt array 
     462        earth = where(mask EQ 0) 
     463        if earth[0] NE -1 then BEGIN 
     464          earth = earth#replicate(1, jpt)+replicate(nx*ny*nz, n_elements(earth))#lindgen(jpt) 
     465        END 
     466      END 
     467    endcase 
     468  ENDIF ELSE earth = -1 
     469; we look for  missing_value 
     470  IF size(missing_value, /type) NE 7 then BEGIN 
     471    IF size(missing_value, /type) EQ 1 THEN BEGIN  
     472      IF isnumber(string(missing_value), tmp) EQ 1 THEN missing_value = tmp 
     473    ENDIF  
     474;    if missing_value NE valmask then begin 
     475    if abs(missing_value) LT 1e6 then missing = where(res EQ missing_value) $ 
     476    ELSE missing = where(abs(res) gt abs(missing_value)/10.) 
     477;    ENDIF ELSE missing = -1 
     478  ENDIF ELSE missing = -1 
    327479; on applique les add_offset, scale_factor et missing_value 
    328    if size(missing_value, /type) NE 7 then BEGIN 
    329       if missing_value NE valmask then begin 
    330          if abs(missing_value) LT 1e6 then missing = where(res EQ missing_value) $ 
    331          ELSE missing = where(abs(res) gt abs(missing_value)/10.) 
    332          if missing[0] NE -1 then res[temporary(missing)] = !values.f_nan 
    333       endif 
    334    ENDIF 
    335    if scale_factor NE 1 then res = res*scale_factor 
    336    if add_offset NE 0 then res = res+add_offset 
    337 ; 
    338    if keyword_set(key_reverse) then res = reverse(temporary(res),  2) 
    339 ;--------------------------------------------------------------------- 
    340    ncdf_close,cdfid 
    341 ;--------------------------------------------------------------------- 
    342    if keyword_set(oldboite) then domdef, oldboite,GRILLE=['T', vargrid] 
    343    if keyword_set(nostruct) then return, res $ 
    344    ELSE return, {tab:res, grille:vargrid, unite:varunit, experience:varexp, nom:varname} 
    345  
    346 end 
    347  
    348  
    349  
    350  
    351  
    352  
     480  if scale_factor NE 1 then res = temporary(res)*scale_factor 
     481  if add_offset NE 0 then res = temporary(res)+add_offset 
     482  if missing[0] NE -1 then res[temporary(missing)] = !values.f_nan 
     483  if earth[0] NE -1 then res[temporary(earth)] = 1e20 
     484;--------------------------------------------------------------------- 
     485  ncdf_close, cdfid 
     486;--------------------------------------------------------------------- 
     487  if keyword_set(savedbox) THEN restoreboxparam, 'boxparam4rdncdf.dat' 
     488  if keyword_set(nostruct) then return, res $ 
     489  ELSE BEGIN  
     490    IF keyword_set(key_forgetold) THEN BEGIN 
     491      return, {arr:res, grid:vargrid, unit:varunit, experiment:varexp, name:varname}  
     492    ENDIF ELSE BEGIN  
     493      return, {tab:res, grille:vargrid, unite:varunit, experience:varexp, nom:varname} 
     494    ENDELSE  
     495  ENDELSE  
     496END 
     497 
     498 
     499 
     500 
     501 
     502 

  • trunk/ToBeReviewed/LECTURE/read_ncdf_varget.pro

    r42 r44  
    44;  attention, qd il y a un shift il faut faire attention a ce que 
    55;  l''on fait et en particulier arriver a positionner la boite 
    6 ;  specifiee par premierx:dernierx par rapport au tableau contenu dans 
     6;  specifiee par firstx:lastx par rapport au tableau contenu dans 
    77;  le fichier NetCdf. 
    88;  pour s''y reperer voici un petit dessin 
    99; 
    10 ;  dans le repere 0:jpi, avec le zoom definit par premierx:dernierx 
     10;  dans le repere 0:jpi, avec le zoom definit par firstx:lastx 
    1111; 
    1212;        0      j      key-1  key          i              jpi-1 
    … …  
    2424;   |,,,,|--------|--------------|      |........|..........|,,,,,,,,,,,| 
    2525; 
    26 ;  apres il suffit de remplacer i et j par premierx ou dernierx qd on 
    27 ;  ne fait pas de stride et par premierx*key_stride[0] ou 
    28 ;  dernierx*key_stride[0] qd on fait un stride (dans ce cas il faut 
    29 ;  aussi remplacer jpi par jpitotal)! 
    30 ; 
     26;  apres il suffit de remplacer i et j par firstx ou lastx qd on 
     27;  ne fait pas de stride. 
     28; 
     29;  When using stride, things are getting more complicated.... 
     30; 
     31;  1) when must replace firstx by firstx*key_stride[0] and 
     32;  lastx by lastx*key_stride[0] 
     33;  2) we must use jpitotal instead of jpi 
     34;  3) finding the position of the offset in the right side part is 
     35;  tricky. look at this illustration of stride = 4 ... 
     36;        0      j      key-1 key          i              jpi-1 
     37;        +...+..|+...+...+  |--+---+---+|--+---+---+---+-| 
     38; 
     39;        0       i-key        jpi-key-1 jpi-key jpi-key+j jpi-1 
     40;        |--+---+---+|--+---+---+---+-|  +...+..|+...+...+ 
     41;  in the ........ part, it is easy the first point is number 0 but in 
     42;  the --------- part the first point is number 3... Ihe position of 
     43;  the first point in the is given by: 
     44;         (key_stride[0]-1)-((key-1) MOD key_stride[0]) 
     45;  4) last point...by default, the number of element read by IDL when 
     46;  using stride is given by n_elements/stride. However, we must read: 
     47;     (n_elements+stride-1)/stride or (n_elements-1)/key_stride+1 
     48;  for example: for n_elements between 10 and 12, with a stride of 3, 
     49;  we must read 4 points instead of 3...:  +--+--+--+-- 
     50;  This problem as an easy solution by using the keyword count with 
     51;  the appropiate value! 
     52; 
     53   ixmin = ixminmesh-ixmindta 
     54   iymin = iyminmesh-iymindta 
     55   izmin = izminmesh-izmindta 
     56   jpitotal = long(ixmaxmesh-ixminmesh+1) 
     57   key = long(key_shift MOD jpitotal) 
     58   if key LT 0 then key = key+jpitotal 
    3159; 
    3260case 1 OF 
    … …  
    4169; on ne shift pas et il n''y a pas de stride 
    4270;,,,,,,,,,,,,,,,,,,,,,,,, 
    43             ncdf_varget,cdfid,name,res,offset=[premierx+ixminmesh-ixmindta $ 
    44                                                ,premiery+iyminmesh-iymindta] $ 
     71            ncdf_varget,cdfid,name,res,offset=[firstx+ixmin,firsty+iymin] $ 
    4572             ,count=[nx,ny], _extra = ex 
    4673         END 
    … …  
    4976; on ne shift pas mais il y a un stride 
    5077;,,,,,,,,,,,,,,,,,,,,,,,, 
    51             ncdf_varget,cdfid,name,res,offset=[premierx*key_stride[0]+ixminmesh-ixmindta $ 
    52                                                ,premiery*key_stride[1]+iyminmesh-iymindta] $ 
     78            ncdf_varget,cdfid,name,res,offset=[firstx*key_stride[0]+ixmin $ 
     79                                               ,firsty*key_stride[1]+iymin] $ 
    5380             ,count=[nx,ny], stride = key_stride[0:1], _extra = ex 
    5481         END 
    … …  
    5885;,,,,,,,,,,,,,,,,,,,,,,,, 
    5986            case 1 of 
    60                premierx GE key:BEGIN ; on peut tout couper d''un coup 
    61                   ncdf_varget,cdfid,name,res,offset=[ixmin-key+premierx $ 
    62                                                      ,premiery+iyminmesh-iymindta] $ 
     87; --------- part, we can directly extract the array in one piece 
     88               firstx GE key:BEGIN 
     89                  ncdf_varget,cdfid,name,res,offset=[ixmin-key+firstx $ 
     90                                                     ,firsty+iymin] $ 
    6391                   ,count=[nx,ny], _extra = ex 
    6492               END 
    65                dernierx LE key-1:BEGIN ; on peut tout couper d''un coup 
    66                   ncdf_varget,cdfid,name,res,offset=[jpi-key+ixmin+premierx $ 
    67                                                      ,premiery+iyminmesh-iymindta] $ 
     93; ......... part, we can directly extract the array in one piece 
     94               lastx LE key-1:BEGIN 
     95                  ncdf_varget,cdfid,name,res,offset=[jpi-key+ixmin+firstx $ 
     96                                                     ,firsty+iymin] $ 
    6897                   ,count=[nx,ny], _extra = ex 
    6998               END 
    70                ELSE:BEGIN       ; Le tableau est separe en 2 morceaux... 
    71                   ncdf_varget,cdfid,name,tab1,offset=[jpi-key+ixmin+premierx $ 
    72                                                       ,premiery+iyminmesh-iymindta] $ 
    73                    ,count=[key-1-premierx+1,ny], _extra = ex 
    74                   ncdf_varget,cdfid,name,tab2,offset=[ixmin $ 
    75                                                       ,premiery+iyminmesh-iymindta] $ 
    76                    ,count=[dernierx-key+1,ny], _extra = ex 
     99               ELSE:BEGIN  ; we have to extract the array in 2 pieces...   
     100; ......... part, first part... 
     101                  ncdf_varget,cdfid,name,tab1,offset=[jpi-key+ixmin+firstx $ 
     102                                                      ,firsty+iymin] $ 
     103                   ,count=[key-1-firstx+1,ny], _extra = ex 
     104; --------- part, second part... 
     105                  ncdf_varget,cdfid,name,tab2,offset=[ixmin,firsty+iymin] $ 
     106                   ,count=[lastx-key+1,ny], _extra = ex 
    77107                  res = [temporary(tab1), temporary(tab2)] 
    78108               END 
    … …  
    84114;,,,,,,,,,,,,,,,,,,,,,,,, 
    85115            case 1 OF           ; case sur la facon de fire le champ 
    86                premierx GE ceil(1.*key/key_stride[0]):BEGIN ; on peut tout coupe d''un coup 
    87                   ncdf_varget,cdfid,name,res,offset=[ixmin-key+premierx*key_stride[0] $ 
    88                                                      ,premiery*key_stride[1]+iyminmesh-iymindta] $ 
     116; --------- part, we can directly extract the array in one piece 
     117               firstx*key_stride[0] GE key:BEGIN  
     118                  ncdf_varget,cdfid,name,res,offset=[ixmin $ 
     119                                                     +(key_stride[0]-1)-((key-1) MOD key_stride[0]) $ 
     120                                                     +(firstx-((key-1)/key_stride[0]+1))*key_stride[0] $ 
     121                                                     ,firsty*key_stride[1]+iymin] $ 
    89122                   ,count=[nx,ny], stride = key_stride[0:1], _extra = ex 
    90                END 
    91                dernierx LE (key-1)/key_stride[0]:BEGIN ; on peut tout coupe d''un coup 
    92                   ncdf_varget,cdfid,name,res,offset=[jpitotal-key+ixmin+premierx*key_stride[0] $ 
    93                                                      ,premiery*key_stride[1]+iyminmesh-iymindta] $ 
     123                END 
     124; ......... part, we can directly extract the array in one piece 
     125               lastx*key_stride[0] LE key-1:BEGIN  
     126                  ncdf_varget,cdfid,name,res,offset=[jpitotal-key+ixmin+firstx*key_stride[0] $ 
     127                                                     ,firsty*key_stride[1]+iymin] $ 
    94128                   ,count=[nx,ny], stride = key_stride[0:1], _extra = ex 
    95129               END 
    96                ELSE:BEGIN       ; le tableau est separe en 2 morceaux... 
    97                   ncdf_varget,cdfid,name,tab1,offset=[jpitotal-key+ixmin+premierx*key_stride[0] $ 
    98                                                       ,premiery*key_stride[1]+iyminmesh-iymindta] $ 
    99                    ,count=[(key-1-premierx*key_stride[0]+1)/key_stride[0],ny] $ 
    100                    , stride = key_stride[0:1], _extra = ex 
    101 ; grosse blague de IDL qui sur un vecteur 0 1 2 3 4 5 6 7 8 
    102 ; extrait pour offset = 0 ou 1 tjs 5 elements pour un stride de 2 par 
    103 ; ex, soit le nombre d''elements/le stride. Il manque donc une colonne 
    104 ; quand nombre d''elements-1 est un multiple de stride...CQFD 
    105 ; il en manque un bout? 
    106                   IF ( (key-1-premierx*key_stride[0]) MOD key_stride[0] ) EQ 0 then begin 
    107                      ncdf_varget,cdfid,name,tab1bis,offset=[jpitotal-1+ixmin $ 
    108                                                             ,premiery*key_stride[1]+iyminmesh-iymindta] $ 
    109                       ,count=[1,ny] $ 
    110                       , stride = [1, key_stride[1]], _extra = ex 
    111                      tab1bis = reform(tab1bis, 1, ny, /over) 
    112                      tab1 = [temporary(tab1), temporary(tab1bis)] 
    113                   ENDIF 
    114 ; 2eme bout... 
    115                   ncdf_varget,cdfid,name,tab2,offset=[ceil(1.*key/key_stride[0])*key_stride[0]-key+ixmin $ 
    116                                                       ,premiery*key_stride[1]+iyminmesh-iymindta] $ 
    117                    ,count=[(dernierx*key_stride[0]-key+1)/key_stride[0],ny] $ 
    118                    , stride = key_stride[0:1], _extra = ex 
    119 ; on recolle le tout 
     130               ELSE:BEGIN ; we have to extract the array in 2 pieces... 
     131; ......... part, first part... 
     132                  ncdf_varget,cdfid,name,tab1,offset=[jpitotal-key+ixmin+firstx*key_stride[0] $ 
     133                                                      ,firsty*key_stride[1]+iymin] $ 
     134                   ,count=[(key-firstx*key_stride[0]-1)/key_stride[0]+1,ny] $ 
     135                    , stride = key_stride[0:1], _extra = ex 
     136; --------- part, second part... 
     137                  ncdf_varget,cdfid,name,tab2,offset=[ixmin $ 
     138                                                      +(key_stride[0]-1)-((key-1) MOD key_stride[0]) $ 
     139                                                      , firsty*key_stride[1]+iymin] $ 
     140                   ,count=[nx-(key-firstx*key_stride[0]-1)/key_stride[0]-1,ny] $ 
     141                    , stride = key_stride[0:1], _extra = ex 
     142; in one unique array... 
    120143                  res = [temporary(tab1), temporary(tab2)] 
    121144               END 
    … …  
    134157; on ne shift pas et il n''y a pas de stride 
    135158;,,,,,,,,,,,,,,,,,,,,,,,, 
    136             ncdf_varget,cdfid,name,res,offset=[premierx+ixminmesh-ixmindta $ 
    137                                                ,premiery+iyminmesh-iymindta $ 
    138                                                ,premierz+izminmesh-izmindta] $ 
     159            ncdf_varget,cdfid,name,res,offset=[firstx+ixmin,firsty+iymin,firstz+izmin] $ 
    139160             ,count=[nx,ny,nz], _extra = ex 
    140161         END 
    … …  
    143164; on ne shift pas mais il y a un stride 
    144165;,,,,,,,,,,,,,,,,,,,,,,,, 
    145             ncdf_varget,cdfid,name,res,offset=[premierx*key_stride[0]+ixminmesh-ixmindta $ 
    146                                                ,premiery*key_stride[1]+iyminmesh-iymindta $ 
    147                                                ,premierz*key_stride[2]+izminmesh-izmindta] $ 
    148              ,count=[nx,ny,nz], stride = key_stride[0:2], _extra = ex 
     166            ncdf_varget,cdfid,name,res,offset=[firstx*key_stride[0]+ixmin $ 
     167                                               ,firsty*key_stride[1]+iymin $ 
     168                                               ,firstz*key_stride[2]+izmin] $ 
     169             ,count=[nx,ny,nz], stride = key_stride, _extra = ex 
    149170         END 
    150171;,,,,,,,,,,,,,,,,,,,,,,,, 
    … …  
    153174;,,,,,,,,,,,,,,,,,,,,,,,, 
    154175            case 1 of 
    155                premierx GE key:BEGIN ; on peut tout couper d''un coup 
    156                   ncdf_varget,cdfid,name,res,offset=[ixmin-key+premierx $ 
    157                                                      ,premiery+iyminmesh-iymindta $ 
    158                                                      ,premierz+izminmesh-izmindta] $ 
     176               firstx GE key:BEGIN ; on peut tout couper d''un coup 
     177                  ncdf_varget,cdfid,name,res,offset=[ixmin-key+firstx $ 
     178                                                     ,firsty+iymin $ 
     179                                                     ,firstz+izmin] $ 
    159180                   ,count=[nx,ny,nz], _extra = ex 
    160181               END 
    161                dernierx LE key-1:BEGIN ; on peut tout couper d''un coup 
    162                   ncdf_varget,cdfid,name,res,offset=[jpi-key+ixmin+premierx $ 
    163                                                      ,premiery+iyminmesh-iymindta $ 
    164                                                      ,premierz+izminmesh-izmindta] $ 
     182               lastx LE key-1:BEGIN ; on peut tout couper d''un coup 
     183                  ncdf_varget,cdfid,name,res,offset=[jpi-key+ixmin+firstx $ 
     184                                                     ,firsty+iymin $ 
     185                                                     ,firstz+izmin] $ 
    165186                   ,count=[nx,ny,nz], _extra = ex 
    166187               END 
    167188               ELSE:BEGIN       ; Le tableau est separe en 2 morceaux... 
    168                   ncdf_varget,cdfid,name,tab1,offset=[jpi-key+ixmin+premierx $ 
    169                                                       ,premiery+iyminmesh-iymindta $ 
    170                                                       ,premierz+izminmesh-izmindta] $ 
    171                    ,count=[key-1-premierx+1,ny,nz], _extra = ex 
    172                   ncdf_varget,cdfid,name,tab2,offset=[ixmin $ 
    173                                                       ,premiery+iyminmesh-iymindta $ 
    174                                                       ,premierz+izminmesh-izmindta] $ 
    175                    ,count=[dernierx-key+1,ny,nz], _extra = ex 
     189                  ncdf_varget,cdfid,name,tab1,offset=[jpi-key+ixmin+firstx $ 
     190                                                      ,firsty+iymin $ 
     191                                                      ,firstz+izmin] $ 
     192                   ,count=[key-1-firstx+1,ny,nz], _extra = ex 
     193                  ncdf_varget,cdfid,name,tab2,offset=[ixmin $ 
     194                                                      ,firsty+iymin $ 
     195                                                      ,firstz+izmin] $ 
     196                   ,count=[lastx-key+1,ny,nz], _extra = ex 
    176197                  res = [temporary(tab1), temporary(tab2)] 
    177198               END 
    … …  
    183204;,,,,,,,,,,,,,,,,,,,,,,,, 
    184205            case 1 OF           ; case sur la facon de fire le champ 
    185                premierx GE ceil(1.*key/key_stride[0]):BEGIN ; on peut tout coupe d''un coup 
    186                   ncdf_varget,cdfid,name,res,offset=[ixmin-key+premierx*key_stride[0] $ 
    187                                                      ,premiery*key_stride[1]+iyminmesh-iymindta $ 
    188                                                      ,premierz*key_stride[2]+izminmesh-izmindta] $ 
    189                    ,count=[nx,ny,nz], stride = key_stride[0:2], _extra = ex 
    190                END 
    191                dernierx LE (key-1)/key_stride[0]:BEGIN ; on peut tout coupe d''un coup 
    192                   ncdf_varget,cdfid,name,res,offset=[jpitotal-key+ixmin+premierx*key_stride[0] $ 
    193                                                      ,premiery*key_stride[1]+iyminmesh-iymindta $ 
    194                                                      ,premierz*key_stride[2]+izminmesh-izmindta] $ 
    195                    ,count=[nx,ny,nz], stride = key_stride[0:2], _extra = ex 
     206               firstx*key_stride[0] GE key:BEGIN ; on peut tout couper d''un coup 
     207                  ncdf_varget,cdfid,name,res,offset=[ixmin $ 
     208                                                     +(key_stride[0]-1)-((key-1) MOD key_stride[0]) $ 
     209                                                     +(firstx-((key-1)/key_stride[0]+1))*key_stride[0] $ 
     210                                                     ,firsty*key_stride[1]+iymin $ 
     211                                                     ,firstz*key_stride[2]+izmin] $ 
     212                   ,count=[nx,ny,nz], stride = key_stride, _extra = ex 
     213               END 
     214               lastx*key_stride[0] LE key-1:BEGIN ; on peut tout couper d''un coup 
     215                  ncdf_varget,cdfid,name,res,offset=[jpitotal-key+ixmin+firstx*key_stride[0] $ 
     216                                                     ,firsty*key_stride[1]+iymin $ 
     217                                                     ,firstz*key_stride[2]+izmin] $ 
     218                   ,count=[nx,ny,nz], stride = key_stride, _extra = ex 
    196219               END 
    197220               ELSE:BEGIN       ; le tableau est separe en 2 morceaux... 
    198                   ncdf_varget,cdfid,name,tab1,offset=[jpitotal-key+ixmin+premierx*key_stride[0] $ 
    199                                                       ,premiery*key_stride[1]+iyminmesh-iymindta $ 
    200                                                       ,premierz*key_stride[2]+izminmesh-izmindta] $ 
    201                    ,count=[(key-1-premierx*key_stride[0]+1)/key_stride[0],ny,nz] $ 
    202                    , stride = key_stride[0:2], _extra = ex 
    203 ; il en manque un bout? 
    204                   IF ( (key-1-premierx*key_stride[0]) MOD key_stride[0] ) EQ 0 then begin 
    205                      ncdf_varget,cdfid,name,tab1bis,offset=[jpitotal-1+ixmin $ 
    206                                                             ,premiery*key_stride[1]+iyminmesh-iymindta $ 
    207                                                             ,premierz*key_stride[2]+izminmesh-izmindta] $ 
    208                       ,count=[1,ny,nz] $ 
    209                       , stride = [1, key_stride[1:2]], _extra = ex 
    210                      tab1bis = reform(tab1bis, 1, ny, nz, /over) 
    211                      tab1 = [temporary(tab1), temporary(tab1bis)] 
    212                   ENDIF 
     221                  ncdf_varget,cdfid,name,tab1,offset=[jpitotal-key+ixmin+firstx*key_stride[0] $ 
     222                                                      ,firsty*key_stride[1]+iymin $ 
     223                                                      ,firstz*key_stride[2]+izmin] $ 
     224                   ,count=[(key-firstx*key_stride[0]-1)/key_stride[0]+1,ny,nz] $ 
     225                    , stride = key_stride, _extra = ex 
    213226; 2eme bout... 
    214                   ncdf_varget,cdfid,name,tab2,offset=[ceil(1.*key/key_stride[0])*key_stride[0]-key+ixmin $ 
    215                                                       ,premiery*key_stride[1]+iyminmesh-iymindta $ 
    216                                                       ,premierz*key_stride[2]+izminmesh-izmindta] $ 
    217                    ,count=[(dernierx*key_stride[0]-key+1)/key_stride[0],ny,nz] $ 
    218                    , stride = key_stride[0:2], _extra = ex 
     227                  ncdf_varget,cdfid,name,tab2,offset=[ixmin $ 
     228                                                      +(key_stride[0]-1)-((key-1) MOD key_stride[0]) $ 
     229                                                      ,firsty*key_stride[1]+iymin $ 
     230                                                      ,firstz*key_stride[2]+izmin] $ 
     231                   ,count=[nx-(key-firstx*key_stride[0]-1)/key_stride[0]-1,ny,nz] $ 
     232                    , stride = key_stride, _extra = ex 
    219233; on recolle le tout 
    220234                  res = [temporary(tab1), temporary(tab2)] 
    … …  
    234248; on ne shift pas et il n''y a pas de stride 
    235249;,,,,,,,,,,,,,,,,,,,,,,,, 
    236             ncdf_varget,cdfid,name,res,offset=[premierx+ixminmesh-ixmindta $ 
    237                                                ,premiery+iyminmesh-iymindta $ 
    238                                                ,premiertps] $ 
     250            ncdf_varget,cdfid,name,res,offset=[firstx+ixmin $ 
     251                                               ,firsty+iymin $ 
     252                                               ,firsttps] $ 
    239253             ,count=[nx,ny,jpt], _extra = ex 
    240254         END 
    … …  
    243257; on ne shift pas mais il y a un stride 
    244258;,,,,,,,,,,,,,,,,,,,,,,,, 
    245             ncdf_varget,cdfid,name,res,offset=[premierx*key_stride[0]+ixminmesh-ixmindta $ 
    246                                                ,premiery*key_stride[1]+iyminmesh-iymindta $ 
    247                                                ,premiertps] $ 
     259            ncdf_varget,cdfid,name,res,offset=[firstx*key_stride[0]+ixmin $ 
     260                                               ,firsty*key_stride[1]+iymin $ 
     261                                               ,firsttps] $ 
    248262             ,count=[nx,ny,jpt], stride = [key_stride[0:1], 1], _extra = ex 
    249263         END 
    … …  
    253267;,,,,,,,,,,,,,,,,,,,,,,,, 
    254268            case 1 of 
    255                premierx GE key:BEGIN ; on peut tout couper d''un coup 
    256                   ncdf_varget,cdfid,name,res,offset=[ixmin-key+premierx $ 
    257                                                      ,premiery+iyminmesh-iymindta $ 
    258                                                      ,premiertps] $ 
     269; --------- part, we can directly extract the array in one piece 
     270               firstx GE key:BEGIN ; on peut tout couper d''un coup 
     271                  ncdf_varget,cdfid,name,res,offset=[ixmin-key+firstx $ 
     272                                                     ,firsty+iymin $ 
     273                                                     ,firsttps] $ 
    259274                   ,count=[nx,ny,jpt], _extra = ex 
    260275               END 
    261                dernierx LE key-1:BEGIN ; on peut tout couper d''un coup 
    262                   ncdf_varget,cdfid,name,res,offset=[jpi-key+ixmin+premierx $ 
    263                                                      ,premiery+iyminmesh-iymindta $ 
    264                                                      ,premiertps] $ 
     276; ......... part, we can directly extract the array in one piece 
     277               lastx LE key-1:BEGIN ; on peut tout couper d''un coup 
     278                  ncdf_varget,cdfid,name,res,offset=[jpi-key+ixmin+firstx $ 
     279                                                     ,firsty+iymin $ 
     280                                                     ,firsttps] $ 
    265281                   ,count=[nx,ny,jpt], _extra = ex 
    266282               END 
    267283               ELSE:BEGIN       ; Le tableau est separe en 2 morceaux... 
    268                   ncdf_varget,cdfid,name,tab1,offset=[jpi-key+ixmin+premierx $ 
    269                                                       ,premiery+iyminmesh-iymindta $ 
    270                                                       ,premiertps] $ 
    271                    ,count=[key-1-premierx+1,ny,jpt], _extra = ex 
    272                   ncdf_varget,cdfid,name,tab2,offset=[ixmin $ 
    273                                                       ,premiery+iyminmesh-iymindta $ 
    274                                                       ,premiertps] $ 
    275                    ,count=[dernierx-key+1,ny,jpt], _extra = ex 
     284; ......... part, first part... 
     285                  ncdf_varget,cdfid,name,tab1,offset=[jpi-key+ixmin+firstx $ 
     286                                                      ,firsty+iymin $ 
     287                                                      ,firsttps] $ 
     288                   ,count=[key-1-firstx+1,ny,jpt], _extra = ex 
     289; --------- part, second part... 
     290                  ncdf_varget,cdfid,name,tab2,offset=[ixmin $ 
     291                                                      ,firsty+iymin $ 
     292                                                      ,firsttps] $ 
     293                   ,count=[lastx-key+1,ny,jpt], _extra = ex 
    276294                  res = [temporary(tab1), temporary(tab2)] 
    277295               END 
    … …  
    282300; on shift et il y a un stride 
    283301;,,,,,,,,,,,,,,,,,,,,,,,, 
     302; --------- part, we can directly extract the array in one piece 
    284303            case 1 OF           ; case sur la facon de fire le champ 
    285                premierx GE ceil(1.*key/key_stride[0]):BEGIN ; on peut tout coupe d''un coup 
    286                   ncdf_varget,cdfid,name,res,offset=[ixmin-key+premierx*key_stride[0] $ 
    287                                                      ,premiery*key_stride[1]+iyminmesh-iymindta $ 
    288                                                      ,premiertps] $ 
     304               firstx*key_stride[0] GE key:BEGIN  
     305                  ncdf_varget,cdfid,name,res,offset=[ixmin $ 
     306                                                     +(key_stride[0]-1)-((key-1) MOD key_stride[0]) $ 
     307                                                     +( firstx-((key-1)/key_stride[0]+1) )*key_stride[0] $ 
     308                                                     ,firsty*key_stride[1]+iymin $ 
     309                                                     ,firsttps] $ 
    289310                   ,count=[nx,ny,jpt], stride = [key_stride[0:1], 1], _extra = ex 
    290311               END 
    291                dernierx LE (key-1)/key_stride[0]:BEGIN ; on peut tout coupe d''un coup 
    292                   ncdf_varget,cdfid,name,res,offset=[jpitotal-key+ixmin+premierx*key_stride[0] $ 
    293                                                      ,premiery*key_stride[1]+iyminmesh-iymindta $ 
    294                                                      ,premiertps] $ 
     312; ......... part, we can directly extract the array in one piece 
     313               lastx*key_stride[0] LE key-1:BEGIN  
     314                  ncdf_varget,cdfid,name,res,offset=[jpitotal-key+ixmin+firstx*key_stride[0] $ 
     315                                                     ,firsty*key_stride[1]+iymin $ 
     316                                                     ,firsttps] $ 
    295317                   ,count=[nx,ny,jpt], stride = [key_stride[0:1], 1], _extra = ex 
    296318               END 
    297                ELSE:BEGIN       ; le tableau est separe en 2 morceaux... 
    298                   ncdf_varget,cdfid,name,tab1,offset=[jpitotal-key+ixmin+premierx*key_stride[0] $ 
    299                                                       ,premiery*key_stride[1]+iyminmesh-iymindta $ 
    300                                                       ,premiertps] $ 
    301                    ,count=[(key-1-premierx*key_stride[0]+1)/key_stride[0],ny,jpt] $ 
    302                    , stride = [key_stride[0:1], 1], _extra = ex 
    303 ; il en manque un bout? 
    304                   IF ( (key-1-premierx*key_stride[0]) MOD key_stride[0] ) EQ 0 then begin 
    305                      ncdf_varget,cdfid,name,tab1bis,offset=[jpitotal-1+ixmin $ 
    306                                                             ,premiery*key_stride[1]+iyminmesh-iymindta $ 
    307                                                             ,premiertps] $ 
    308                       ,count=[1,ny,jpt] $ 
    309                       , stride = [1, key_stride[1:1], 1], _extra = ex 
    310                      tab1bis = reform(tab1bis, 1, ny, jpt, /over) 
    311                      tab1 = [temporary(tab1), temporary(tab1bis)] 
    312                   ENDIF 
    313 ; 2eme bout... 
    314                   ncdf_varget,cdfid,name,tab2,offset=[ceil(1.*key/key_stride[0])*key_stride[0]-key+ixmin $ 
    315                                                       ,premiery*key_stride[1]+iyminmesh-iymindta $ 
    316                                                       ,premiertps] $ 
    317                    ,count=[(dernierx*key_stride[0]-key+1)/key_stride[0],ny,jpt] $ 
    318                    , stride = [key_stride[0:1], 1], _extra = ex 
     319               ELSE:BEGIN       
     320; ......... part, first part... 
     321                  ncdf_varget,cdfid,name,tab1,offset=[jpitotal-key+ixmin+firstx*key_stride[0] $ 
     322                                                      ,firsty*key_stride[1]+iymin $ 
     323                                                      ,firsttps] $ 
     324                   ,count=[(key-firstx*key_stride[0]-1)/key_stride[0]+1,ny,jpt] $ 
     325                    , stride = [key_stride[0:1], 1], _extra = ex 
     326; --------- part, second part... 
     327                  ncdf_varget,cdfid,name,tab2,offset=[ixmin $ 
     328                                                      +(key_stride[0]-1)-((key-1) MOD key_stride[0]) $ 
     329                                                      , firsty*key_stride[1]+iymin $ 
     330                                                      ,firsttps] $ 
     331                   ,count=[nx-(key-firstx*key_stride[0]-1)/key_stride[0]-1,ny,jpt] $ 
     332                    , stride = [key_stride[0:1], 1], _extra = ex 
    319333; on recolle le tout 
    320334                  res = [temporary(tab1), temporary(tab2)] 
    … …  
    334348; on ne shift pas et il n''y a pas de stride 
    335349;,,,,,,,,,,,,,,,,,,,,,,,, 
    336             ncdf_varget,cdfid,name,res,offset=[premierx+ixminmesh-ixmindta $ 
    337                                                ,premiery+iyminmesh-iymindta $ 
    338                                                       ,premierz+izminmesh-izmindta $ 
    339                                                ,premiertps] $ 
     350            ncdf_varget,cdfid,name,res,offset=[firstx+ixmin $ 
     351                                               ,firsty+iymin $ 
     352                                                      ,firstz+izmin $ 
     353                                               ,firsttps] $ 
    340354             ,count=[nx,ny,nz,jpt], _extra = ex 
    341355         END 
    … …  
    344358; on ne shift pas mais il y a un stride 
    345359;,,,,,,,,,,,,,,,,,,,,,,,, 
    346             ncdf_varget,cdfid,name,res,offset=[premierx*key_stride[0]+ixminmesh-ixmindta $ 
    347                                                ,premiery*key_stride[1]+iyminmesh-iymindta $ 
    348                                                ,premierz*key_stride[2]+izminmesh-izmindta $ 
    349                                                ,premiertps] $ 
    350              ,count=[nx,ny,nz,jpt], stride = [key_stride[0:2], 1], _extra = ex 
     360            ncdf_varget,cdfid,name,res,offset=[firstx*key_stride[0]+ixmin $ 
     361                                               ,firsty*key_stride[1]+iymin $ 
     362                                               ,firstz*key_stride[2]+izmin $ 
     363                                               ,firsttps] $ 
     364             ,count=[nx,ny,nz,jpt], stride = [key_stride, 1], _extra = ex 
    351365         END 
    352366;,,,,,,,,,,,,,,,,,,,,,,,, 
    … …  
    355369;,,,,,,,,,,,,,,,,,,,,,,,, 
    356370            case 1 of 
    357                premierx GE key:BEGIN ; on peut tout couper d''un coup 
    358                   ncdf_varget,cdfid,name,res,offset=[ixmin-key+premierx $ 
    359                                                      ,premiery+iyminmesh-iymindta $ 
    360                                                      ,premierz+izminmesh-izmindta $ 
    361                                                      ,premiertps] $ 
     371               firstx GE key:BEGIN ; on peut tout couper d''un coup 
     372                  ncdf_varget,cdfid,name,res,offset=[ixmin-key+firstx $ 
     373                                                     ,firsty+iymin $ 
     374                                                     ,firstz+izmin $ 
     375                                                     ,firsttps] $ 
    362376                   ,count=[nx,ny,nz,jpt], _extra = ex 
    363377               END 
    364                dernierx LE key-1:BEGIN ; on peut tout couper d''un coup 
    365                   ncdf_varget,cdfid,name,res,offset=[jpi-key+ixmin+premierx $ 
    366                                                      ,premiery+iyminmesh-iymindta $ 
    367                                                      ,premierz+izminmesh-izmindta $ 
    368                                                      ,premiertps] $ 
     378               lastx LE key-1:BEGIN ; on peut tout couper d''un coup 
     379                  ncdf_varget,cdfid,name,res,offset=[jpi-key+ixmin+firstx $ 
     380                                                     ,firsty+iymin $ 
     381                                                     ,firstz+izmin $ 
     382                                                     ,firsttps] $ 
    369383                   ,count=[nx,ny,nz,jpt], _extra = ex 
    370384               END 
    371385               ELSE:BEGIN       ; Le tableau est separe en 2 morceaux... 
    372                   ncdf_varget,cdfid,name,tab1,offset=[jpi-key+ixmin+premierx $ 
    373                                                       ,premiery+iyminmesh-iymindta $ 
    374                                                       ,premierz+izminmesh-izmindta $ 
    375                                                       ,premiertps] $ 
    376                    ,count=[key-1-premierx+1,ny,nz,jpt], _extra = ex 
    377                   ncdf_varget,cdfid,name,tab2,offset=[ixmin $ 
    378                                                       ,premiery+iyminmesh-iymindta $ 
    379                                                       ,premierz+izminmesh-izmindta $ 
    380                                                       ,premiertps] $ 
    381                    ,count=[dernierx-key+1,ny,nz,jpt], _extra = ex 
     386                  ncdf_varget,cdfid,name,tab1,offset=[jpi-key+ixmin+firstx $ 
     387                                                      ,firsty+iymin $ 
     388                                                      ,firstz+izmin $ 
     389                                                      ,firsttps] $ 
     390                   ,count=[key-1-firstx+1,ny,nz,jpt], _extra = ex 
     391                  ncdf_varget,cdfid,name,tab2,offset=[ixmin $ 
     392                                                      ,firsty+iymin $ 
     393                                                      ,firstz+izmin $ 
     394                                                      ,firsttps] $ 
     395                   ,count=[lastx-key+1,ny,nz,jpt], _extra = ex 
    382396                  res = [temporary(tab1), temporary(tab2)] 
    383397               END 
    … …  
    389403;,,,,,,,,,,,,,,,,,,,,,,,, 
    390404            case 1 OF           ; case sur la facon de fire le champ 
    391                premierx GE ceil(1.*key/key_stride[0]):BEGIN ; on peut tout coupe d''un coup 
    392                   ncdf_varget,cdfid,name,res,offset=[ixmin-key+premierx*key_stride[0] $ 
    393                                                      ,premiery*key_stride[1]+iyminmesh-iymindta $ 
    394                                                      ,premierz*key_stride[2]+izminmesh-izmindta $ 
    395                                                      ,premiertps] $ 
    396                    ,count=[nx,ny,nz,jpt], stride = [key_stride[0:2], 1], _extra = ex 
    397                END 
    398                dernierx LE (key-1)/key_stride[0]:BEGIN ; on peut tout coupe d''un coup 
    399                   ncdf_varget,cdfid,name,res,offset=[jpitotal-key+ixmin+premierx*key_stride[0] $ 
    400                                                      ,premiery*key_stride[1]+iyminmesh-iymindta $ 
    401                                                      ,premierz*key_stride[2]+izminmesh-izmindta $ 
    402                                                      ,premiertps] $ 
    403                    ,count=[nx,ny,nz,jpt], stride = [key_stride[0:2], 1], _extra = ex 
     405               firstx*key_stride[0] GE key:BEGIN ; on peut tout coupe d''un coup 
     406                  ncdf_varget,cdfid,name,res,offset=[ixmin $ 
     407                                                     +(key_stride[0]-1)-((key-1) MOD key_stride[0]) $ 
     408                                                     +(firstx-((key-1)/key_stride[0]+1))*key_stride[0] $ 
     409                                                     ,firsty*key_stride[1]+iymin $ 
     410                                                     ,firstz*key_stride[2]+izmin $ 
     411                                                     ,firsttps] $ 
     412                   ,count=[nx,ny,nz,jpt], stride = [key_stride, 1], _extra = ex 
     413               END 
     414               lastx*key_stride[0] LE key-1:BEGIN ; on peut tout coupe d''un coup 
     415                  ncdf_varget,cdfid,name,res,offset=[jpitotal-key+ixmin+firstx*key_stride[0] $ 
     416                                                     ,firsty*key_stride[1]+iymin $ 
     417                                                     ,firstz*key_stride[2]+izmin $ 
     418                                                     ,firsttps] $ 
     419                   ,count=[nx,ny,nz,jpt], stride = [key_stride, 1], _extra = ex 
    404420               END 
    405421               ELSE:BEGIN       ; le tableau est separe en 2 morceaux... 
    406                   ncdf_varget,cdfid,name,tab1,offset=[jpitotal-key+ixmin+premierx*key_stride[0] $ 
    407                                                       ,premiery*key_stride[1]+iyminmesh-iymindta $ 
    408                                                      ,premierz*key_stride[2]+izminmesh-izmindta $ 
    409                                                       ,premiertps] $ 
    410                    ,count=[(key-1-premierx*key_stride[0]+1)/key_stride[0],ny,nz,jpt] $ 
    411                    , stride = [key_stride[0:2], 1], _extra = ex 
    412 ; il en manque un bout? 
    413                   IF ( (key-1-premierx*key_stride[0]) MOD key_stride[0] ) EQ 0 then begin 
    414                      ncdf_varget,cdfid,name,tab1bis,offset=[jpitotal-1+ixmin $ 
    415                                                             ,premiery*key_stride[1]+iyminmesh-iymindta $ 
    416                                                      ,premierz*key_stride[2]+izminmesh-izmindta $ 
    417                                                             ,premiertps] $ 
    418                       ,count=[1,ny,nz,jpt] $ 
    419                       , stride = [1, key_stride[1:2], 1], _extra = ex 
    420                      tab1bis = reform(tab1bis, 1, ny, nz, jpt, /over) 
    421                      tab1 = [temporary(tab1), temporary(tab1bis)] 
    422                   ENDIF 
     422                  ncdf_varget,cdfid,name,tab1,offset=[jpitotal-key+ixmin+firstx*key_stride[0] $ 
     423                                                      ,firsty*key_stride[1]+iymin $ 
     424                                                     ,firstz*key_stride[2]+izmin $ 
     425                                                      ,firsttps] $ 
     426                   ,count=[(key-firstx*key_stride[0]-1)/key_stride[0]+1,ny,nz,jpt] $ 
     427                   , stride = [key_stride, 1], _extra = ex 
    423428; 2eme bout... 
    424                   ncdf_varget,cdfid,name,tab2,offset=[ceil(1.*key/key_stride[0])*key_stride[0]-key+ixmin $ 
    425                                                       ,premiery*key_stride[1]+iyminmesh-iymindta $ 
    426                                                      ,premierz*key_stride[2]+izminmesh-izmindta $ 
    427                                                       ,premiertps] $ 
    428                    ,count=[(dernierx*key_stride[0]-key+1)/key_stride[0],ny,nz,jpt] $ 
    429                    , stride = [key_stride[0:2], 1], _extra = ex 
     429                  ncdf_varget,cdfid,name,tab2,offset=[ixmin $ 
     430                                                      +(key_stride[0]-1)-((key-1) MOD key_stride[0]) $ 
     431                                                      , firsty*key_stride[1]+iymin $ 
     432                                                      , firstz*key_stride[2]+izmin $ 
     433                                                      , firsttps] $ 
     434                   ,count=[nx-(key-firstx*key_stride[0]-1)/key_stride[0]-1,ny,nz,jpt] $ 
     435                   , stride = [key_stride, 1], _extra = ex 
    430436; on recolle le tout 
    431437                  res = [temporary(tab1), temporary(tab2)] 
Note: See TracChangeset for help on using the changeset viewer.