Changeset 219 for trunk/SRC/ToBeReviewed/INIT/initncdf.pro

Timestamp:

03/14/07 16:44:26 (17 years ago)

Author:

smasson

Message:

get back to changeset:217

File:

• trunk/SRC/ToBeReviewed/INIT/initncdf.pro (modified) (6 diffs)

trunk/SRC/ToBeReviewed/INIT/initncdf.pro

@@ -28 +28 @@
 ; @keyword START1 {default=0}{type=scalar: 0 or 1}
 ; Index the axis from 1 instead of 0 when using
-; /xyindex and/or /zindex
+; /xyindex and/or zindex
 ;
 ; @keyword USEASMASK {type=scalar string}
@@ -39 +39 @@
 ;    3 NaN values if existing
 ;
+; @keyword XAXISNAME {default='x', 'longitude', 'nav_lon', 'lon', 'lon_rho' or 'NbLongitudes'}{type=scalar string}
+; A string giving the name of the variable in the file 
+; that contains the [xyz]axis. 
+;      
+; @keyword YAXISNAME {default='y', 'latitude', 'nav_lat','lat', 'lat_rho' or 'NbLatitudes'}{type=scalar string}
+; A string giving the name of the variable in the file 
+; that contains the [xyz]axis. 
+;
 ; @keyword ZAXISNAME {default='z', 'level', 'lev', 'depth...'}{type=scalar string}
 ; A string giving the name of the variable in the file 
@@ -55 +63 @@
 ; 
 ; @keyword _EXTRA
-; Used to pass keywords to computegrid and ncdf_getaxis
+; Used to pass keywords to computegrid
 ;
 ; @uses
@@ -80 +88 @@
 ;------------------------------------------------------------
 ;------------------------------------------------------------
-PRO initncdf, ncfilein $
+PRO initncdf, ncfilein, XAXISNAME = xaxisname, YAXISNAME = yaxisname $
               , ZAXISNAME = zaxisname, MASKNAME = maskname $
               , INVMASK = invmask, USEASMASK = useasmask $
@@ -105 +113 @@
 ; what is inside the file
   inside = ncdf_inquire(cdfid)
-;------------------------------------------------------------
-; name of all dimensions
-  namedim = strarr(inside.ndims)
-  for dimiq = 0, inside.ndims-1 do begin
-    ncdf_diminq, cdfid, dimiq, tmpname, value 
-    namedim[dimiq] = strlowcase(tmpname)
-  ENDFOR
 ;----------------------------------------------------------
 ; name of the variables
@@ -120 +121 @@
   ENDFOR
 ;----------------------------------------------------------
-; find the x/yaxis
- ncdf_getaxis, cdfid, dimidx, dimidy, xaxis, yaxis $
-               , START1 = start1, XYINDEX = xyindex, ROMSGRID = romsgrid, _extra = ex
+; find the xaxis
+  if keyword_set(xaxisname) then xaxisname = strlowcase(xaxisname) ELSE xaxisname = 'x'
+  xvarid = (where(namevar EQ xaxisname OR namevar EQ 'longitude' $
+                  OR namevar EQ 'nav_lon' OR namevar EQ 'lon' $
+                  OR namevar EQ 'lon_rho' OR namevar EQ 'nblongitudes'))[0]
+  if xvarid EQ -1 then begin
+    print, 'the xaxis was not found, check the use of XAXISNAME keyword'
+    stop
+  endif
+  romsgrid = strmid(namevar[xvarid], 0, 4) EQ 'lon_'
+; get the size of xaxis
+  xinq = ncdf_varinq(cdfid, xvarid)
+  ncdf_diminq, cdfid, xinq.dim[0], blabla, jpifromx
+; should we read or compute the xaxis?
+  IF NOT keyword_set(xyindex) THEN BEGIN
+; read the xaxis
+    ncdf_varget, cdfid, xvarid, xaxis
+; make sure of the shape of xaxis
+    IF xinq.ndims GE 2 THEN BEGIN
+      ncdf_diminq, cdfid, xinq.dim[1], blabla, jpjfromx
+      xaxis = reform(xaxis, jpifromx, jpjfromx, /over)
+    ENDIF
+  ENDIF ELSE xaxis = keyword_set(start1) + findgen(jpifromx)
+;----------------------------------------------------------
+; find the yaxis
+  if keyword_set(yaxisname) then yaxisname = strlowcase(yaxisname) ELSE yaxisname = 'y'
+  yvarid = (where(namevar EQ yaxisname OR namevar EQ 'latitude' $
+                  OR namevar EQ 'nav_lat' OR namevar EQ 'lat' $
+                  OR namevar EQ 'lat_rho' OR namevar EQ 'nblatitudes'))[0]
+  yvarid = yvarid[0]
+  if yvarid EQ -1 then begin
+    print, 'the yaxis was not found, check the use of YAXISNAME keyword'
+    stop
+  endif
+; get the size of yaxis and check it is ok with the values found for x
+  yinq = ncdf_varinq(cdfid, yvarid)
+  IF xinq.ndims GE 2 THEN BEGIN
+    ncdf_diminq, cdfid, yinq.dim[0], blabla, jpifromy
+    ncdf_diminq, cdfid, yinq.dim[1], blabla, jpjfromy
+    IF jpifromy NE jpifromx THEN BEGIN 
+      print, 'xaxis and y axis do not have the same x dimension...'
+    ENDIF
+  ENDIF ELSE ncdf_diminq, cdfid, yinq.dim[0], blabla, jpjfromy
+  IF n_elements(jpjfromx) NE 0 THEN BEGIN
+    IF jpjfromy NE jpjfromx THEN BEGIN 
+      print, 'xaxis and y axis do not have the same y dimension...'
+    ENDIF    
+  ENDIF
+; should we read or compute the xaxis?
+  IF NOT keyword_set(xyindex) THEN BEGIN
+; read the yaxis
+    ncdf_varget, cdfid, yvarid, yaxis
+; make sure of the shape of xaxis
+    IF xinq.ndims GE 2 THEN yaxis = reform(yaxis, jpifromy, jpjfromy, /over)
+  ENDIF ELSE yaxis = keyword_set(start1) + findgen(jpjfromy)
 ;----------------------------------------------------------
 ; find the zaxis