MODULE flincom !- !$Id: flincom.f90 1932 2012-11-28 09:56:17Z jgipsl $ !- ! This software is governed by the CeCILL license ! See IOIPSL/IOIPSL_License_CeCILL.txt !--------------------------------------------------------------------- USE netcdf !- USE calendar, ONLY : ju2ymds, ymds2ju, ioconf_calendar USE errioipsl, ONLY : histerr, ipslout,ipslerr,ipsldbg USE stringop, ONLY : strlowercase !- IMPLICIT NONE !- PRIVATE PUBLIC :: flinput, flincre, flinget, flinclo, & flinopen, flininfo, flininspect, flinquery_var !- INTERFACE flinopen !--------------------------------------------------------------------- !- The "flinopen" routines will open an input file !- !- INPUT !- !- filename : Name of the netCDF file to be opened !- !- iideb : index i for zoom ! !- iilen : length of zoom ! for !- jjdeb : index j for zoom ! zoom !- jjlen : length of zoom ! !- !- do_test : A flag that enables the testing of the content !- of the file against the input from the model !- !- INPUT if do_test=TRUE OUTPUT else !- !- iim : size in the x direction in the file (longitude) !- jjm : size in the y direction !- llm : number of levels !- (llm = 0 means no axis to be expected) !- lon : array of (iilen,jjlen) (zoom), or (iim,jjm) (no zoom), !- that contains the longitude of each point !- lat : same for latitude !- lev : An array of llm for the latitude !- !- WARNING : !- In the case of do_test=FALSE it is for the user to check !- that the dimensions of lon lat and lev are correct when passed to !- flinopen. This can be done after the call when iim and jjm have !- been retrieved from the netCDF file. In F90 this problem will !- be solved with an internal assign !- IF iim, jjm, llm or ttm are parameters in the calling program and !- you use the option do_test=FALSE it will create a segmentation fault !- !- OUTPUT !- !- ttm : size of time axis !- itaus : Time steps within this file !- date0 : Julian date at which itau = 0 !- dt : length of the time steps of the data !- fid : returned file ID which is later used to read the data !--------------------------------------------------------------------- MODULE PROCEDURE flinopen_zoom2d, flinopen_nozoom END INTERFACE !- INTERFACE flinput !--------------------------------------------------------------------- !- The "flinput" routines will put a variable !- on the netCDF file created by flincre. !- If the sizes of the axis do not match the one of the IDs !- then a new axis is created. !- That is we loose the possibility of writting hyperslabs of data. !- !- Again here if iim = jjm = llm = ttm = 0 !- then a global attribute is added to the file. !- !- INPUT !- !- fid : Identification of the file in which we will write !- varname : Name of variable to be written !- iim : size in x of variable !- nlonid : ID of x axis which could fit for this axis !- jjm : size in y of variable !- nlatid : ID of y axis which could fit for this axis !- llm : size in z of variable !- zdimid : ID of z axis which could fit for this axis !- ttm : size in t of variable !- tdimid : ID of t axis which could fit for this axis !- !- OUTPUT !- !- NONE !--------------------------------------------------------------------- MODULE PROCEDURE flinput_r4d, flinput_r3d, flinput_r2d, & flinput_r1d, flinput_scal END INTERFACE !- INTERFACE flinget MODULE PROCEDURE flinget_r4d, flinget_r3d, flinget_r2d, & flinget_r1d, flinget_scal, & flinget_r4d_zoom2d, flinget_r3d_zoom2d, & flinget_r2d_zoom2d END INTERFACE !- ! This is the data we keep on each file we open !- INTEGER, PARAMETER :: nbfile_max = 200 INTEGER, SAVE :: nbfiles = 0 INTEGER, SAVE :: ncids(nbfile_max), ncnbd(nbfile_max), & ncfunli(nbfile_max), ncnba(nbfile_max) INTEGER, SAVE :: ncnbva(nbfile_max), ncdims(nbfile_max,4) LOGICAL, SAVE :: ncfileopen(nbfile_max)=.FALSE. !- INTEGER, SAVE :: cind_vid, cind_fid, cind_len INTEGER,DIMENSION(:),ALLOCATABLE,SAVE :: cindex !- INTEGER,DIMENSION(4) :: w_sta, w_len, w_dim !- CONTAINS !- !=== !- SUBROUTINE flincre & (filename, iim1, jjm1, lon1, lat1, llm1, lev1, ttm1, itaus, & time0, dt, fid_out, nlonid1, nlatid1, zdimid1, tdimid1) !--------------------------------------------------------------------- !- This is a "low level" subroutine for opening netCDF files wich !- contain the major coordinate system of the model. !- Other coordinates needed for other variables !- will be added as they are needed. !- !- INPUT !- !- filename : Name of the file to be created !- iim1, jjm1 : Horizontal size of the grid !- which will be stored in the file !- lon1, lat1 : Horizontal grids !- llm1 : Size of the vertical grid !- lev1 : Vertical grid !- ttm1 : Size of time axis !- itaus : time steps on the time axis !- time0 : Time in julian days at which itau = 0 !- dt : time step in seconds between itaus !- (one step of itau) !- !- OUTPUT !- !- fid : File identification !- nlonid1 : Identification of longitudinal axis !- nlatid1 : Identification of latitudinal axis !- zdimid1 : ID of vertical axis !- tdimid1 : ID of time axis !--------------------------------------------------------------------- IMPLICIT NONE !- ! ARGUMENTS !- CHARACTER(LEN=*) :: filename INTEGER :: iim1, jjm1, llm1, ttm1 REAL :: lon1(iim1,jjm1) REAL :: lat1(iim1,jjm1) REAL :: lev1(llm1) INTEGER :: itaus(ttm1) REAL :: time0 REAL :: dt INTEGER :: fid_out, zdimid1, nlonid1, nlatid1, tdimid1 !- ! LOCAL !- INTEGER :: iret, lll, fid INTEGER :: lonid, latid, levid, timeid INTEGER :: year, month, day REAL :: sec CHARACTER(LEN=250):: name !- LOGICAL :: l_dbg !--------------------------------------------------------------------- CALL ipsldbg (old_status=l_dbg) lll = LEN_TRIM(filename) IF (filename(lll-2:lll) /= '.nc') THEN name=filename(1:lll)//'.nc' ELSE name=filename(1:lll) ENDIF !- iret = NF90_CREATE (name, NF90_CLOBBER, fid) !- iret = NF90_DEF_DIM (fid, 'x', iim1, nlonid1) iret = NF90_DEF_DIM (fid, 'y', jjm1, nlatid1) iret = NF90_DEF_DIM (fid, 'lev', llm1, zdimid1) iret = NF90_DEF_DIM (fid, 'tstep', ttm1, tdimid1) !- ! Vertical axis !- IF (l_dbg) WRITE(ipslout,*) 'flincre Vertical axis' !- iret = NF90_DEF_VAR (fid, 'lev', NF90_FLOAT, zdimid1, levid) iret = NF90_PUT_ATT (fid, levid, 'units', '-') iret = NF90_PUT_ATT (fid, levid, 'title', 'levels') iret = NF90_PUT_ATT (fid, levid, 'long_name', 'Sigma Levels') !- ! Time axis !- IF (l_dbg) WRITE(ipslout,*) 'flincre time axis' !- iret = NF90_DEF_VAR (fid, 'tstep', NF90_FLOAT, tdimid1, timeid) iret = NF90_PUT_ATT (fid, timeid, 'units', '-') iret = NF90_PUT_ATT (fid, timeid, 'title', 'time') iret = NF90_PUT_ATT (fid, timeid, 'long_name', 'time steps') !- ! The longitude !- IF (l_dbg) WRITE(ipslout,*) 'flincre Longitude axis' !- iret = NF90_DEF_VAR (fid, "nav_lon", NF90_FLOAT, & (/ nlonid1, nlatid1 /), lonid) iret = NF90_PUT_ATT (fid, lonid, 'units', "degrees_east") iret = NF90_PUT_ATT (fid, lonid, 'title', "Longitude") iret = NF90_PUT_ATT (fid, lonid, 'nav_model', & "Lambert projection of PROMES") iret = NF90_PUT_ATT (fid, lonid, 'valid_min', & REAL(MINVAL(lon1),KIND=4)) iret = NF90_PUT_ATT (fid, lonid, 'valid_max', & REAL(MAXVAL(lon1),KIND=4)) !- ! The Latitude !- IF (l_dbg) WRITE(ipslout,*) 'flincre Latitude axis' !- iret = NF90_DEF_VAR (fid, "nav_lat", NF90_FLOAT, & (/ nlonid1, nlatid1 /), latid) iret = NF90_PUT_ATT (fid, latid, 'units', "degrees_north") iret = NF90_PUT_ATT (fid, latid, 'title', "Latitude") iret = NF90_PUT_ATT (fid, latid, 'nav_model', & "Lambert projection of PROMES") iret = NF90_PUT_ATT (fid, latid, 'valid_min', & REAL(MINVAL(lat1),KIND=4)) iret = NF90_PUT_ATT (fid, latid, 'valid_max', & REAL(MAXVAL(lat1),KIND=4)) !- ! The time coordinates !- iret = NF90_PUT_ATT (fid, NF90_GLOBAL, 'delta_tstep_sec', & REAL(dt,KIND=4)) !- CALL ju2ymds (time0, year, month, day, sec) !- iret = NF90_PUT_ATT (fid, NF90_GLOBAL, 'year0', REAL(year,KIND=4)) iret = NF90_PUT_ATT (fid, NF90_GLOBAL, 'month0', REAL(month,KIND=4)) iret = NF90_PUT_ATT (fid, NF90_GLOBAL, 'day0', REAL(day,KIND=4)) iret = NF90_PUT_ATT (fid, NF90_GLOBAL, 'sec0', REAL(sec,KIND=4)) !- iret = NF90_ENDDEF (fid) !- IF (l_dbg) WRITE(ipslout,*) 'flincre Variable' !- iret = NF90_PUT_VAR (fid, levid, lev1(1:llm1)) !- IF (l_dbg) WRITE(ipslout,*) 'flincre Time Variable' !- iret = NF90_PUT_VAR (fid, timeid, REAL(itaus(1:ttm1))) !- IF (l_dbg) WRITE(ipslout,*) 'flincre Longitude' !- iret = NF90_PUT_VAR (fid, lonid, lon1(1:iim1,1:jjm1)) !- IF (l_dbg) WRITE(ipslout,*) 'flincre Latitude' !- iret = NF90_PUT_VAR (fid, latid, lat1(1:iim1,1:jjm1)) !- ! Keep all this information !- nbfiles = nbfiles+1 !- IF (nbfiles > nbfile_max) THEN CALL histerr (3,'flincre', & 'Too many files. Please increase nbfil_max', & 'in program flincom.F90.',' ') ENDIF !- ncids(nbfiles) = fid ncnbd(nbfiles) = 4 !- ncdims(nbfiles,1:4) = (/ iim1, jjm1, llm1, ttm1 /) !- ncfunli(nbfiles) = -1 ncnba(nbfiles) = 4 ncnbva(nbfiles) = 0 ncfileopen(nbfiles) = .TRUE. !- fid_out = nbfiles !--------------------- END SUBROUTINE flincre !- !=== !- SUBROUTINE flinopen_zoom2d & (filename, iideb, iilen, jjdeb, jjlen, do_test, & iim, jjm, llm, lon, lat, lev, ttm, itaus, date0, dt, fid_out) !--------------------------------------------------------------------- IMPLICIT NONE !- ! ARGUMENTS !- CHARACTER(LEN=*) :: filename LOGICAL :: do_test INTEGER :: iim, jjm, llm, ttm, iideb, iilen, jjdeb, jjlen REAL :: lon(iilen,jjlen), lat(iilen,jjlen), lev(llm) INTEGER :: itaus(ttm) REAL :: date0, dt INTEGER :: fid_out !- LOGICAL :: l_dbg !--------------------------------------------------------------------- CALL ipsldbg (old_status=l_dbg) IF (l_dbg) WRITE (*,*) ' iideb, iilen, jjdeb, jjlen, iim, jjm ', & iideb, iilen, jjdeb, jjlen, iim, jjm IF (l_dbg) WRITE (*,*) ' lon ', lon(1,1), lon(iilen,jjlen) IF (l_dbg) WRITE (*,*) ' lat ', lat(1,1), lat(iilen,jjlen) !- CALL flinopen_work & (filename, iideb, iilen, jjdeb, jjlen, do_test, & iim, jjm, llm, lon, lat, lev, ttm, itaus, date0, dt, fid_out) !----------------------------- END SUBROUTINE flinopen_zoom2d !- !=== !- SUBROUTINE flinopen_nozoom & (filename, do_test, iim, jjm, llm, lon, lat, lev, ttm, & itaus, date0, dt, fid_out) !--------------------------------------------------------------------- IMPLICIT NONE !- ! ARGUMENTS !- CHARACTER(LEN=*) :: filename LOGICAL :: do_test INTEGER :: iim, jjm, llm, ttm REAL :: lon(iim,jjm), lat(iim,jjm), lev(llm) INTEGER :: itaus(ttm) REAL :: date0, dt INTEGER :: fid_out INTEGER :: iimc, jjmc !--------------------------------------------------------------------- iimc=iim jjmc=jjm CALL flinopen_work & (filename, 1, iimc, 1, jjmc, do_test, & iim, jjm, llm, lon, lat, lev, ttm, itaus, date0, dt, fid_out) !------------------------- END SUBROUTINE flinopen_nozoom !- !=== !- SUBROUTINE flinopen_work & (filename, iideb, iilen, jjdeb, jjlen, do_test, & iim, jjm, llm, lon, lat, lev, ttm, itaus, date0, dt, fid_out) !--------------------------------------------------------------------- IMPLICIT NONE !- ! ARGUMENTS !- CHARACTER(LEN=*) :: filename LOGICAL :: do_test INTEGER :: iim, jjm, llm, ttm, iideb, iilen, jjdeb, jjlen REAL :: lon(iilen,jjlen), lat(iilen,jjlen), lev(llm) INTEGER :: itaus(ttm) REAL :: date0, dt INTEGER :: fid_out !- ! LOCAL !- REAL, PARAMETER :: eps = 1.e-4 !- INTEGER :: iret, vid, fid, nbdim, i, iilast, jjlast INTEGER :: gdtt_id, old_id, iv, gdtmaf_id CHARACTER(LEN=250) :: name CHARACTER(LEN=80) :: units, calendar INTEGER :: tmp_iim, tmp_jjm, tmp_llm, tmp_ttm REAL :: x_first, x_last INTEGER :: year, month, day REAL :: r_year, r_month, r_day INTEGER :: year0, month0, day0, hours0, minutes0, seci REAL :: sec, sec0 CHARACTER :: strc !- REAL,DIMENSION(:),ALLOCATABLE :: vec_tmp !- LOGICAL :: open_file LOGICAL :: l_dbg !--------------------------------------------------------------------- CALL ipsldbg (old_status=l_dbg) iilast = iideb+iilen-1 jjlast = jjdeb+jjlen-1 IF (l_dbg) WRITE (*,*) & ' flinopen_work zoom 2D information '// & ' iideb, iilen, iilast, jjdeb, jjlen, jjlast ', & iideb, iilen, iilast, jjdeb, jjlen, jjlast !- ! 1.0 get all infos on the file !- ! Either the fid_out has not been initialized (0 or very large) ! then we have to open anyway. Else we only need to open the file ! if it has not been opened before. !- IF ( (fid_out < 1).OR.(fid_out > nbfile_max) ) THEN open_file = .TRUE. ELSE IF (.NOT.ncfileopen(fid_out)) THEN open_file = .TRUE. ELSE open_file = .FALSE. ENDIF !- IF (open_file) THEN CALL flininfo (filename,tmp_iim,tmp_jjm,tmp_llm,tmp_ttm,fid_out) ELSE !-- The user has already opened the file !-- and we trust that he knows the dimensions tmp_iim = iim tmp_jjm = jjm tmp_llm = llm tmp_ttm = ttm ENDIF !- IF (l_dbg) & WRITE(ipslout,*) 'OUT OF flininfo :',tmp_iim,tmp_jjm,tmp_llm,tmp_ttm !- fid = ncids(fid_out) !- ! 2.0 get the sizes and names of the different coordinates ! and do a first set of verification. !- ! 2.2 We test the axis if we have to. !- IF (l_dbg) & WRITE(ipslout,*) 'flininfo 2.2 We test if we have to test : ',do_test !- IF (do_test) THEN IF (iim /= tmp_iim) THEN CALL histerr (3,'flinopen', & 'file '//filename//' does not have the ', & 'required dimension in x direction (longitude)',' ') ELSE IF (jjm /= tmp_jjm) THEN CALL histerr (3,'flinopen', & 'file '//filename//' does not have the ', & 'required dimension in y direction (latitude)',' ') ELSE IF ( llm /= tmp_llm .AND. llm > 0 ) THEN CALL histerr (3,'flinopen', & 'file '//filename//' does not have the ', & 'required dimension in the vertical',' ') ENDIF ELSE !--- !-- 2.3 Else the sizes of the axes are returned to the user !--- IF (l_dbg) WRITE(ipslout,*) 'flinopen 2.3 Else sizes are returned' !--- iim = tmp_iim jjm = tmp_jjm llm = tmp_llm ENDIF !- ttm = tmp_ttm !- ! 3.0 Check if we are realy talking about the same coodinate system ! if not then we get the lon, lat and lev variables from the file !- IF (l_dbg) WRITE(ipslout,*) 'flinopen 3.0 we are realy talking' !- IF (do_test) THEN !--- CALL flinfindcood (fid_out, 'lon', vid, nbdim) iret = NF90_GET_VAR (fid, vid, x_first, start=(/ iideb, jjdeb /)) iret = NF90_GET_VAR (fid, vid, x_last, start=(/ iilast, jjlast /)) !--- IF (l_dbg) & WRITE(ipslout,*) 'from file lon first and last, modulo 360. ', & x_first, x_last, MODULO(x_first,360.), MODULO(x_last,360.) IF (l_dbg) & WRITE(ipslout,*) 'from model lon first and last, modulo 360. ', & lon(1,1),lon(iilen,jjlen), & MODULO(lon(1,1),360.), MODULO(lon(iilen,jjlen),360.) IF ( (ABS( MODULO(x_first,360.) & -MODULO(lon(1,1),360.)) > eps) & .OR.(ABS( MODULO(x_last,360.) & -MODULO(lon(iilen ,jjlen),360.)) > eps ) ) THEN CALL histerr (3,'flinopen', & 'file '//filename//' and the model do not', & 'share the same longitude coordinate', & 'Obtained by comparing the first and last values ') ENDIF !--- CALL flinfindcood (fid_out, 'lat', vid, nbdim) iret = NF90_GET_VAR (fid, vid, x_first, start=(/ iideb, jjdeb /)) iret = NF90_GET_VAR (fid, vid, x_last, start=(/ iilast, jjlast /)) !--- IF (l_dbg) WRITE(ipslout,*) & 'from file lat first and last ',x_first,x_last IF (l_dbg) WRITE(ipslout,*) & 'from model lat first and last ',lat(1,1),lat(iilen,jjlen) !--- IF ( (ABS(x_first-lat(1,1)) > eps) & .OR.(ABS(x_last-lat(iilen,jjlen)) > eps) ) THEN CALL histerr (3,'flinopen', & 'file '//filename//' and the model do not', & 'share the same latitude coordinate', & 'Obtained by comparing the first and last values ') ENDIF !--- IF (llm > 0) THEN CALL flinfindcood (fid_out, 'lev', vid, nbdim) iret = NF90_GET_VAR (fid, vid, x_first, start=(/ 1 /)) iret = NF90_GET_VAR (fid, vid, x_last, start=(/ llm /)) !----- IF (l_dbg) WRITE(ipslout,*) & 'from file lev first and last ',x_first ,x_last IF (l_dbg) WRITE(ipslout,*) & 'from model lev first and last ',lev(1),lev(llm) !----- IF ( (ABS(x_first-lev(1)) > eps) & .OR.(ABS(x_last-lev(llm)) > eps) ) THEN CALL histerr (3,'flinopen', & 'file '//filename//' and the model do not', & 'share the same vertical coordinate', & 'Obtained by comparing the first and last values') ENDIF ENDIF !--- ELSE !--- !-- 4.0 extracting the coordinates if we do not check !--- IF (l_dbg) WRITE(ipslout,*) 'flinopen 4.0 extracting the coordinates' !--- CALL flinfindcood (fid_out, 'lon', vid, nbdim) IF (nbdim == 2) THEN iret = NF90_GET_VAR (fid, vid, lon, & start=(/ iideb, jjdeb /), count=(/ iilen, jjlen /)) ELSE ALLOCATE(vec_tmp(iilen)) iret = NF90_GET_VAR (fid, vid, vec_tmp, & start=(/ iideb /), count=(/ iilen /)) DO i=1,jjlen lon(:,i) = vec_tmp(:) ENDDO DEALLOCATE(vec_tmp) ENDIF !--- CALL flinfindcood (fid_out, 'lat', vid, nbdim) IF (nbdim == 2) THEN iret = NF90_GET_VAR (fid, vid, lat, & start=(/ iideb, jjdeb /), count=(/ iilen, jjlen /)) ELSE ALLOCATE(vec_tmp(jjlen)) iret = NF90_GET_VAR (fid, vid, vec_tmp, & start=(/ jjdeb /), count=(/ jjlen /)) DO i=1,iilen lat(i,:) = vec_tmp(:) ENDDO DEALLOCATE(vec_tmp) ENDIF !--- IF (llm > 0) THEN CALL flinfindcood (fid_out, 'lev', vid, nbdim) IF (nbdim == 1) THEN iret = NF90_GET_VAR (fid, vid, lev, & start=(/ 1 /), count=(/ llm /)) ELSE CALL histerr (3,'flinopen', & 'Can not handle vertical coordinates that have more',& 'than 1 dimension',' ') ENDIF ENDIF ENDIF !- ! 5.0 Get all the details for the time if possible needed !- IF (l_dbg) WRITE(ipslout,*) 'flinopen 5.0 Get time' !- IF (ttm > 0) THEN !--- !-- 5.1 Find the time axis. Prefered method is the 'timestep since' !--- gdtmaf_id = -1 gdtt_id = -1 old_id = -1 DO iv=1,ncnbva(fid_out) name='' iret = NF90_INQUIRE_VARIABLE (fid, iv, name=name) units='' iret = NF90_GET_ATT (fid, iv, 'units', units) IF (INDEX(units,'seconds since') > 0) gdtmaf_id = iv IF (INDEX(units,'timesteps since') > 0) gdtt_id = iv IF (INDEX(name, 'tstep') > 0 .OR. INDEX(name,'time') > 0 ) old_id = iv ENDDO !--- IF (gdtt_id > 0) THEN vid = gdtt_id ELSE IF (gdtmaf_id > 0) THEN vid = gdtmaf_id ELSE IF (old_id > 0) THEN vid = old_id ELSE CALL histerr (3, 'flinopen', 'No time axis found',' ',' ') ENDIF !--- ALLOCATE(vec_tmp(ttm)) iret = NF90_GET_VAR (fid,vid,vec_tmp,start=(/ 1 /),count=(/ ttm /)) itaus(1:ttm) = NINT(vec_tmp(1:ttm)) DEALLOCATE(vec_tmp) !--- IF (l_dbg) WRITE(ipslout,*) 'flinopen 5.1 Times ',itaus !--- !-- Getting all the details for the time axis !--- !-- Find the calendar calendar = '' iret = NF90_GET_ATT (fid,gdtmaf_id,'calendar',calendar) IF (iret == NF90_NOERR) THEN CALL ioconf_calendar(calendar) ENDIF !-- units = '' iret = NF90_GET_ATT (fid,vid,'units',units) IF (gdtt_id > 0) THEN units = units(INDEX(units,'since')+6:LEN_TRIM(units)) READ (units,'(I4.4,5(a,I2.2))') & year0, strc, month0, strc, day0, & strc, hours0, strc, minutes0, strc, seci sec0 = hours0*3600. + minutes0*60. + seci CALL ymds2ju (year0, month0, day0, sec0, date0) IF (l_dbg) & WRITE(ipslout,*) 'flinopen 5.1 gdtt_id year0 ... date0 ', & year0, month0, day0, sec0, date0 !----- iret = NF90_GET_ATT (fid, gdtt_id, 'tstep_sec', dt) ELSE IF (gdtmaf_id > 0) THEN units = units(INDEX(units,'since')+6:LEN_TRIM(units)) READ (units,'(I4.4,5(a,I2.2))') & year0, strc, month0, strc, day0, & strc, hours0, strc, minutes0, strc, seci sec0 = hours0*3600. + minutes0*60. + seci CALL ymds2ju (year0, month0, day0, sec0, date0) !----- IF (l_dbg) & WRITE(ipslout,*) 'flinopen 5.1 gdtmaf_id year0 ... date0 ', & year0, month0, day0, sec0, date0 ELSE IF (old_id > 0) THEN iret = NF90_GET_ATT (fid, NF90_GLOBAL, 'delta_tstep_sec', dt) iret = NF90_GET_ATT (fid, NF90_GLOBAL, 'day0', r_day) iret = NF90_GET_ATT (fid, NF90_GLOBAL, 'sec0', sec) iret = NF90_GET_ATT (fid, NF90_GLOBAL, 'year0', r_year) iret = NF90_GET_ATT (fid, NF90_GLOBAL, 'month0', r_month) !----- day = INT(r_day) month = INT(r_month) year = INT(r_year) !----- CALL ymds2ju (year, month, day, sec, date0) ENDIF ENDIF !- IF (l_dbg) WRITE(ipslout,*) 'flinopen 6.0 File opened', date0, dt !--------------------------- END SUBROUTINE flinopen_work !- !=== !- SUBROUTINE flininfo (filename, iim, jjm, llm, ttm, fid_out) !--------------------------------------------------------------------- !- This subroutine allows to get some information. !- It is usualy done within flinopen but the user may want to call !- it before in order to allocate the space needed to extract the !- data from the file. !--------------------------------------------------------------------- IMPLICIT NONE !- ! ARGUMENTS !- CHARACTER(LEN=*) :: filename INTEGER :: iim, jjm, llm, ttm, fid_out !- ! LOCAL !- INTEGER :: iret, fid, ndims, nvars, nb_atts, id_unlim INTEGER :: iv, lll INTEGER :: xid, yid, zid, tid CHARACTER(LEN=80) :: name CHARACTER(LEN=30) :: axname !- LOGICAL :: l_dbg !--------------------------------------------------------------------- CALL ipsldbg (old_status=l_dbg) lll = LEN_TRIM(filename) IF (filename(lll-2:lll) /= '.nc') THEN name = filename(1:lll)//'.nc' ELSE name = filename(1:lll) ENDIF !- iret = NF90_OPEN (name, NF90_NOWRITE, fid) IF (iret /= NF90_NOERR) THEN CALL histerr(3, 'flininfo','Could not open file :',TRIM(name),' ') ENDIF !- iret = NF90_INQUIRE (fid, nDimensions=ndims, nVariables=nvars, & nAttributes=nb_atts, unlimitedDimId=id_unlim) !- xid = -1; iim = 0; yid = -1; jjm = 0; zid = -1; llm = 0; tid = -1; ttm = 0; !- DO iv=1,ndims !--- iret = NF90_INQUIRE_DIMENSION (fid, iv, name=axname, len=lll) CALL strlowercase (axname) axname = ADJUSTL(axname) !--- IF (l_dbg) WRITE(ipslout,*) & 'flininfo - getting axname',iv,axname,lll !--- IF ( (INDEX(axname,'x') == 1) & .OR.(INDEX(axname,'lon') == 1) ) THEN xid = iv; iim = lll; ELSE IF ( (INDEX(axname,'y') == 1) & .OR.(INDEX(axname,'lat') == 1) ) THEN yid = iv; jjm = lll; ELSE IF ( (INDEX(axname,'lev') == 1) & .OR.(INDEX(axname,'plev') == 1) & .OR.(INDEX(axname,'z') == 1) & .OR.(INDEX(axname,'depth') == 1) ) THEN zid = iv; llm = lll; ELSE IF ( (INDEX(axname,'tstep') == 1) & .OR.(INDEX(axname,'time') == 1) & .OR.(INDEX(axname,'time_counter') == 1) ) THEN !---- For the time we certainly need to allow for other names tid = iv; ttm = lll; ELSE IF (ndims == 1) THEN !---- Nothing was found and ndims=1 then we have a vector of data xid = 1; iim = lll; ENDIF !--- ENDDO !- ! Keep all this information !- nbfiles = nbfiles+1 !- IF (nbfiles > nbfile_max) THEN CALL histerr (3,'flininfo', & 'Too many files. Please increase nbfil_max', & 'in program flincom.F90.',' ') ENDIF !- ncids(nbfiles) = fid ncnbd(nbfiles) = ndims !- ncdims(nbfiles,1:4) = (/ iim, jjm, llm, ttm /) !- ncfunli(nbfiles) = id_unlim ncnba(nbfiles) = nb_atts ncnbva(nbfiles) = nvars ncfileopen(nbfiles) = .TRUE. !- fid_out = nbfiles !---------------------- END SUBROUTINE flininfo !- !=== !- SUBROUTINE flinput_r1d & (fid_in,varname,iim,nlonid,jjm,nlatid,llm,zdimid,ttm,tdimid,var) !--------------------------------------------------------------------- IMPLICIT NONE !- INTEGER :: fid_in CHARACTER(LEN=*) :: varname INTEGER :: iim, nlonid, jjm, nlatid, llm, zdimid, ttm, tdimid REAL :: var(:) !- INTEGER :: fid, ncvarid, ndim, iret LOGICAL :: l_dbg !--------------------------------------------------------------------- CALL ipsldbg (old_status=l_dbg) IF (l_dbg) WRITE(ipslout,*) & "flinput_r1d : SIZE(var) = ",SIZE(var) !- CALL flinput_mat & (fid_in,varname,iim,nlonid,jjm,nlatid,llm,zdimid,ttm,tdimid, & fid,ncvarid,ndim) !- iret = NF90_PUT_VAR (fid, ncvarid, var, & start=w_sta(1:ndim), count=w_len(1:ndim)) !------------------------- END SUBROUTINE flinput_r1d !- !=== !- SUBROUTINE flinput_r2d & (fid_in,varname,iim,nlonid,jjm,nlatid,llm,zdimid,ttm,tdimid,var) !--------------------------------------------------------------------- IMPLICIT NONE !- INTEGER :: fid_in CHARACTER(LEN=*) :: varname INTEGER :: iim, nlonid, jjm, nlatid, llm, zdimid, ttm, tdimid REAL :: var(:,:) !- INTEGER :: fid, ncvarid, ndim, iret LOGICAL :: l_dbg !--------------------------------------------------------------------- CALL ipsldbg (old_status=l_dbg) IF (l_dbg) WRITE(ipslout,*) & "flinput_r2d : SIZE(var) = ",SIZE(var) !- CALL flinput_mat & (fid_in,varname,iim,nlonid,jjm,nlatid,llm,zdimid,ttm,tdimid, & fid,ncvarid,ndim) !- iret = NF90_PUT_VAR (fid, ncvarid, var, & start=w_sta(1:ndim), count=w_len(1:ndim)) !------------------------- END SUBROUTINE flinput_r2d !- !=== !- SUBROUTINE flinput_r3d & (fid_in,varname,iim,nlonid,jjm,nlatid,llm,zdimid,ttm,tdimid,var) !--------------------------------------------------------------------- IMPLICIT NONE !- INTEGER :: fid_in CHARACTER(LEN=*) :: varname INTEGER :: iim, nlonid, jjm, nlatid, llm, zdimid, ttm, tdimid REAL :: var(:,:,:) !- INTEGER :: fid, ncvarid, ndim, iret LOGICAL :: l_dbg !--------------------------------------------------------------------- CALL ipsldbg (old_status=l_dbg) IF (l_dbg) WRITE(ipslout,*) & "flinput_r3d : SIZE(var) = ",SIZE(var) !- CALL flinput_mat & (fid_in,varname,iim,nlonid,jjm,nlatid,llm,zdimid,ttm,tdimid, & fid,ncvarid,ndim) !- iret = NF90_PUT_VAR (fid, ncvarid, var, & start=w_sta(1:ndim), count=w_len(1:ndim)) !------------------------- END SUBROUTINE flinput_r3d !- !=== !- SUBROUTINE flinput_r4d & (fid_in,varname,iim,nlonid,jjm,nlatid,llm,zdimid,ttm,tdimid,var) !--------------------------------------------------------------------- IMPLICIT NONE !- INTEGER :: fid_in CHARACTER(LEN=*) :: varname INTEGER :: iim, nlonid, jjm, nlatid, llm, zdimid, ttm, tdimid REAL :: var(:,:,:,:) !- INTEGER :: fid, ncvarid, ndim, iret LOGICAL :: l_dbg !--------------------------------------------------------------------- CALL ipsldbg (old_status=l_dbg) IF (l_dbg) WRITE(ipslout,*) & "flinput_r4d : SIZE(var) = ",SIZE(var) !- CALL flinput_mat & (fid_in,varname,iim,nlonid,jjm,nlatid,llm,zdimid,ttm,tdimid, & fid,ncvarid,ndim) !- iret = NF90_PUT_VAR (fid, ncvarid, var, & start=w_sta(1:ndim), count=w_len(1:ndim)) !------------------------- END SUBROUTINE flinput_r4d !- !=== !- SUBROUTINE flinput_mat & (fid_in,varname,iim,nlonid,jjm,nlatid, & llm,zdimid,ttm,tdimid,fid,ncvarid,ndim) !--------------------------------------------------------------------- IMPLICIT NONE !- INTEGER :: fid_in CHARACTER(LEN=*) :: varname INTEGER :: iim, nlonid, jjm, nlatid, llm, zdimid, ttm, tdimid INTEGER :: fid, ncvarid, ndim !- ! LOCAL !- INTEGER :: iret !--------------------------------------------------------------------- fid = ncids(fid_in) !- w_sta(1:4) = (/ 1, 1, 1, 1 /) w_len(1:2) = (/ iim, jjm /) w_dim(1:2) = (/ nlonid, nlatid /) !- IF ( (llm > 0).AND.(ttm > 0) ) THEN ndim = 4 w_len(3:4) = (/ llm, ttm /) w_dim(3:4) = (/ zdimid, tdimid /) ELSE IF (llm > 0) THEN ndim = 3 w_dim(3) = zdimid w_len(3) = llm ELSE IF (ttm > 0) THEN ndim = 3 w_dim(3) = tdimid w_len(3) = ttm ELSE ndim = 2 ENDIF !- iret = NF90_REDEF (fid) iret = NF90_DEF_VAR (fid,varname,NF90_FLOAT,w_dim(1:ndim),ncvarid) iret = NF90_PUT_ATT (fid,ncvarid,'short_name',TRIM(varname)) iret = NF90_ENDDEF (fid) !-------------------------- END SUBROUTINE flinput_mat !- !=== !- SUBROUTINE flinput_scal & (fid_in, varname, iim, nlonid, jjm, nlatid, & llm, zdimid, ttm, tdimid, var) !--------------------------------------------------------------------- IMPLICIT NONE !- INTEGER :: fid_in CHARACTER(LEN=*) :: varname INTEGER :: iim, nlonid, jjm, nlatid, llm, zdimid, ttm, tdimid REAL :: var !- ! LOCAL !- INTEGER :: fid, iret !--------------------------------------------------------------------- fid = ncids(fid_in) !- iret = NF90_REDEF (fid) iret = NF90_PUT_ATT (fid, NF90_GLOBAL, varname, REAL(var,KIND=4)) iret = NF90_ENDDEF (fid) !--------------------------- END SUBROUTINE flinput_scal !- !=== !- SUBROUTINE flinget_r1d & (fid_in,varname,iim,jjm,llm,ttm,itau_dep,itau_fin,var) !--------------------------------------------------------------------- IMPLICIT NONE !- INTEGER :: fid_in CHARACTER(LEN=*) :: varname INTEGER :: iim, jjm, llm, ttm, itau_dep, itau_fin REAL :: var(:) !- INTEGER :: jl, ji REAL,DIMENSION(:),ALLOCATABLE,SAVE :: buff_tmp LOGICAL :: l_dbg !--------------------------------------------------------------------- CALL ipsldbg (old_status=l_dbg) IF (.NOT.ALLOCATED(buff_tmp)) THEN IF (l_dbg) WRITE(ipslout,*) & "flinget_r1d : allocate buff_tmp for buff_sz = ",SIZE(var) ALLOCATE (buff_tmp(SIZE(var))) ELSE IF (SIZE(var) > SIZE(buff_tmp)) THEN IF (l_dbg) WRITE(ipslout,*) & "flinget_r1d : re-allocate buff_tmp for buff_sz = ",SIZE(var) DEALLOCATE (buff_tmp) ALLOCATE (buff_tmp(SIZE(var))) ENDIF !- CALL flinget_mat (fid_in,varname,iim,jjm,llm,ttm, & itau_dep,itau_fin,1,iim,1,jjm,buff_tmp) !- jl=0 DO ji=1,SIZE(var,1) jl=jl+1 var(ji) = buff_tmp(jl) ENDDO !------------------------- END SUBROUTINE flinget_r1d !- !=== !- SUBROUTINE flinget_r2d & (fid_in,varname,iim,jjm,llm,ttm,itau_dep,itau_fin,var) !--------------------------------------------------------------------- IMPLICIT NONE !- INTEGER :: fid_in CHARACTER(LEN=*) :: varname INTEGER :: iim, jjm, llm, ttm, itau_dep, itau_fin REAL :: var(:,:) !- INTEGER :: jl, jj, ji REAL,DIMENSION(:),ALLOCATABLE,SAVE :: buff_tmp LOGICAL :: l_dbg !--------------------------------------------------------------------- CALL ipsldbg (old_status=l_dbg) IF (.NOT.ALLOCATED(buff_tmp)) THEN IF (l_dbg) WRITE(ipslout,*) & "flinget_r2d : allocate buff_tmp for buff_sz = ",SIZE(var) ALLOCATE (buff_tmp(SIZE(var))) ELSE IF (SIZE(var) > SIZE(buff_tmp)) THEN IF (l_dbg) WRITE(ipslout,*) & "flinget_r2d : re-allocate buff_tmp for buff_sz = ",SIZE(var) DEALLOCATE (buff_tmp) ALLOCATE (buff_tmp(SIZE(var))) ENDIF !- CALL flinget_mat (fid_in,varname,iim,jjm,llm,ttm, & itau_dep,itau_fin,1,iim,1,jjm,buff_tmp) !- jl=0 DO jj=1,SIZE(var,2) DO ji=1,SIZE(var,1) jl=jl+1 var(ji,jj) = buff_tmp(jl) ENDDO ENDDO !------------------------- END SUBROUTINE flinget_r2d !- !=== !- SUBROUTINE flinget_r2d_zoom2d & (fid_in,varname,iim,jjm,llm,ttm, & itau_dep,itau_fin,iideb,iilen,jjdeb,jjlen,var) !--------------------------------------------------------------------- IMPLICIT NONE !- INTEGER :: fid_in CHARACTER(LEN=*) :: varname INTEGER :: iim,jjm,llm,ttm,itau_dep,itau_fin,iideb,jjdeb,iilen,jjlen REAL :: var(:,:) !- INTEGER :: jl, jj, ji REAL,DIMENSION(:),ALLOCATABLE,SAVE :: buff_tmp LOGICAL :: l_dbg !--------------------------------------------------------------------- CALL ipsldbg (old_status=l_dbg) IF (.NOT.ALLOCATED(buff_tmp)) THEN IF (l_dbg) WRITE(ipslout,*) & "flinget_r2d_zoom : allocate buff_tmp for buff_sz = ",SIZE(var) ALLOCATE (buff_tmp(SIZE(var))) ELSE IF (SIZE(var) > SIZE(buff_tmp)) THEN IF (l_dbg) WRITE(ipslout,*) & "flinget_r2d_zoom : re-allocate buff_tmp for buff_sz = ",SIZE(var) DEALLOCATE (buff_tmp) ALLOCATE (buff_tmp(SIZE(var))) ENDIF !- CALL flinget_mat (fid_in,varname,iim,jjm,llm,ttm, & itau_dep,itau_fin,iideb,iilen,jjdeb,jjlen,buff_tmp) !- jl=0 DO jj=1,SIZE(var,2) DO ji=1,SIZE(var,1) jl=jl+1 var(ji,jj) = buff_tmp(jl) ENDDO ENDDO !-------------------------------- END SUBROUTINE flinget_r2d_zoom2d !- !=== !- SUBROUTINE flinget_r3d & (fid_in,varname,iim,jjm,llm,ttm,itau_dep,itau_fin,var) !--------------------------------------------------------------------- IMPLICIT NONE !- INTEGER :: fid_in CHARACTER(LEN=*) :: varname INTEGER :: iim, jjm, llm, ttm, itau_dep, itau_fin REAL :: var(:,:,:) !- INTEGER :: jl, jk, jj, ji REAL,DIMENSION(:),ALLOCATABLE,SAVE :: buff_tmp LOGICAL :: l_dbg !--------------------------------------------------------------------- CALL ipsldbg (old_status=l_dbg) IF (.NOT.ALLOCATED(buff_tmp)) THEN IF (l_dbg) WRITE(ipslout,*) & "flinget_r3d : allocate buff_tmp for buff_sz = ",SIZE(var) ALLOCATE (buff_tmp(SIZE(var))) ELSE IF (SIZE(var) > SIZE(buff_tmp)) THEN IF (l_dbg) WRITE(ipslout,*) & "flinget_r3d : re-allocate buff_tmp for buff_sz = ",SIZE(var) DEALLOCATE (buff_tmp) ALLOCATE (buff_tmp(SIZE(var))) ENDIF !- CALL flinget_mat (fid_in,varname,iim,jjm,llm,ttm, & itau_dep,itau_fin,1,iim,1,jjm,buff_tmp) !- jl=0 DO jk=1,SIZE(var,3) DO jj=1,SIZE(var,2) DO ji=1,SIZE(var,1) jl=jl+1 var(ji,jj,jk) = buff_tmp(jl) ENDDO ENDDO ENDDO !------------------------- END SUBROUTINE flinget_r3d !- !=== !- SUBROUTINE flinget_r3d_zoom2d & (fid_in,varname,iim,jjm,llm,ttm, & itau_dep,itau_fin,iideb,iilen,jjdeb,jjlen,var) !--------------------------------------------------------------------- IMPLICIT NONE !- INTEGER :: fid_in CHARACTER(LEN=*) :: varname INTEGER :: iim,jjm,llm,ttm,itau_dep,itau_fin,iideb,jjdeb,iilen,jjlen REAL :: var(:,:,:) !- INTEGER :: jl, jk, jj, ji REAL,DIMENSION(:),ALLOCATABLE,SAVE :: buff_tmp LOGICAL :: l_dbg !--------------------------------------------------------------------- CALL ipsldbg (old_status=l_dbg) IF (.NOT.ALLOCATED(buff_tmp)) THEN IF (l_dbg) WRITE(ipslout,*) & "flinget_r3d_zoom : allocate buff_tmp for buff_sz = ",SIZE(var) ALLOCATE (buff_tmp(SIZE(var))) ELSE IF (SIZE(var) > SIZE(buff_tmp)) THEN IF (l_dbg) WRITE(ipslout,*) & "flinget_r3d_zoom : re-allocate buff_tmp for buff_sz = ",SIZE(var) DEALLOCATE (buff_tmp) ALLOCATE (buff_tmp(SIZE(var))) ENDIF !- CALL flinget_mat (fid_in,varname,iim,jjm,llm,ttm, & itau_dep,itau_fin,iideb,iilen,jjdeb,jjlen,buff_tmp) !- jl=0 DO jk=1,SIZE(var,3) DO jj=1,SIZE(var,2) DO ji=1,SIZE(var,1) jl=jl+1 var(ji,jj,jk) = buff_tmp(jl) ENDDO ENDDO ENDDO !-------------------------------- END SUBROUTINE flinget_r3d_zoom2d !- !=== !- SUBROUTINE flinget_r4d & (fid_in,varname,iim,jjm,llm,ttm,itau_dep,itau_fin,var) !--------------------------------------------------------------------- IMPLICIT NONE !- INTEGER :: fid_in CHARACTER(LEN=*) :: varname INTEGER :: iim, jjm, llm, ttm, itau_dep, itau_fin REAL :: var(:,:,:,:) !- INTEGER :: jl, jk, jj, ji, jm REAL,DIMENSION(:),ALLOCATABLE,SAVE :: buff_tmp LOGICAL :: l_dbg !--------------------------------------------------------------------- CALL ipsldbg (old_status=l_dbg) IF (.NOT.ALLOCATED(buff_tmp)) THEN IF (l_dbg) WRITE(ipslout,*) & "flinget_r4d : allocate buff_tmp for buff_sz = ",SIZE(var) ALLOCATE (buff_tmp(SIZE(var))) ELSE IF (SIZE(var) > SIZE(buff_tmp)) THEN IF (l_dbg) WRITE(ipslout,*) & "flinget_r4d : re-allocate buff_tmp for buff_sz = ",SIZE(var) DEALLOCATE (buff_tmp) ALLOCATE (buff_tmp(SIZE(var))) ENDIF !- CALL flinget_mat (fid_in,varname,iim,jjm,llm,ttm, & itau_dep,itau_fin,1,iim,1,jjm,buff_tmp) !- jl=0 DO jm=1,SIZE(var,4) DO jk=1,SIZE(var,3) DO jj=1,SIZE(var,2) DO ji=1,SIZE(var,1) jl=jl+1 var(ji,jj,jk,jm) = buff_tmp(jl) ENDDO ENDDO ENDDO ENDDO !------------------------- END SUBROUTINE flinget_r4d !- !=== !- SUBROUTINE flinget_r4d_zoom2d & (fid_in,varname,iim,jjm,llm,ttm, & itau_dep,itau_fin,iideb,iilen,jjdeb,jjlen,var) !--------------------------------------------------------------------- IMPLICIT NONE !- INTEGER :: fid_in CHARACTER(LEN=*) :: varname INTEGER :: iim,jjm,llm,ttm,itau_dep,itau_fin,iideb,jjdeb,iilen,jjlen REAL :: var(:,:,:,:) !- INTEGER :: jl, jk, jj, ji, jm REAL,DIMENSION(:),ALLOCATABLE,SAVE :: buff_tmp LOGICAL :: l_dbg !--------------------------------------------------------------------- CALL ipsldbg (old_status=l_dbg) IF (.NOT.ALLOCATED(buff_tmp)) THEN IF (l_dbg) WRITE(ipslout,*) & "flinget_r4d_zoom : allocate buff_tmp for buff_sz = ",SIZE(var) ALLOCATE (buff_tmp(SIZE(var))) ELSE IF (SIZE(var) > SIZE(buff_tmp)) THEN IF (l_dbg) WRITE(ipslout,*) & "flinget_r4d_zoom : re-allocate buff_tmp for buff_sz = ",SIZE(var) DEALLOCATE (buff_tmp) ALLOCATE (buff_tmp(SIZE(var))) ENDIF !- CALL flinget_mat (fid_in,varname,iim,jjm,llm,ttm, & itau_dep,itau_fin,iideb,iilen,jjdeb,jjlen,buff_tmp) !- jl=0 DO jm=1,SIZE(var,4) DO jk=1,SIZE(var,3) DO jj=1,SIZE(var,2) DO ji=1,SIZE(var,1) jl=jl+1 var(ji,jj,jk,jm) = buff_tmp(jl) ENDDO ENDDO ENDDO ENDDO !-------------------------------- END SUBROUTINE flinget_r4d_zoom2d !- !=== !- SUBROUTINE flinget_mat & (fid_in, varname, iim, jjm, llm, ttm, itau_dep, & itau_fin, iideb, iilen, jjdeb, jjlen, var) !--------------------------------------------------------------------- !- This subroutine will read the variable named varname from !- the file previously opened by flinopen and identified by fid !- !- It is checked that the dimensions of the variable to be read !- correspond to what the user requested when he specified !- iim, jjm and llm. The only exception which is allowed is !- for compressed data where the horizontal grid is not expected !- to be iim x jjm. !- !- If variable is of size zero a global attribute is read. !- This global attribute will be of type real !- !- INPUT !- !- fid : File ID returned by flinopen !- varname : Name of the variable to be read from the file !- iim : | These three variables give the size of the variables !- jjm : | to be read. It will be verified that the variables !- llm : | fits in there. !- ttm : | !- itau_dep : Time step at which we will start to read !- itau_fin : Time step until which we are going to read !- For the moment this is done on indexes !- but it should be in the physical space. !- If there is no time-axis in the file then use a !- itau_fin < itau_dep, this will tell flinget not to !- expect a time-axis in the file. !- iideb : index i for zoom !- iilen : length of zoom !- jjdeb : index j for zoom !- jjlen : length of zoom !- !- OUTPUT !- !- var : array that will contain the data !--------------------------------------------------------------------- IMPLICIT NONE !- ! ARGUMENTS !- INTEGER, INTENT(IN) :: fid_in CHARACTER(LEN=*), INTENT(IN) :: varname INTEGER, INTENT(IN) :: iim, jjm, llm, ttm INTEGER, INTENT(IN) :: itau_dep, itau_fin, iideb, iilen, jjdeb, jjlen REAL, INTENT(OUT) :: var(:) !- ! LOCAL !- INTEGER :: iret, fid INTEGER :: vid, cvid, clen CHARACTER(LEN=70) :: str1 CHARACTER(LEN=250) :: att_n, tmp_n CHARACTER(LEN=5) :: axs_l INTEGER :: tmp_i REAL,SAVE :: mis_v=0. REAL :: tmp_r INTEGER :: ndims, x_typ, nb_atts INTEGER,DIMENSION(NF90_MAX_VAR_DIMS) :: dimids INTEGER :: i, nvars, i2d, cnd REAL,DIMENSION(:),ALLOCATABLE,SAVE :: var_tmp INTEGER :: itau_len LOGICAL :: uncompress = .FALSE. INTEGER :: il, ip, i2p, it !- LOGICAL :: l_dbg !--------------------------------------------------------------------- CALL ipsldbg (old_status=l_dbg) !- fid = ncids(fid_in) !- IF (l_dbg) THEN WRITE(ipslout,*) & 'flinget_mat : fid_in, fid, varname :', fid_in, fid, TRIM(varname) WRITE(ipslout,*) & 'flinget_mat : iim, jjm, llm, ttm, itau_dep, itau_fin :', & iim, jjm, llm, ttm, itau_dep, itau_fin WRITE(ipslout,*) & 'flinget_mat : iideb, iilen, jjdeb, jjlen :', & iideb, iilen, jjdeb, jjlen ENDIF !- uncompress = .FALSE. !- ! 1.0 We get first all the details on this variable from the file !- nvars = ncnbva(fid_in) !- vid = -1 iret = NF90_INQ_VARID (fid, varname, vid) !- IF (vid < 0 .OR. iret /= NF90_NOERR) THEN CALL histerr (3,'flinget', & 'Variable '//TRIM(varname)//' not found in file',' ',' ') ENDIF !- iret = NF90_INQUIRE_VARIABLE (fid, vid, & ndims=ndims, dimids=dimids, nAtts=nb_atts) IF (l_dbg) THEN WRITE(ipslout,*) & 'flinget_mat : fid, vid :', fid, vid WRITE(ipslout,*) & 'flinget_mat : ndims, dimids(1:ndims), nb_atts :', & ndims, dimids(1:ndims), nb_atts ENDIF !- w_dim(:) = 0 DO i=1,ndims iret = NF90_INQUIRE_DIMENSION (fid, dimids(i), len=w_dim(i)) ENDDO IF (l_dbg) WRITE(ipslout,*) & 'flinget_mat : w_dim :', w_dim(1:ndims) !- mis_v = 0.0; axs_l = ' '; !- IF (nb_atts > 0) THEN IF (l_dbg) THEN WRITE(ipslout,*) 'flinget_mat : attributes for variable :' ENDIF ENDIF DO i=1,nb_atts iret = NF90_INQ_ATTNAME (fid, vid, i, att_n) iret = NF90_INQUIRE_ATTRIBUTE (fid, vid, att_n, xtype=x_typ) CALL strlowercase (att_n) IF ( (x_typ == NF90_INT).OR.(x_typ == NF90_SHORT) & .OR.(x_typ == NF90_BYTE) ) THEN iret = NF90_GET_ATT (fid, vid, att_n, tmp_i) IF (l_dbg) THEN WRITE(ipslout,*) ' ',TRIM(att_n),' : ',tmp_i ENDIF ELSE IF ( (x_typ == NF90_FLOAT).OR.(x_typ == NF90_DOUBLE) ) THEN iret = NF90_GET_ATT (fid, vid, att_n, tmp_r) IF (l_dbg) THEN WRITE(ipslout,*) ' ',TRIM(att_n),' : ',tmp_r ENDIF IF (index(att_n,'missing_value') > 0) THEN mis_v = tmp_r ENDIF ELSE tmp_n = '' iret = NF90_GET_ATT (fid, vid, att_n, tmp_n) IF (l_dbg) THEN WRITE(ipslout,*) ' ',TRIM(att_n),' : ',TRIM(tmp_n) ENDIF IF (index(att_n,'axis') > 0) THEN axs_l = tmp_n ENDIF ENDIF ENDDO !? !!!!!!!!!! We will need a verification on the type of the variable !? !- ! 2.0 The dimensions are analysed to determine what is to be read !- ! 2.1 the longitudes !- IF ( w_dim(1) /= iim .OR. w_dim(2) /= jjm) THEN !--- !-- There is a possibility that we have to deal with a compressed axis ! !--- iret = NF90_INQUIRE_DIMENSION (fid, dimids(1), & name=tmp_n, len=clen) iret = NF90_INQ_VARID (fid, tmp_n, cvid) !--- IF (l_dbg) WRITE(ipslout,*) & 'Dimname, iret , NF90_NOERR : ',TRIM(tmp_n),iret,NF90_NOERR !--- !-- If we have an axis which has the same name !-- as the dimension we can see if it is compressed !--- !-- TODO TODO for zoom2d !--- IF (iret == NF90_NOERR) THEN iret = NF90_GET_ATT (fid, cvid, 'compress', str1) !----- IF (iret == NF90_NOERR) THEN iret = NF90_INQUIRE_VARIABLE (fid,cvid,xtype=x_typ,ndims=cnd) !------- IF ( cnd /= 1 .AND. x_typ /= NF90_INT) THEN CALL histerr (3,'flinget', & 'Variable '//TRIM(tmp_n)//' can not be a compressed axis', & 'Either it has too many dimensions'// & ' or it is not of type integer', ' ') ELSE !--------- !-------- Let us see if we already have that index table !--------- IF ( (cind_len /= clen).OR.(cind_vid /= cvid) & .OR.(cind_fid /= fid) ) THEN IF (ALLOCATED(cindex)) DEALLOCATE(cindex) ALLOCATE(cindex(clen)) cind_len = clen cind_vid = cvid cind_fid = fid iret = NF90_GET_VAR (fid, cvid, cindex) ENDIF !--------- !-------- In any case we need to set the slab of data to be read !--------- uncompress = .TRUE. w_sta(1) = 1 w_len(1) = clen i2d = 1 ENDIF ELSE str1 = 'The horizontal dimensions of '//varname CALL histerr (3,'flinget',str1, & 'is not compressed and does not'// & ' correspond to the requested size',' ') ENDIF ELSE IF (w_dim(1) /= iim) THEN str1 = 'The longitude dimension of '//varname CALL histerr (3,'flinget',str1, & 'in the file is not equal to the dimension', & 'that should be read') ENDIF IF (w_dim(2) /= jjm) THEN str1 = 'The latitude dimension of '//varname CALL histerr (3,'flinget',str1, & 'in the file is not equal to the dimension', & 'that should be read') ENDIF ENDIF ELSE w_sta(1:2) = (/ iideb, jjdeb /) w_len(1:2) = (/ iilen, jjlen /) i2d = 2 ENDIF !- ! 2.3 Now the difficult part, the 3rd dimension which can be ! time or levels. !- ! Priority is given to the time axis if only three axes are present. !- IF (ndims > i2d) THEN !--- !-- 2.3.1 We have a vertical axis !--- IF (llm == 1 .AND. ndims == i2d+2 .OR. llm == w_dim(i2d+1)) THEN !----- IF (w_dim(i2d+1) /= llm) THEN CALL histerr (3,'flinget', & 'The vertical dimension of '//varname, & 'in the file is not equal to the dimension', & 'that should be read') ELSE w_sta(i2d+1) = 1 IF (llm > 0) THEN w_len(i2d+1) = llm ELSE w_len(i2d+1) = w_sta(i2d+1) ENDIF ENDIF !----- IF ((itau_fin-itau_dep) >= 0) THEN IF (ndims /= i2d+2) THEN CALL histerr (3,'flinget', & 'You attempt to read a time slab', & 'but there is no time axis on this variable', varname) ELSE IF ((itau_fin - itau_dep) <= w_dim(i2d+2)) THEN w_sta(i2d+2) = itau_dep w_len(i2d+2) = itau_fin-itau_dep+1 ELSE CALL histerr (3,'flinget', & 'The time step you try to read is not', & 'in the file (1)', varname) ENDIF ELSE IF (ndims == i2d+2 .AND. w_dim(i2d+2) > 1) THEN CALL histerr (3,'flinget', & 'There is a time axis in the file but no', & 'time step give in the call', varname) ELSE w_sta(i2d+2) = 1 w_len(i2d+2) = 1 ENDIF ELSE !----- !---- 2.3.2 We do not have any vertical axis !----- IF (ndims == i2d+2) THEN CALL histerr (3,'flinget', & 'The file contains 4 dimensions', & 'but only 3 are requestes for variable ', varname) ENDIF IF ((itau_fin-itau_dep) >= 0) THEN IF (ndims == i2d+1) THEN IF ((itau_fin-itau_dep) < w_dim(i2d+1) ) THEN w_sta(i2d+1) = itau_dep w_len(i2d+1) = itau_fin-itau_dep+1 ELSE CALL histerr (3,'flinget', & 'The time step you try to read is not', & 'in the file (2)', varname) ENDIF ELSE CALL histerr (3,'flinget', & 'From your input you sould have 3 dimensions', & 'in the file but there are 4', varname) ENDIF ELSE IF (ndims == i2d+1 .AND. w_dim(i2d+1) > 1) THEN CALL histerr (3,'flinget', & 'There is a time axis in the file but no', & 'time step given in the call', varname) ELSE w_sta(i2d+1) = 1 w_len(i2d+1) = 1 ENDIF ENDIF ENDIF ELSE !--- !-- 2.3.3 We do not have any vertical axis !--- w_sta(i2d+1:i2d+2) = (/ 0, 0 /) w_len(i2d+1:i2d+2) = (/ 0, 0 /) ENDIF !- ! 3.0 Reading the data !- IF (l_dbg) WRITE(ipslout,*) & 'flinget_mat 3.0 : ', uncompress, w_sta, w_len !--- var(:) = mis_v IF (uncompress) THEN !--- IF (ALLOCATED(var_tmp)) THEN IF (SIZE(var_tmp) < PRODUCT(w_len(:),mask=(w_len>1))) THEN DEALLOCATE(var_tmp) ALLOCATE(var_tmp(PRODUCT(w_len(:),mask=(w_len>1)))) ENDIF ELSE ALLOCATE(var_tmp(PRODUCT(w_len(:),mask=(w_len>1)))) ENDIF !--- iret = NF90_GET_VAR (fid, vid, var_tmp, & start=w_sta(:), count=w_len(:)) !--- itau_len=itau_fin-itau_dep+1 IF (l_dbg) WRITE(ipslout,*) 'flinget_mat 3.0 : clen, itau_len ',clen,itau_len var(:) = mis_v IF (itau_len > 0) THEN DO it=1,itau_len DO il=1,clen ip = il + (it-1)*clen i2p = cindex(il)+(it-1)*iim*jjm var(i2p) = var_tmp(ip) ENDDO ENDDO ELSE var(cindex(:)) = var_tmp(:) ENDIF !--- ELSE iret = NF90_GET_VAR (fid, vid, var, & start=w_sta(:), count=w_len(:)) ENDIF !- IF (l_dbg) WRITE(ipslout,*) 'flinget_mat 3.1 : ',NF90_STRERROR (iret) !-------------------------- END SUBROUTINE flinget_mat !- !=== !- SUBROUTINE flinget_scal & (fid_in, varname, iim, jjm, llm, ttm, itau_dep, itau_fin, var) !--------------------------------------------------------------------- !- This subroutine will read the variable named varname from !- the file previously opened by flinopen and identified by fid !- !- If variable is of size zero a global attribute is read. This !- global attribute will be of type real !- !- INPUT !- !- fid : File ID returned by flinopen !- varname : Name of the variable to be read from the file !- iim : | These three variables give the size of the variables !- jjm : | to be read. It will be verified that the variables !- llm : | fits in there. !- ttm : | !- itau_dep : Time step at which we will start to read !- itau_fin : Time step until which we are going to read !- For the moment this is done on indeces but it should be !- in the physical space !- If there is no time-axis in the file then use a !- itau_fin < itau_dep, this will tell flinget not to !- expect a time-axis in the file. !- !- OUTPUT !- !- var : scalar that will contain the data !--------------------------------------------------------------------- IMPLICIT NONE !- ! ARGUMENTS !- INTEGER :: fid_in CHARACTER(LEN=*) :: varname INTEGER :: iim, jjm, llm, ttm, itau_dep, itau_fin REAL :: var !- ! LOCAL !- INTEGER :: iret, fid, vid INTEGER :: attlen, attnum INTEGER :: ndims, nb_atts INTEGER,DIMENSION(NF90_MAX_VAR_DIMS) :: dimids LOGICAL :: var_exists !- LOGICAL :: l_dbg INTEGER :: lll !--------------------------------------------------------------------- CALL ipsldbg (old_status=l_dbg) IF (l_dbg) THEN WRITE (ipslout,*) 'flinget_scal in file with id ',fid_in ENDIF !- fid = ncids(fid_in) iret = NF90_INQUIRE_ATTRIBUTE(fid, NF90_GLOBAL, varname, len=attlen, attnum=attnum) !- ! 1.0 Reading a global attribute !- IF ( iret == nf90_noerr ) THEN ! ! This seems to be a Global attribute ! iret = NF90_GET_ATT (fid, NF90_GLOBAL, varname, var) ELSE ! ! If there was an error on the test for a global attribute it ! is perhaps a scalar variable. ! vid = -1 iret = NF90_INQ_VARID (fid, varname, vid) ! IF ( (vid >= 0).AND.(iret == NF90_NOERR) ) THEN iret = NF90_INQUIRE_VARIABLE (fid, vid, & ndims=ndims, dimids=dimids, nAtts=nb_atts) IF (ndims == 1) THEN iret = NF90_INQUIRE_DIMENSION (fid, dimids(1), len=lll) ENDIF IF ( ((ndims == 0) .OR. ((ndims == 1).AND.(lll == 1))) .AND. (nb_atts >= 0) ) THEN iret = NF90_GET_VAR(fid, vid, var) ELSE CALL histerr (3,'flinget_scal', & 'The variable has coordinates and thus is probably not a scalar.', & 'Check your netCDF file.', " ") ENDIF ENDIF IF (l_dbg) THEN WRITE(ipslout,*) "Reading a Scalar value for varibale ", varname," It has value ", var ENDIF ENDIF !- !--------------------------- END SUBROUTINE flinget_scal !- !=== !- SUBROUTINE flinfindcood (fid_in, axtype, vid, ndim) !--------------------------------------------------------------------- !- This subroutine explores the file in order to find !- the coordinate according to a number of rules !--------------------------------------------------------------------- IMPLICIT NONE !- ! ARGUMENTS !- INTEGER :: fid_in, vid, ndim CHARACTER(LEN=3) :: axtype !- ! LOCAL !- INTEGER :: iv, iret, dimnb CHARACTER(LEN=40) :: dimname, dimuni1, dimuni2, dimuni3 CHARACTER(LEN=80) :: str1 LOGICAL :: found_rule = .FALSE. !--------------------------------------------------------------------- vid = -1 !- ! Make sure all strings are invalid !- dimname = '?-?' dimuni1 = '?-?' dimuni2 = '?-?' dimuni3 = '?-?' !- ! First rule : we look for the correct units ! lon : east ! lat : north ! We make an exact check as it would be too easy to mistake ! some units by just comparing the substrings. !- SELECTCASE(axtype) CASE ('lon') dimuni1 = 'degree_e' dimuni2 = 'degrees_e' found_rule = .TRUE. CASE('lat') dimuni1 = 'degree_n' dimuni2 = 'degrees_n' found_rule = .TRUE. CASE('lev') dimuni1 = 'm' dimuni2 = 'km' dimuni3 = 'hpa' found_rule = .TRUE. CASE DEFAULT found_rule = .FALSE. END SELECT !- IF (found_rule) THEN iv = 0 DO WHILE ( (vid < 0).AND.(iv < ncnbva(fid_in)) ) iv = iv+1 str1 = '' iret = NF90_GET_ATT (ncids(fid_in), iv, 'units', str1) IF (iret == NF90_NOERR) THEN CALL strlowercase (str1) IF ( (INDEX(str1, TRIM(dimuni1)) == 1) & .OR.(INDEX(str1, TRIM(dimuni2)) == 1) & .OR.(INDEX(str1, TRIM(dimuni3)) == 1) ) THEN vid = iv iret = NF90_INQUIRE_VARIABLE (ncids(fid_in), iv, ndims=ndim) ENDIF ENDIF ENDDO ENDIF !- ! Second rule : we find specific names : ! lon : nav_lon ! lat : nav_lat ! Here we can check if we find the substring as the ! names are more specific. !- SELECTCASE(axtype) CASE ('lon') dimname = 'nav_lon lon longitude' found_rule = .TRUE. CASE('lat') dimname = 'nav_lat lat latitude' found_rule = .TRUE. CASE('lev') dimname = 'plev level depth deptht' found_rule = .TRUE. CASE DEFAULT found_rule = .FALSE. END SELECT !- IF (found_rule) THEN iv = 0 DO WHILE ( (vid < 0).AND.(iv < ncnbva(fid_in)) ) iv = iv+1 str1='' iret = NF90_INQUIRE_VARIABLE (ncids(fid_in), iv, & name=str1, ndims=ndim) IF (INDEX(dimname,TRIM(str1)) >= 1) THEN vid = iv ENDIF ENDDO ENDIF !- ! Third rule : we find a variable with the same name as the dimension ! lon = 1 ! lat = 2 ! lev = 3 !- IF (vid < 0) THEN SELECTCASE(axtype) CASE ('lon') dimnb = 1 found_rule = .TRUE. CASE('lat') dimnb = 2 found_rule = .TRUE. CASE('lev') dimnb = 3 found_rule = .TRUE. CASE DEFAULT found_rule = .FALSE. END SELECT !--- IF (found_rule) THEN iret = NF90_INQUIRE_DIMENSION (ncids(fid_in), dimnb, name=dimname) iv = 0 DO WHILE ( (vid < 0).AND.(iv < ncnbva(fid_in)) ) iv = iv+1 str1='' iret = NF90_INQUIRE_VARIABLE (ncids(fid_in), iv, & name=str1, ndims=ndim) IF (INDEX(dimname,TRIM(str1)) == 1) THEN vid = iv ENDIF ENDDO ENDIF ENDIF !- ! Stop the program if no coordinate was found !- IF (vid < 0) THEN CALL histerr (3,'flinfindcood', & 'No coordinate axis was found in the file', & 'The data in this file can not be used', axtype) ENDIF !-------------------------- END SUBROUTINE flinfindcood !- !=== !- SUBROUTINE flinclo (fid_in) !--------------------------------------------------------------------- IMPLICIT NONE !- INTEGER :: fid_in !- INTEGER :: iret !--------------------------------------------------------------------- iret = NF90_CLOSE (ncids(fid_in)) ncfileopen(fid_in) = .FALSE. !--------------------- END SUBROUTINE flinclo !- !=== !- SUBROUTINE flinquery_var(fid_in, varname, exists) !--------------------------------------------------------------------- !- Queries the existance of a variable in the file. !--------------------------------------------------------------------- IMPLICIT NONE !- INTEGER :: fid_in CHARACTER(LEN=*) varname LOGICAL :: exists !- INTEGER :: iret, fid, vid !--------------------------------------------------------------------- fid = ncids(fid_in) vid = -1 iret = NF90_INQ_VARID (fid, varname, vid) !- exists = ( (vid >= 0).AND.(iret == NF90_NOERR) ) !--------------------------- END SUBROUTINE flinquery_var !- !=== !- SUBROUTINE flininspect (fid_in) !--------------------------------------------------------------------- IMPLICIT NONE !- ! fid : File id to inspect !- INTEGER :: fid_in !- !- LOCAL !- INTEGER :: iim, jjm, llm, ttm, fid_out INTEGER :: iret, fid, ndims, nvars, nb_atts, id_unlim INTEGER :: iv, in, lll INTEGER :: xid, yid, zid, tid INTEGER,DIMENSION(NF90_MAX_VAR_DIMS) :: idimid CHARACTER(LEN=80) :: name CHARACTER(LEN=30) :: axname !--------------------------------------------------------------------- fid = ncids(fid_in) !- iret = NF90_INQUIRE (fid, nDimensions=ndims, nVariables=nvars, & nAttributes=nb_atts, unlimitedDimId=id_unlim) !- WRITE (*,*) 'IOIPSL ID : ',fid_in WRITE (*,*) 'NetCDF ID : ',fid WRITE (*,*) 'Number of dimensions : ',ndims WRITE (*,*) 'Number of variables : ',nvars WRITE (*,*) 'Number of global attributes : ',nb_atts WRITE (*,*) 'ID unlimited : ',id_unlim !- xid = -1; iim = 0; yid = -1; jjm = 0; zid = -1; llm = 0; tid = -1; ttm = 0; !- DO iv=1,ndims !--- iret = NF90_INQUIRE_DIMENSION (fid, iv, name=axname, len=lll) CALL strlowercase (axname) axname = ADJUSTL(axname) !--- WRITE (*,*) 'Dimension number : ',iv WRITE (*,*) 'Dimension name : ',TRIM(axname) !--- IF ( (INDEX(axname,'x') == 1) & .OR.(INDEX(axname,'lon') == 1)) THEN xid = iv; iim = lll; WRITE (*,*) 'Dimension X size : ',iim ELSE IF ( (INDEX(axname,'y') == 1) & .OR.(INDEX(axname,'lat') == 1)) THEN yid = iv; jjm = lll; WRITE (*,*) 'Dimension Y size : ',jjm ELSE IF ( (INDEX(axname,'lev') == 1) & .OR.(INDEX(axname,'plev') == 1) & .OR.(INDEX(axname,'z') == 1) & .OR.(INDEX(axname,'depth') == 1)) THEN zid = iv; llm = lll; WRITE (*,*) 'Dimension Z size : ',llm ELSE IF ( (INDEX(axname,'tstep') == 1) & .OR.(INDEX(axname,'time') == 1) & .OR.(INDEX(axname,'time_counter') == 1)) THEN !---- For the time we certainly need to allow for other names tid = iv; ttm = lll; ELSE IF (ndims == 1) THEN !---- Nothing was found and ndims=1 then we have a vector of data xid = 1; iim = lll; ENDIF !--- ENDDO !- ! Keep all this information !- nbfiles = nbfiles+1 !- IF (nbfiles > nbfile_max) THEN CALL histerr(3,'flininspect', & 'Too many files. Please increase nbfil_max', & 'in program flincom.F90.',' ') ENDIF !- ncids(nbfiles) = fid ncnbd(nbfiles) = ndims !- ncdims(nbfiles,1:4) = (/ iim, jjm, llm, ttm /) !- ncfunli(nbfiles) = id_unlim ncnba(nbfiles) = nb_atts ncnbva(nbfiles) = nvars ncfileopen(nbfiles) = .TRUE. !- fid_out = nbfiles !- DO in=1,nvars iret = NF90_INQUIRE_VARIABLE (fid, in, & name=name, ndims=ndims, dimids=idimid, nAtts=nb_atts) WRITE (*,*) 'Variable number ------------ > ', in WRITE (*,*) 'Variable name : ', TRIM(name) WRITE (*,*) 'Number of dimensions : ', ndims WRITE (*,*) 'Dimensions ID''s : ', idimid(1:ndims) WRITE (*,*) 'Number of attributes : ', nb_atts ENDDO !------------------------- END SUBROUTINE flininspect !- !=== !- END MODULE flincom