Changeset 6 for trunk


Ignore:
Timestamp:
11/23/06 10:01:24 (17 years ago)
Author:
pinsard
Message:

improvement of interpolation of condmag

File:
1 moved

Legend:

Unmodified
Added
Removed
  • trunk/condmag_on_orca.pro

    r5 r6  
    33; @file_comments 
    44; interpolate condmag.nc file (sediment and magnetic fields) on ORCA grid 
    5 ; and produce one file cond_sed_<orcares>.nc 
     5; and produce cond_sed_<orcares>.nc and Br_<orcares>.nc 
    66; 
    77; @categories 
     
    2727;  - must have meshmask in ${GEOMAG_ID}/ 
    2828;  - must not have cond_sed_*.nc in ${GEOMAG_OD}/ 
     29;  - must not have Br_*.nc in ${GEOMAG_OD}/ 
    2930; 
    3031; @todo 
     
    4243; @post 
    4344; cond_sed_<I>orcasres</I>.nc is now present in ${GEOMAG_OD}/ 
     45; Br_<I>orcasres</I>.nc is now present in ${GEOMAG_OD}/ 
    4446; see geomag_env.sh 
    4547; 
    4648; @examples 
    47 ; IDL> .run condmag_cond_sed_on_orca 
    48 ; IDL> condmag_cond_sed_on_orca, 'ORCA2', 'bilinear','T' 
     49; IDL> .run condmag_on_orca 
     50; IDL> condmag_on_orca, 'ORCA2', 'bilinear','T' 
    4951; 
    5052; @history 
     53; fplod 2006-11-23T08:57:10Z aedon.locean-ipsl.upmc.fr (Darwin) 
     54; rename to condmag_on_orca.pro 
     55; fplod 2006-11-22T15:10:31Z aedon.locean-ipsl.upmc.fr (Darwin) 
     56; add creation of Br_<I>orcasres</I>.nc 
     57; ++ not very beautifuly implemented 
     58; fplod 2006-11-22T10:38:51Z aedon.locean-ipsl.upmc.fr (Darwin) 
     59; always use msg = 'iii : ...' or msg = 'eee : ...'  
     60; information on open io 
    5161; fplod 2006-11-20T10:32:02Z aedon.locean-ipsl.upmc.fr (Darwin) 
    5262; cleaning before introduction into svn repository 
     
    6676; 
    6777;- 
    68 PRO condmag_cond_sed_on_orca, orcares, method, gridtype 
     78PRO condmag_on_orca, orcares, method, gridtype 
    6979; 
    7080  compile_opt idl2, strictarrsubs 
     
    7888  CASE orcares OF 
    7989     'ORCA2': BEGIN 
     90                 msg = 'iii : valid orcares parameter = ' +  orcares 
     91                 PRINT, msg 
    8092                 filename_oce = 'meshmask_bab.nc' 
    81                  PRINT, 'valid orcares parameter' 
    8293              END 
    8394      ELSE  : BEGIN 
    84                  msg = 'orcares must be ORCA2' 
     95                 msg = 'eee : invalid orcares parameter = ' + orcares 
     96                 PRINT, msg 
     97                 msg = 'eee : orcares must be ORCA2' 
    8598                 PRINT, msg 
    8699                 RETURN 
     
    91104  CASE method OF 
    92105     'bilinear': BEGIN 
    93                     PRINT, 'valid method parameter' 
     106                    msg = 'iii : valid method parameter = ' + method 
     107                    PRINT, msg 
    94108                 END 
    95109     'imoms3'  : BEGIN 
    96                     PRINT, 'valid method parameter' 
     110                    msg = 'iii : valid method parameter = ' + method 
     111                    PRINT, msg 
    97112                 END 
    98113     ELSE      : BEGIN 
    99                     msg = 'method must be bilinear or imoms3' 
     114                    msg = 'eee : invalid method parameter = ' + method  
     115                    PRINT, msg 
     116                    msg = 'eee : method must be bilinear or imoms3' 
    100117                    PRINT, msg 
    101118                    RETURN 
     
    107124  CASE gridtype OF 
    108125     'T' : BEGIN 
    109               PRINT, 'valid gridtype parameter' 
     126              msg = 'iii : valid gridtype parameter = ' + gridtype 
     127              PRINT, msg 
    110128           END 
    111129     'U' : BEGIN 
    112               PRINT, 'valid gridtype parameter' 
     130              msg = 'iii : valid gridtype parameter = ' + gridtype 
     131              PRINT, msg 
    113132           END 
    114133     'V' : BEGIN 
    115               PRINT, 'valid gridtype parameter' 
     134              msg = 'iii : valid gridtype parameter = ' + gridtype 
     135              PRINT, msg 
    116136           END 
    117137     ELSE: BEGIN 
    118               msg = 'gridtype must be T, U or V' 
     138              msg = 'eee : invalid gridtype parameter = ' + gridtype  
     139              PRINT, msg 
     140              msg = 'eee : gridtype must be T, U or V' 
    119141              PRINT, msg 
    120142              RETURN 
     
    128150  CASE geomag_id_env OF 
    129151     ''  :  BEGIN 
    130               msg = '${GEOMAG_ID} is not defined' 
     152              msg = 'eee : ${GEOMAG_ID} is not defined' 
    131153              PRINT, msg 
    132154              RETURN 
    133155            END 
    134156     ELSE: BEGIN 
    135              msg = '${GEOMAG_ID} is ' + geomag_id_env 
     157             msg = 'iii : ${GEOMAG_ID} is ' + geomag_id_env 
    136158             PRINT, msg 
    137159           END  
     
    142164; existence and protection of ${GEOMAG_ID} 
    143165  IF (FILE_TEST(iodirin, /DIRECTORY,/EXECUTABLE , /READ) EQ 0) THEN BEGIN 
    144      msg = 'the directory' + iodirin  + ' is not accessible.' 
     166     msg = 'eee : the directory' + iodirin  + ' is not accessible.' 
    145167     PRINT, msg 
    146168     RETURN 
     
    151173; 
    152174; check if this file exists 
    153   fullfilename_condmag = isafile(iodirin + filename_condmag) 
     175  fullfilename_condmag = isafile(iodirin + filename_condmag, NEW=0,/MUST_EXIST) 
    154176  IF fullfilename_condmag[0] EQ '' THEN BEGIN 
    155      msg = 'the file ' + fullfilename_condmag + ' was not found.' 
     177     msg = 'eee : the file ' + fullfilename_condmag + ' was not found.' 
    156178     PRINT, msg 
    157179     RETURN 
     
    160182; protection 
    161183  IF (FILE_TEST(fullfilename_condmag[0], /READ) EQ 0) THEN BEGIN 
    162      msg = 'the file ' + fullfilename_condmag[0] + ' is not readable.' 
     184     msg = 'eee : the file ' + fullfilename_condmag[0] + ' is not readable.' 
    163185     PRINT, msg 
    164186     RETURN 
     
    167189; mesh mask 
    168190; check if this file exists 
    169   fullfilename_oce = isafile(iodirin + filename_oce) 
     191  fullfilename_oce = isafile(iodirin + filename_oce, NEW=0,/MUST_EXIST, RECURSIVE=0) 
    170192  IF fullfilename_oce[0] EQ '' THEN BEGIN 
    171      msg = 'the file ' + fullfilename_oce + ' was not found.' 
     193     msg = 'eee : the file ' + fullfilename_oce + ' was not found.' 
    172194     PRINT, msg 
    173195     RETURN 
     
    176198; protection 
    177199  IF (FILE_TEST(fullfilename_oce[0], /READ) EQ 0) THEN BEGIN 
    178      msg = 'the file ' + fullfileoce_condmag[0] + ' is not readable.' 
     200     msg = 'eee : the file ' + fullfileoce_condmag[0] + ' is not readable.' 
    179201     PRINT, msg 
    180202     RETURN 
     
    185207  CASE geomag_od_env OF 
    186208     '' : BEGIN 
    187              msg = '${GEOMAG_OD} is not defined' 
     209             msg = 'eee : ${GEOMAG_OD} is not defined' 
    188210             PRINT, msg 
    189211             RETURN 
    190212          END 
    191213     ELSE: BEGIN 
    192              msg = '${GEOMAG_OD} is ' + geomag_od_env 
     214             msg = 'iii : ${GEOMAG_OD} is ' + geomag_od_env 
    193215             PRINT, msg 
    194216           END  
     
    200222; existence and protection 
    201223  IF (FILE_TEST(iodirout, /DIRECTORY,/WRITE) EQ 0) THEN BEGIN 
    202      msg = 'the directory' + iodirout  + ' was not found.' 
    203      PRINT, msg 
    204      RETURN 
    205   ENDIF 
    206 ; 
    207 ; build output filename 
     224     msg = 'eee : the directory' + iodirout  + ' was not found.' 
     225     PRINT, msg 
     226     RETURN 
     227  ENDIF 
     228; 
     229; build output filenames 
    208230  filename_cond_sed = 'cond_sed' + '_' + orcares +'.nc' 
    209231  fullfilename_cond_sed = iodirout + filename_cond_sed 
     
    211233; in order to avoid unexpected overwritten 
    212234  IF (FILE_TEST(fullfilename_cond_sed) EQ 1) THEN BEGIN 
    213      msg = 'the file ' + fullfilename_cond_sed  + ' already exists.' 
     235     msg = 'eee : the file ' + fullfilename_cond_sed  + ' already exists.' 
     236     PRINT, msg 
     237     RETURN 
     238  ENDIF 
     239; 
     240  filename_br = 'Br' + '_' + orcares +'.nc' 
     241  fullfilename_br = iodirout + filename_br 
     242; 
     243; in order to avoid unexpected overwritten 
     244  IF (FILE_TEST(fullfilename_br) EQ 1) THEN BEGIN 
     245     msg = 'eee : the file ' + fullfilename_br  + ' already exists.' 
    214246     PRINT, msg 
    215247     RETURN 
     
    220252  condmaglatname = 'la' 
    221253  varname_cond_sed = 'cond_sed' 
     254  varname_br = 'Br' 
    222255; 
    223256;---- 
     
    251284; 
    252285; reading condmag file 
    253   netcdf_id_condmag = NCDF_OPEN(fullfilename_condmag[0]) 
     286  netcdf_id_condmag = NCDF_OPEN(fullfilename_condmag[0], /NOWRITE) 
     287  msg = 'iii : ' + fullfilename_condmag[0] + ' opened for read' 
     288  PRINT, msg 
     289   
    254290  varinq_cond_sed = NCDF_VARINQ(netcdf_id_condmag, varname_cond_sed) 
    255 ; 
    256 ;--------------------------- 
    257 ; Creation of the NetCdf file for cond_sed 
     291  varinq_br = NCDF_VARINQ(netcdf_id_condmag, varname_br) 
     292; 
     293;--------------------------- 
     294; Creation of the NetCdf file for cond_sed and Br 
    258295;--------------------------- 
    259296; 
    260297  netcdf_id_cond_sed = NCDF_CREATE(fullfilename_cond_sed, /clobber) 
    261298  NCDF_CONTROL, netcdf_id_cond_sed, /NOFILL 
     299  netcdf_id_br = NCDF_CREATE(fullfilename_br, /clobber) 
     300  NCDF_CONTROL, netcdf_id_br, /NOFILL 
     301; 
    262302; dimension 
    263303  dimidx = NCDF_DIMDEF(netcdf_id_cond_sed, 'x' ,  jpio) 
    264304  dimidy = NCDF_DIMDEF(netcdf_id_cond_sed, 'y' ,  jpjo) 
    265305  dimidt = NCDF_DIMDEF(netcdf_id_cond_sed, 'lo', /UNLIMITED) 
     306  dimidx = NCDF_DIMDEF(netcdf_id_br, 'x' ,  jpio) 
     307  dimidy = NCDF_DIMDEF(netcdf_id_br, 'y' ,  jpjo) 
     308  dimidt = NCDF_DIMDEF(netcdf_id_br, 'lo', /UNLIMITED) 
    266309; 
    267310; global attributes 
    268   NCDF_ATTPUT, netcdf_id_cond_sed, 'Conventions', 'GDT 1.2', /GLOBAL 
     311  NCDF_ATTPUT, netcdf_id_cond_sed, 'Conventions', 'GDT 1.2', /GLOBAL ; ++ conformité 
    269312  NCDF_ATTPUT, netcdf_id_cond_sed, 'file_name'  , filename_cond_sed, /GLOBAL 
    270   NCDF_ATTPUT, netcdf_id_cond_sed, 'Title'      , 'Monthly Levitus Sea Salinity corrected', /GLOBAL 
     313  NCDF_ATTPUT, netcdf_id_cond_sed, 'Title'      , 'Conductance', /GLOBAL 
     314  NCDF_ATTPUT, netcdf_id_br, 'Conventions', 'GDT 1.2', /GLOBAL ; ++ conformité 
     315  NCDF_ATTPUT, netcdf_id_br, 'file_name'  , filename_Br, /GLOBAL 
     316  NCDF_ATTPUT, netcdf_id_br, 'Title'      , 'Magnetic field', /GLOBAL 
    271317; 
    272318; declaration of variables 
     
    279325  NCDF_ATTPUT, netcdf_id_cond_sed, varid[0], 'valid_max', omax,/FLOAT 
    280326  NCDF_ATTPUT, netcdf_id_cond_sed, varid[0], 'long_name', 'Longitude at t-point' 
     327  varid[0] = NCDF_VARDEF(netcdf_id_br, 'nav_lon'  , [dimidx, dimidy], /FLOAT) 
     328  NCDF_ATTPUT, netcdf_id_br, varid[0], 'units'    , 'degrees_east' 
     329  NCDF_ATTPUT, netcdf_id_br, varid[0], 'valid_min', min(olon, max = omax),/FLOAT 
     330  NCDF_ATTPUT, netcdf_id_br, varid[0], 'valid_max', omax,/FLOAT 
     331  NCDF_ATTPUT, netcdf_id_br, varid[0], 'long_name', 'Longitude at t-point' 
    281332; 
    282333  varid[1] = NCDF_VARDEF(netcdf_id_cond_sed, 'nav_lat'  , [dimidx, dimidy], /FLOAT) 
     
    285336  NCDF_ATTPUT, netcdf_id_cond_sed, varid[1], 'valid_max', omax,/FLOAT 
    286337  NCDF_ATTPUT, netcdf_id_cond_sed, varid[1], 'long_name', 'Latitude at t-point' 
     338  varid[1] = NCDF_VARDEF(netcdf_id_br, 'nav_lat'  , [dimidx, dimidy], /FLOAT) 
     339  NCDF_ATTPUT, netcdf_id_br, varid[1], 'units'    , 'degrees_north' 
     340  NCDF_ATTPUT, netcdf_id_br, varid[1], 'valid_min', min(olat, max = omax),/FLOAT 
     341  NCDF_ATTPUT, netcdf_id_br, varid[1], 'valid_max', omax,/FLOAT 
     342  NCDF_ATTPUT, netcdf_id_br, varid[1], 'long_name', 'Latitude at t-point' 
    287343; 
    288344  varid[2] = NCDF_VARDEF(netcdf_id_cond_sed, 'time'     , [dimidt], /FLOAT) 
    289 ; 
    290 ; each of the the two files to produce contains a specific variable 
     345  varid[2] = NCDF_VARDEF(netcdf_id_br, 'time'     , [dimidt], /FLOAT) 
     346; 
     347; each of the two files to produce contains a specific variable 
     348; 
    291349  varid_cond_sed = lonarr(1) 
    292350  varid_cond_sed[0] = NCDF_VARDEF(netcdf_id_cond_sed, varname_cond_sed, [dimidx, dimidy, dimidt], /FLOAT) 
    293351; 
    294352; variable 3 
    295   NCDF_ATTPUT, netcdf_id_cond_sed, varid_cond_sed[0], 'long_name', 'Salinity' ; ++ faux et non cf 
    296   NCDF_ATTPUT, netcdf_id_cond_sed, varid_cond_sed[0], 'units', 'siemens' 
     353  NCDF_ATTPUT, netcdf_id_cond_sed, varid_cond_sed[0], 'long_name', 'Conductance' ; ++ non cf  
     354  NCDF_ATTPUT, netcdf_id_cond_sed, varid_cond_sed[0], 'units', 'siemens' ; ++ récupérer l'unite de cond_sed dans condmag.nc 
    297355; pour min  et max, il faut avoir lu la variable ... cf. plus bas ++ 
    298356; donc pour l'instant on les met à valeur manquante 
    299357  NCDF_ATTPUT, netcdf_id_cond_sed, varid_cond_sed[0], 'valid_min', !VALUES.F_NAN ,/FLOAT 
    300358  NCDF_ATTPUT, netcdf_id_cond_sed, varid_cond_sed[0], 'valid_max', !VALUES.F_NAN ,/FLOAT 
    301 ;--------------------------- 
    302 ; end of header definition, writing of the NetCdf file 
     359; 
     360  varid_br = lonarr(1) 
     361  varid_br[0] = NCDF_VARDEF(netcdf_id_br, varname_br, [dimidx, dimidy, dimidt], /FLOAT) 
     362; 
     363; variable 3 
     364  NCDF_ATTPUT, netcdf_id_br, varid_br[0], 'long_name', 'magnetic field' ; ++ non cf  
     365  NCDF_ATTPUT, netcdf_id_br, varid_br[0], 'units', 'tesla' ; ++ récupérer l'unité de Br dans condmag.nc 
     366; pour min  et max, il faut avoir lu la variable ... cf. plus bas ++ 
     367; donc pour l'instant on les met à valeur manquante 
     368  NCDF_ATTPUT, netcdf_id_br, varid_br[0], 'valid_min', !VALUES.F_NAN ,/FLOAT 
     369  NCDF_ATTPUT, netcdf_id_cond_sed, varid_cond_sed[0], 'valid_max', !VALUES.F_NAN ,/FLOAT 
     370; 
     371;--------------------------- 
     372; end of header definition, writing of the NetCdf files 
    303373;--------------------------- 
    304374  NCDF_CONTROL, netcdf_id_cond_sed, /ENDEF 
    305   NCDF_CLOSE, netcdf_id_condmag ; ++ si le close n'est pa ici ça merde mais je ne sais pas pourquoi 
     375  NCDF_CONTROL, netcdf_id_br, /ENDEF 
     376  NCDF_CLOSE, netcdf_id_condmag ; ++ si le close n'est pas ici ça merde mais je ne sais pas pourquoi 
     377; 
    306378; grid 
    307 ;--------------------------- 
    308   NCDF_VARPUT, netcdf_id_cond_sed, 'nav_lon', TEMPORARY(olon) 
    309   NCDF_VARPUT, netcdf_id_cond_sed, 'nav_lat', TEMPORARY(olat) 
     379  NCDF_VARPUT, netcdf_id_cond_sed, 'nav_lon', olon 
     380  NCDF_VARPUT, netcdf_id_cond_sed, 'nav_lat', olat 
    310381  NCDF_VARPUT, netcdf_id_cond_sed, varid[2], FLOAT(0.5) ; ++ c'est quoi cette valeur 
    311 ;--------------------------- 
    312 ; data 
     382  NCDF_VARPUT, netcdf_id_br, 'nav_lon', olon 
     383  NCDF_VARPUT, netcdf_id_br, 'nav_lat', olat 
     384  NCDF_VARPUT, netcdf_id_br, varid[2], FLOAT(0.5) ; ++ c'est quoi cette valeur 
     385 
    313386;--------------------------- 
    314387; réouverture du fichier ... mais pourquoi on le relit !!++ 
    315       netcdf_id_condmag = NCDF_OPEN(fullfilename_condmag[0]) 
    316 ; 
    317       NCDF_VARGET, netcdf_id_condmag, varname_cond_sed, varin_cond_sed ; , COUNT = [jpia, jpja, 1], OFFSET = [0, 0, 0] 
     388  netcdf_id_condmag = NCDF_OPEN(fullfilename_condmag[0],/NOWRITE) 
     389  msg = 'iii : ' + fullfilename_condmag[0] + ' reopened for read' 
     390  PRINT, msg 
     391; 
     392  NCDF_VARGET, netcdf_id_condmag, varname_cond_sed, varin_cond_sed ; , COUNT = [jpia, jpja, 1], OFFSET = [0, 0, 0] 
     393  NCDF_VARGET, netcdf_id_condmag, varname_br, varin_br ; , COUNT = [jpia, jpja, 1], OFFSET = [0, 0, 0] 
    318394; 
    319395; do the interpolation 
    320396      varin_cond_sed = TOTAL(weig*varin_cond_sed[addr], 1) 
    321397      varin_cond_sed = REFORM(varin_cond_sed, jpio, jpjo, /OVER) 
     398      varin_br = TOTAL(weig*varin_br[addr], 1) 
     399      varin_br = REFORM(varin_br, jpio, jpjo, /OVER) 
    322400; 
    323401      NCDF_CLOSE, netcdf_id_condmag 
    324402      varout_cond_sed = TEMPORARY(varin_cond_sed) 
     403      varout_br = TEMPORARY(varin_br) 
    325404; compute min max 
    326       minarr = min(varout_cond_sed, max = maxarr) 
     405      minarr_cond_sed = min(varout_cond_sed, max = maxarr_cond_sed) 
     406      minarr_br = min(varout_br, max = maxarr_br) 
    327407; 
    328408; put back the masked value 
    329409;++      IF bad[0] NE -1 THEN varout_cond_sed[TEMPORARY(bad)] = 32767 
     410;++      IF bad[0] NE -1 THEN varout_br[TEMPORARY(bad)] = 32767 
    330411; 
    331412; write the data 
    332413      NCDF_VARPUT, netcdf_id_cond_sed, varname_cond_sed, varout_cond_sed, COUNT = [jpio, jpjo, 1], OFFSET = [0, 0, 0] 
     414      NCDF_VARPUT, netcdf_id_br, varname_br, varout_br, COUNT = [jpio, jpjo, 1], OFFSET = [0, 0, 0] 
    333415; 
    334416; update min max attributes 
    335417  NCDF_CONTROL, netcdf_id_cond_sed, /REDEF 
    336   NCDF_ATTPUT, netcdf_id_cond_sed, varid_cond_sed[0], 'valid_min', minarr,/FLOAT 
    337   NCDF_ATTPUT, netcdf_id_cond_sed, varid_cond_sed[0], 'valid_max', maxarr,/FLOAT 
     418  NCDF_ATTPUT, netcdf_id_cond_sed, varid_cond_sed[0], 'valid_min', minarr_cond_sed,/FLOAT 
     419  NCDF_ATTPUT, netcdf_id_cond_sed, varid_cond_sed[0], 'valid_max', maxarr_cond_sed,/FLOAT 
    338420  NCDF_CONTROL, netcdf_id_cond_sed, /ENDEF 
    339 ; 
    340 ;--------------------------- 
    341 ; close the netcdf file 
     421  NCDF_CONTROL, netcdf_id_br, /REDEF 
     422  NCDF_ATTPUT, netcdf_id_br, varid_br[0], 'valid_min', minarr_br,/FLOAT 
     423  NCDF_ATTPUT, netcdf_id_br, varid_br[0], 'valid_max', maxarr_br,/FLOAT 
     424  NCDF_CONTROL, netcdf_id_br, /ENDEF 
     425; 
     426;--------------------------- 
     427; close the netcdf files 
    342428  NCDF_CLOSE, netcdf_id_cond_sed 
    343   msg =' done' 
     429  NCDF_CLOSE, netcdf_id_br 
     430; 
     431  msg = 'iii : ' + fullfilename_cond_sed + ' created' 
    344432  PRINT, msg 
     433  msg = 'iii : ' + fullfilename_br + ' created' 
     434  PRINT, msg 
    345435; 
    346436END 
Note: See TracChangeset for help on using the changeset viewer.