Changeset 1073 for trunk/LMDZ.MARS/util/concatnc.F90

Timestamp:

Oct 17, 2013, 4:19:26 PM (11 years ago)

Author:

tnavarro

Message:

Better handling of the first date of the file : Read and write the controle field, if possible, for zrecast, localtime and concatnc + cosmetic change in lslin. This is retrocompatible with previous versions.

File:

• trunk/LMDZ.MARS/util/concatnc.F90 (modified) (13 diffs)

trunk/LMDZ.MARS/util/concatnc.F90

@@ -18 +18 @@
 ! + handle absence of ap() and bp() if aps and bps are available 
 !    (case of stats file) FF, November 2011
+! + read and write controle field, if available. TN, October 2013
 ! ********************************************************
 
@@ -54 +55 @@
 ! memolonlen: # of elements of lon(), read from the first input file
 ! memoaltlen: # of elements of alt(), read from the first input file
-real, dimension(:), allocatable:: lat,lon,alt,time
+real, dimension(:), allocatable:: lat,lon,alt,ctl,time
 ! lat(): array, stores latitude coordinates
 ! lon(): array, stores longitude coordinates
 ! alt(): array, stores altitude coordinates
+! ctl(): array, stores controle coordinates
 ! time(): array, stores time coordinates
 integer :: nbvar,nbfile,nbvarfile,ndim
@@ -64 +66 @@
 ! nbvarfile: total # of variables in an input file
 ! ndim: [netcdf] # (3 or 4) of dimensions (for variables)
-integer :: latdim,londim,altdim,timedim
+integer :: latdim,londim,altdim,ctldim,timedim
 ! latdim: [netcdf] "latitude" dim ID
 ! londim: [netcdf] "longitude" dim ID
 ! altdim: [netcdf] "altdim" dim ID
+! ctldim: [netcdf] "ctldim" dim ID
 ! timedim: [netcdf] "timedim" dim ID
 integer :: gcmlayerdim ! NetCDF dimension ID for # of layers in GCM
-integer :: latvar,lonvar,altvar,timevar
+integer :: latvar,lonvar,altvar,ctlvar,timevar
 ! latvar: [netcdf] ID of "latitude" variable
 ! lonvar: [netcdf] ID of "longitude" variable
 ! altvar: [netcdf] ID of "altitude" variable
+! ctlvar: [netcdf] ID of "controle" variable
 ! timevar: [netcdf] ID of "Time" variable
-integer :: latlen,lonlen,altlen,timelen
+integer :: latlen,lonlen,altlen,ctllen,timelen
 ! latlen: # of elements of lat() array
 ! lonlen: # of elements of lon() array
