Changeset 999 for trunk/LMDZ.MARS/util/expandstartfi.F90

Timestamp:

Jul 2, 2013, 9:40:28 AM (11 years ago)

Author:

tnavarro

Message:

Possibility to store multiple initial states in one start/startfi. This is RETROCOMPATIBLE. New option ecrithist in run.def to write data in start/startfi every ecrithist dynamical timestep. New option timestart in run.def to initialize the GCM with the time timestart stored in start

File:

• trunk/LMDZ.MARS/util/expandstartfi.F90 (modified) (8 diffs)

trunk/LMDZ.MARS/util/expandstartfi.F90

@@ -20 +20 @@
 integer :: varid ! to store the ID of a variable
 integer :: datashape(4) ! to store dimension IDs of a given dataset
+integer :: corner(3),edges(3) ! to read data with a time axis
 character(len=90) :: varname ! name of a variable
 character(len=90) :: varatt ! name of attribute of a variable
@@ -28 +29 @@
 integer :: nlayer_plus_1_dimid
 integer :: number_of_advected_fields_dimid
+integer :: time_dimid
 integer :: nbindims ! number of dimensions in input file
 integer :: nbinvars ! number of variables in input file
@@ -36 +38 @@
 integer :: nlayer_plus_1
 integer :: number_of_advected_fields
+integer :: timelen
 real,allocatable :: surf_field(:) ! to store a 1D field of physical_points elements
 real,allocatable :: subsurf_field(:,:) ! to store subsurface (2D field)
@@ -155 +158 @@
 endif
 
+status=nf90_inq_dimid(inid,"Time",time_dimid)
+if (status.ne.nf90_noerr) then
+  write(*,*)"Failed to find Time dimension"
+  write(*,*)trim(nf90_strerror(status))
+  timelen = 0
+else
+  status=nf90_inquire_dimension(inid,time_dimid,len=timelen)
+  if (status.ne.nf90_noerr) then
+    write(*,*)"Failed to read Time dimension"
+    write(*,*)trim(nf90_strerror(status))
+   stop
+  else
+    write(*,*) " time length = ",timelen
+  endif
+endif
+
 ! 1.3 Allocate memory for input fields
 allocate(surf_field(physical_points))
@@ -361 +380 @@
 allocate(out_surf_field(lonlen,latlen))
 
+shape(:) = 0
 do ivar=1,nbinvars ! loop on all input variables
   ! find out what dimensions are linked to this variable
   status=nf90_inquire_variable(inid,ivar,name=varname,ndims=nbdim,&
                                dimids=shape,natts=nbatt)
-  if ((nbdim==1).and.(shape(1)==physical_points_dimid)) then
+  if (((nbdim==1).and.(shape(1)==physical_points_dimid))&
+  .or.((nbdim==2).and.(shape(1)==physical_points_dimid)&
+                 .and.(shape(2)==time_dimid))) then
+  
+    corner(1) = 1
+    corner(2) = timelen
+    edges(1)  = physical_points
+    edges(2)  = 1
     
     ! skip "longitude" and "latitude"
@@ -375 +402 @@
     ! load input data:
     status=nf90_inq_varid(inid,varname,invarid)
-    status=nf90_get_var(inid,invarid,surf_field)
+    status=nf90_get_var(inid,invarid,surf_field,corner,edges)
     
     ! switch output file to to define mode
@@ -455 +482 @@
   status=nf90_inquire_variable(inid,ivar,name=varname,ndims=nbdim,&
                                dimids=shape,natts=nbatt)
-  if ((nbdim==2).and.(shape(1)==physical_points_dimid) &
-                .and.(shape(2)==subsurface_layers_dimid)) then
+  if (((nbdim==2).and.(shape(1)==physical_points_dimid)&
+                 .and.(shape(2)==subsurface_layers_dimid))&
+  .or.((nbdim==3).and.(shape(1)==physical_points_dimid)&
+                 .and.(shape(2)==subsurface_layers_dimid)&
+                 .and.(shape(3)==time_dimid))) then
+    
+    corner(1) = 1
+    corner(2) = 1
+    corner(3) = timelen
+    edges(1)  = physical_points
+    edges(2)  = subsurface_layers
+    edges(3)  = 1
     
     write(*,*) " processing: ",trim(varname)
@@ -462 +499 @@
     ! load input data:
     status=nf90_inq_varid(inid,varname,invarid)
-    status=nf90_get_var(inid,invarid,subsurf_field)
+    status=nf90_get_var(inid,invarid,subsurf_field,corner,edges)
     
     ! switch output file to to define mode