Changeset 191 for trunk/src/interp_erai_msl.pro

Timestamp:

07/04/12 17:52:46 (12 years ago)

Author:

pinsard

Message:

introduce new erai dataset : new terminology and new time unit and origin

File:

• trunk/src/interp_erai_msl.pro (modified) (15 diffs)

trunk/src/interp_erai_msl.pro

@@ -5 +5 @@
 ; ===================
 ;
-; .. function:: interp_erai_msl(yyyymmddb,yyyymmdde)
+; .. function:: interp_erai_msl(yyyymmddb,yyyymmdde, eraitype)
 ;
 ; DESCRIPTION
@@ -12 +12 @@
 ; Interpolation of msl from ERA-I grid to OAFLUX grid
 ;
+; If eraitype is set to 1,
 ; :file:`${PROJECT_ID}/20c3m_erai_msl_TROP_1989_2009.nc`
 ; containing
@@ -18 +19 @@
 ; :ref:`compute_erai_daily_region_2d.sh`.
 ;
+; If eraitype is set to 2,
+; :file:`${PROJECT_ID}/erai_msl_{yyyyb}_{yyyye}.nc`
+; containing
+; msl from ERA-I
+; has been produced by
+; :ref:`concat_eraiy.sh`.
+;
 ; :file:`${PROJECT_ID}/mask_oaflux_30N30S.nc`
 ; containing
@@ -44 +52 @@
 ;        digraph interp_erai_msl {
 ;
-;           file_in [shape=ellipse,fontname=Courier,label="${PROJECT_ID}/20c3m_erai_msl_TROP_1989_2009.nc"];
+;           file_in_erai1 [shape=ellipse,fontname=Courier,label="${PROJECT_ID}/20c3m_erai_msl_TROP_1989_2009.nc"];
+;           file_in_erai2 [shape=ellipse,fontname=Courier,label="${PROJECT_ID}/erai_msl_{yyyyb}_{yyyye}.nc"];
 ;           mask [shape=ellipse,fontname=Courier,label="${PROJECT_ID}/mask_oaflux_30N30S.nc"];
 ;           file_out [shape=ellipse,fontname=Courier,label="${PROJECT_OD}/erai_msl_{yyyymmdd}_{yyyymmdd}_oafluxgrid.nc"];
@@ -54 +63 @@
 ;           label="${PROJECT}/src/interp_erai_msl.pro"];
 ;
-;           {file_in mask} -> {interp_erai_msl} -> {file_out}
+;           {file_in_erai1 file_in_erai2 mask} -> {interp_erai_msl} -> {file_out}
 ;
 ;        }
@@ -66 +75 @@
 ;
 ; Used by :ref:`tropflux.sh`
+;
+; :ref:`compute_erai_daily_region_2d.sh`
+; :ref:`concat_eraiy.sh`
 ;
 ; :func:`report <saxo:report>`
@@ -81 +93 @@
 ; ========
 ;
-; ::
+; With ERA-I type 1::
 ;
 ;  .compile file_interp
 ;  yyyymmddb = 20000101L
 ;  yyyymmdde = 20001231L
-;  result = interp_erai_msl(yyyymmddb, yyyymmdde)
+;  eraitype = 1
+;  result = interp_erai_msl(yyyymmddb, yyyymmdde, eraitype)
+;  print, result
+;
+; With ERA-I type 2::
+;
+;  .compile file_interp
+;  yyyymmddb = 20000101L
+;  yyyymmdde = 20001231L
+;  eraitype = 2
+;  result = interp_erai_msl(yyyymmddb, yyyymmdde, eraitype)
 ;  print, result
 ;
@@ -92 +114 @@
 ; ====
 ;
-; use real output of :ref:`compute_erai_daily_region_2d.sh`
+; use real output of :ref:`compute_erai_daily_region_2d.sh` if eraitype is set to 1
 ;
 ; get_erai.py + grib to netcdf conversion
@@ -116 +138 @@
 ;
 ; $URL$
+;
+; - fplod 20120704T113504Z cratos (Linux)
+;
+;   * add eraitype parameter
+;   * conversion mPa to hPa
+;   * deal with time units (days vs hours) and time origin
 ;
 ; - fplod 20120320
@@ -150 +178 @@
 function interp_erai_msl $
          , yyyymmddb $
-         , yyyymmdde
+         , yyyymmdde $
+         , eraitype
 ;
 ;++compile_opt idl2, strictarrsubs, logical_predicate
@@ -165 +194 @@
 result = -1
 ;
-usage = 'result = interp_erai_msl(yyyymmddb, yyyymmdde)'
+usage = 'result = interp_erai_msl(yyyymmddb, yyyymmdde, eraitype)'
 nparam = N_PARAMS()
-IF (nparam NE 2) THEN BEGIN
+IF (nparam NE 3) THEN BEGIN
    ras = report(['Incorrect number of arguments.' $
    + '!C' $
@@ -212 +241 @@
 ;
 ; build data filename
-filename='20c3m_erai_msl_TROP_1989_2009.nc'
+IF eraitype EQ 1 THEN BEGIN
+    erai_var = 'msl'
+    filename='20c3m_erai_msl_TROP_1989_2009.nc'
+    tutoday = 24.d
+    tmporig = date2jul(19570101L,month=immo,day=iddo,year=iyyyyo)
+endif
+IF (eraitype eq 2) THEN BEGIN
+    erai_var = 'msl'
+    tmpb = date2jul(yyyymmddb,month=immb,day=iddb,year=iyyyyb)
+    tmpe = date2jul(yyyymmdde,month=imme,day=idde,year=iyyyye)
+    filename = 'erai_' + erai_var + '_' + string(iyyyyb,format='(I4.4)') + '_' + string(iyyyye,format='(I4.4)') + '.nc'
+    tutoday = 1.d
+    tmporig = date2jul(19500101L,month=immo,day=iddo,year=iyyyyo)
+ENDIF
 ;
 ; check if this file exists
@@ -222 +264 @@
    ras = report(msg)
    return, result
-ENDIF
+ENDIF ELSE BEGIN
+   msg = 'iii : the file ' + fullfilename + ' found.'
+   ras = report(msg)
+ENDELSE
 ;
 ; test if ${PROJECT_OD} defined
@@ -263 +308 @@
 print, 'lat grid ',min(latin),max(latin),latin[1]-latin[0]
 print, 'lon grid ',min(lonin),max(lonin),lonin[1]-lonin[0]
-;++mslin=read_ncdf("msl",yyyymmddb-.5d,yyymmdde,file=fullfilename,/nostr)
-mslin=read_ncdf("msl",19880101,20100930,file=fullfilename,/nostr)
-;
-timein=24.d*(time-julday(1,1,1957,0,0,0))
+mslin=read_ncdf(erai_var,yyyymmddb,yyyymmdde,file=fullfilename,/nostr)
+;
+help, time
+timein=tutoday*(time-julday(immo,iddo,iyyyyo,0,0,0))
 jptin=n_elements(timein)
 da=jul2date(time[0])
@@ -287 +332 @@
 ;
 ; conversion to hPa
-ncdf_getatt, fullfilename, 'msl', units=units
+ncdf_getatt, fullfilename, erai_var, units=units
 CASE units OF
     'Pascal': BEGIN
         mslin=mslin/100
+    END
+    'mPa': BEGIN
+        mslin=mslin*10
     END
     ELSE: BEGIN