-! altvar: # of elements of alt() array
+! altlen: # of elements of alt() array
+! ctllen: # of elements of ctl() array
 ! timelen: # of elemnets of time() array
 integer :: GCM_layers ! number of GCM atmospheric layers (may not be
@@ -141 +146 @@
 !write(*,*) "Beginning day of the first specified file?"
 write(*,*) "Starting day of the run stored in the first input file?"
+write(*,*) " (Obsolete if the controle field is present, answer any number)"
 write(*,*) "(e.g.: 100 if that run started at time=100 sols)"
 read(*,*) memotime
@@ -320 +326 @@
 !  write(*,*) "altlen: ",altlen
 
+   ierr=NF_INQ_DIMID(nid,"index",ctldim)
+   ierr=NF_INQ_VARID(nid,"controle",ctlvar)
+   if (ierr.NE.NF_NOERR) then
+      write(*,*) 'Field <controle> is missing in file'//file
+      ctllen=0
+      !stop ""
+   else
+      ierr=NF_INQ_DIMLEN(nid,ctldim,ctllen)
+   endif
+!  write(*,*) "controle: ",controle
+
 ! load size of aps() or sigma() (in case it is not altlen)
    ! default is that GCM_layers=altlen
@@ -345 +362 @@
       allocate(lon(lonlen))
       allocate(alt(altlen))
+      allocate(ctl(ctllen))
 #ifdef NC_DOUBLE
       ierr = NF_GET_VAR_DOUBLE(nid,latvar,lat)
       ierr = NF_GET_VAR_DOUBLE(nid,lonvar,lon)
       ierr = NF_GET_VAR_DOUBLE(nid,altvar,alt)
+      if (ctllen .ne. -1) ierr = NF_GET_VAR_DOUBLE(nid,ctlvar,ctl)
 #else
       ierr = NF_GET_VAR_REAL(nid,latvar,lat)
       ierr = NF_GET_VAR_REAL(nid,lonvar,lon)
       ierr = NF_GET_VAR_REAL(nid,altvar,alt)
-#endif
+      if (ctllen .ne. -1) ierr = NF_GET_VAR_REAL(nid,ctlvar,ctl)
+#endif
+      if (ctllen .ne. -1) then
+         if (modulo(int(memotime),669)/=modulo(int(ctl(4)),669)) then
+           write(*,'(2(A,I4),A)') "WARNING: Starting day of the run is not ",&
+                                modulo(int(memotime),669)," but ",modulo(int(ctl(4)),669),"!!"
+           write(*,*) "Starting day of the run has been corrected."
+           memotime=float(modulo(int(ctl(4)),669))
+         endif
+      endif
    ! Initialize output file's lat,lon,alt and time dimensions
-      call initiate (filename,lat,lon,alt,GCM_layers,nout,&
+      call initiate (filename,lat,lon,alt,ctl,GCM_layers,nout,&
        latdimout,londimout,altdimout,timedimout,&
        layerdimout,interlayerdimout,timevarout)
@@ -570 +598 @@
 
 !******************************************************************************
-subroutine initiate (filename,lat,lon,alt,GCM_layers,nout,&
+subroutine initiate (filename,lat,lon,alt,ctl,GCM_layers,nout,&
          latdimout,londimout,altdimout,timedimout,&
          layerdimout,interlayerdimout,timevarout)
@@ -596 +624 @@
 real, dimension(:), intent(in):: alt
 ! alt(): altitude
+real, dimension(:), intent(in):: ctl
+! ctl(): controle
 integer,intent(in) :: GCM_layers ! number of GCM layers
 integer, intent(out):: nout
@@ -619 +649 @@
 !integer :: latdim,londim,altdim,timedim
 integer :: nvarid,ierr
+integer :: ctldimout
 ! nvarid: [netcdf] ID of a variable
 ! ierr: [netcdf]  return error code (from called subroutines)
@@ -640 +671 @@
 ierr = NF_DEF_DIM(nout, "longitude", size(lon), londimout)
 ierr = NF_DEF_DIM(nout, "altitude", size(alt), altdimout)
+if (size(ctl).ne.0) ierr = NF_DEF_DIM(nout, "index", size(ctl), ctldimout)
 ierr = NF_DEF_DIM(nout, "Time", NF_UNLIMITED, timedimout)
 ierr = NF_DEF_DIM(nout, "GCM_layers", GCM_layers, layerdimout)
@@ -681 +713 @@
 
 !==============================================================================
-! 4. Write "altitude" (data and attributes)
+! 5. Write "altitude" (data and attributes)
 !==============================================================================
 
@@ -704 +736 @@
 #else
 ierr = NF_PUT_VAR_REAL (nout,nvarid,alt)
-#endif 
+#endif
+
+!==============================================================================
+! 6. Write "controle" (data and attributes)
+!==============================================================================
+
+if (size(ctl).ne.0) then
+   ! Switch to netcdf define mode
+   ierr = NF_REDEF (nout)
+
+   #ifdef NC_DOUBLE
+   ierr = NF_DEF_VAR (nout,"controle",NF_DOUBLE,1,ctldimout,nvarid)
+   #else
+   ierr = NF_DEF_VAR (nout,"controle",NF_FLOAT,1,ctldimout,nvarid)
+   #endif
+
+   ierr = NF_PUT_ATT_TEXT (nout,nvarid,"title",18,"Control parameters")
+
+   ! End netcdf define mode
+   ierr = NF_ENDDEF(nout)
+
+   #ifdef NC_DOUBLE
+   ierr = NF_PUT_VAR_DOUBLE (nout,nvarid,ctl)
+   #else
+   ierr = NF_PUT_VAR_REAL (nout,nvarid,ctl)
+   #endif
+endif
 
 end Subroutine initiate
@@ -1162 +1220 @@
 do i=1,timelen-1
    if ((ls(i+1)-ls(i)) > 350) then
-       write(*,*) "+ 360° Ls jump solved:", ls(i), ls(i+1), "at timestep", i
+       write(*,*) "+ 360 deg Ls jump solved:", ls(i), ls(i+1), "at timestep", i
       ls(i+1) = ls(i+1) - 360
        write(*,*) " corrected to now be:", ls(i), ls(i+1)
    else if ((ls(i)-ls(i+1)) > 350) then
-       write(*,*) "- 360° Ls jump solved:", ls(i), ls(i+1), "at timestep", i
+       write(*,*) "- 360 deg Ls jump solved:", ls(i), ls(i+1), "at timestep", i
       ls(i+1) = ls(i+1) + 360
        write(*,*) " corrected to now be:", ls(i), ls(i+1)