program caliporbit use netcdf implicit none character(len=200) :: filein,month2, file2, filemod, modfileout, modfileout2,daytime integer :: err, imonth, id, ifile,day, month character*2 :: dayc,monthc real :: t,hour integer :: ilon, ilat, ialt, nprof, iz, i, icat, itemp,n, lontmp, lattmp real,dimension(:),allocatable :: lat,lon,time ! instant lat/lon integer,parameter :: lonmax=95,latmax=95,altmax=40, daymax2=120, daymax=30 real,dimension(latmax+1) :: latmod real,dimension(lonmax+1) :: lonmod real,dimension(latmax) :: latmid real,dimension(lonmax) :: lonmid real,dimension(altmax) :: altmid real,dimension(altmax+1) :: alt real,dimension(lonmax+1,latmax+1,altmax,daymax2) :: clcalipso real,dimension(lonmax+1,latmax+1,daymax2,4) :: cltcalipso real,dimension(lonmax+1,latmax+1,altmax,daymax) :: cf3d,cf3dind real,dimension(lonmax+1,latmax+1,daymax,4) :: cf2d,cf2dind 100 format(A200) !filemod='/bdd/CFMIP/MODEL_LOCAL/LMDZ/LMDZg9520121016.testing/modipsl/modeles/LMDZ5/BENCH96x95x39/200701/histhfCOSP.nc' cf3d(:,:,:,:)=0. cf3dind(:,:,:,:)=0. cf2d(:,:,:,:)=0. cf2dind(:,:,:,:)=0. cltcalipso(:,:,:,:)=0. clcalipso(:,:,:,:)=0. ! Reading the daily fime list do write (*,'(a)',advance='no') 'input file name = ' read *, filein open(unit=1,file=filein,iostat=err,status='OLD') if (err==0) exit print *,'--- input file not found' enddo do write (*,'(a)',advance='no') 'input month = ' read *, month2 if (err==0) exit enddo do write (*,'(a)',advance='no') 'input day/night = ' read *, daytime if (err==0) exit enddo filemod='/bdd/CFMIP/MODEL_LOCAL/LMDZ/LMDZg9520121016.testing/modipsl/modeles/LMDZ5/BENCH96x95x39/2007'//trim(month2)//'/histhfCOSP.nc' ! definition of lonmod and latmod latmod(1)=90 do ilat=1,latmax latmod(ilat+1)=latmod(ilat)-1.875 enddo latmid=latmod(1:latmax)-0.9375 lonmod(1)=-180 do ilon=1,lonmax lonmod(ilon+1)=lonmod(ilon)+3.75 enddo lonmid=lonmod(1:lonmax)+1.875 alt(1)=0 do ilat=1,altmax alt(ilat+1)=alt(ilat)+0.48 enddo altmid=alt(1:altmax)+0.24 !print *, 'reading file2d' call rdnc3(trim(filemod),lonmax+1,latmax+1,daymax2,cltcalipso(:,:,:,1),'cltcalipso') call rdnc3(trim(filemod),lonmax+1,latmax+1,daymax2,cltcalipso(:,:,:,2),'clhcalipso') call rdnc3(trim(filemod),lonmax+1,latmax+1,daymax2,cltcalipso(:,:,:,3),'clmcalipso') call rdnc3(trim(filemod),lonmax+1,latmax+1,daymax2,cltcalipso(:,:,:,4),'cllcalipso') !print *, 'reading file3d' call rdnc4(trim(filemod),lonmax+1,latmax+1,altmax,daymax2,clcalipso,'clcalipso') ! Reading the files 888 read(1,100,end=999)file2 !print *, 'file :',trim(file2) dayc=file2(65:66) monthc=file2(62:63) call getdim(trim(file2),'it',nprof) allocate(lon(nprof),lat(nprof),time(nprof)) !print *, 'allocation ok' !print *, 'prof:',nprof call rdnc(trim(file2),nprof,lon,'longitude') call rdnc(trim(file2),nprof,lat,'latitude') call rdnc(trim(file2),nprof,time,'time') read(dayc,'(I2)')day read(monthc,'(I2)')month print *, 'month, day ',month,' ',day hour=time(1) hour=nint(hour) t=(day-1)*4+(hour/6+1) t=nint(t) !print *, 'day, hour, ind',day,hour,t lontmp=1 lattmp=1 if(t.le.120)then do i=1,nprof ! Cloud Fraction 3D do ilon=1,lonmax !longitude if( (lon(i).ge.lonmod(ilon)) .and. (lon(i).lt.lonmod(ilon+1)) )then do ilat=1,latmax !latitude if ( (lat(i).le.latmod(ilat)) .and. (lat(i).gt.latmod(ilat+1)) )then if( (i.gt.1).and. ((lontmp.ne.ilon).or.(lattmp.ne.ilat)) )then lontmp=ilon lattmp=ilat !print *, ilon,ilat do iz=1,40 if (clcalipso(lontmp,lattmp,iz,t).lt.105) then cf3d(lontmp,lattmp,iz,day)=cf3d(lontmp,lattmp,iz,day)+clcalipso(lontmp,lattmp,iz,t) cf3dind(lontmp,lattmp,iz,day)=cf3dind(lontmp,lattmp,iz,day)+1 endif enddo do icat=1,4 if (cltcalipso(lontmp,lattmp,t,icat).lt.105) then cf2d(lontmp,lattmp,day,icat)=cf2d(lontmp,lattmp,day,icat)+cltcalipso(lontmp,lattmp,t,icat) cf2dind(lontmp,lattmp,day,icat)=cf2dind(lontmp,lattmp,day,icat)+1 !print *, cltcalipso(lontmp,lattmp,t,icat) endif enddo endif endif enddo endif enddo enddo endif deallocate(lat,lon,time) !print *, 'Next file' goto 888 999 continue do day=1,30 do ilon=1,lonmax+1 do ilat=1,latmax+1 do icat=1,4 if(cf2dind(ilon,ilat,day,icat)>0.)then cf2d(ilon,ilat,day,icat)=cf2d(ilon,ilat,day,icat)/cf2dind(ilon,ilat,day,icat) else cf2d(ilon,ilat,day,icat)=-9999. endif enddo do iz = 1, altmax ! loop by vertical bins if(cf3dind(ilon,ilat,iz,day)>0.)then cf3d(ilon,ilat,iz,day)=cf3d(ilon,ilat,iz,day)/cf3dind(ilon,ilat,iz,day) else cf3d(ilon,ilat,iz,day)=-9999. endif enddo enddo enddo enddo ! Recording daily files print *, 'Recording files' modfileout='/bdd/CFMIP/workspace/gcesana/CF3D330m_2007'//monthc//'_'//trim(daytime)//'_LMDZ.nc' modfileout2='/bdd/CFMIP/workspace/gcesana/CF2D330m_2007'//monthc//'_'//trim(daytime)//'_LMDZ.nc' call create_3d(trim(modfileout),'clcalipso','lon','lat','alt', & lonmax+1,latmax+1,altmax,cf3d,lonmod,latmod,altmid) call create_2d(trim(modfileout2),'clcalipso','lon','lat', & lonmax+1,latmax+1,cf2d,lonmod,latmod) contains ! Check the status of the nf90 command subroutine check(status) implicit none include 'netcdf.inc' integer, intent ( in) :: status if(status /= nf90_noerr) then print *, trim(nf90_strerror(status)) stop "Stopped" end if end subroutine check subroutine getdim(fname,vname,fdim) use netcdf implicit none integer :: ncid, dimid, fdim character(len=*) :: fname,vname !print *, trim(fname) call check(NF90_OPEN(trim(fname),NF90_NOWRITE,ncid)) call check(NF90_inq_dimid(ncid, trim(vname), dimid)) call check(NF90_inquire_dimension(ncid, dimid, len = fdim)) end subroutine getdim ! Read 1D netcdf variable subroutine rdnc(fname,fdim,var,nvar) use netcdf implicit none character(len=*) :: nvar, fname integer :: fdim integer :: ncid,varid real,dimension(fdim) :: var call check(NF90_OPEN(trim(fname),NF90_NOWRITE,ncid)) call check(NF90_inq_varid(ncid,nvar,varid)) call check(NF90_get_var(ncid,varid,var)) call check(NF90_CLOSE(ncid)) end subroutine rdnc subroutine rdnc2(fname,prof,alt,var,nvar) use netcdf implicit none character(len=*) :: nvar, fname integer :: alt,prof integer :: ncid,varid real,dimension(alt,prof) :: var !print *, trim(fname) call check(NF90_OPEN(trim(fname),NF90_NOWRITE,ncid)) call check(NF90_inq_varid(ncid,nvar,varid)) call check(NF90_get_var(ncid,varid,var)) call check(NF90_CLOSE(ncid)) end subroutine rdnc2 subroutine rdnc3(fname,lon,lat,alt,var,nvar) use netcdf implicit none character(len=*) :: nvar, fname integer :: lon,lat,alt integer :: ncid,varid real,dimension(lon,lat,alt) :: var print *, trim(fname) call check(NF90_OPEN(trim(fname),NF90_NOWRITE,ncid)) call check(NF90_inq_varid(ncid,nvar,varid)) call check(NF90_get_var(ncid,varid,var)) call check(NF90_CLOSE(ncid)) end subroutine rdnc3 subroutine rdnc4(fname,lon,lat,alt,day,var,nvar) use netcdf implicit none character(len=*) :: nvar, fname integer :: lon,lat,alt,day integer :: ncid,varid real,dimension(lon,lat,alt,day) :: var call check(NF90_OPEN(trim(fname),NF90_NOWRITE,ncid)) call check(NF90_inq_varid(ncid,nvar,varid)) call check(NF90_get_var(ncid,varid,var)) call check(NF90_CLOSE(ncid)) end subroutine rdnc4 subroutine create_3d(fname,vname,ndim1,ndim2,ndim3,dim1,dim2,dim3,var5,lonmod,latmod,tempmod) use netcdf implicit none character(LEN=*) :: fname,ndim1,ndim2,ndim3,vname character*3 :: ndim4='day' integer,parameter :: ndims=4 real,parameter :: nan=-9999 integer,parameter :: dim4=30 integer :: i,dim(ndims),nc_id,dim1,dim2,dim3 integer :: varid1, varid2,tiime_varid,varid3,varid4,varid5 integer :: dimid1,dimid2,dimid3,dimid4, time_dimid real*4,dimension(dim1,dim2,dim3,dim4) :: var5 real,dimension(dim1) :: lonmod real,dimension(dim2) :: latmod real,dimension(dim3) :: altmod real,dimension(dim4) :: tempmod tempmod(:)=0. do i=1,30 tempmod(i)=i enddo ! CREATION DU FICHIER nf90_create call check(nf90_create(fname, NF90_CLOBBER, & nc_id)) call check(nf90_put_att(nc_id, NF90_GLOBAL, 'Description', & 'LMDZ CALIOP Orbit files')) call check(nf90_put_att(nc_id, NF90_GLOBAL, 'Version','caliporbit2')) call check(nf90_put_att(nc_id, NF90_GLOBAL, 'Author', & 'Gregory CESANA; Cesana et al., 2015, Evaluation of 3 CALIPSO & Cloud Phase Products Using In-Situ Airborne Measurements')) ! DEFINITION DES DIMENSIONS call check(nf90_def_dim(nc_id,ndim1, dim1, dimid1)) call check(nf90_def_dim(nc_id,ndim2, dim2, dimid2 )) call check(nf90_def_dim(nc_id,ndim3, dim3, dimid3)) call check(nf90_def_dim(nc_id,ndim4, NF90_UNLIMITED, dimid4)) ! call check(nf90_def_dim(nc_id, 'time', NF90_UNLIMITED, time_dimid)) ! DIMENSION DES MATRICES dim = (/dimid1, dimid2, dimid3, dimid4/) ! VARIABLE 4D +time call check(nf90_def_var(nc_id, trim(vname), NF90_FLOAT, dim, varid5)) call check(nf90_put_att(nc_id, varid5, '_FillValue',nan)) call check(nf90_def_var(nc_id, ndim1, NF90_FLOAT, dimid1, varid1)) call check(nf90_def_var(nc_id, ndim2, NF90_FLOAT, dimid2, varid2)) call check(nf90_def_var(nc_id, ndim3, NF90_FLOAT, dimid3, varid3)) call check(nf90_def_var(nc_id, ndim4, NF90_FLOAT, dimid4, varid4)) ! FIN DU MODE DEFINITION call check(nf90_enddef(nc_id)) ! ENREGISTREMENT DES VARIABLES/VECTEURS DIMENSIONS call check(nf90_put_var(nc_id, varid1, lonmod)) call check(nf90_put_var(nc_id, varid2, latmod)) call check(nf90_put_var(nc_id, varid3, altmod)) call check(nf90_put_var(nc_id, varid4, tempmod)) call check(nf90_put_var(nc_id, varid5, var5)) call check(nf90_close(nc_id)) end subroutine create_3d subroutine create_2d(fname,vname,ndim1,ndim2,dim1,dim2,var5,lonmod,latmod) use netcdf implicit none character(LEN=*) :: fname,ndim1,ndim2,vname integer,parameter :: ndims=3 real,parameter :: nan=-9999 character*3 :: ndim3='day' integer,parameter :: dim3=30,dim4=4 integer :: i,dim(ndims),nc_id,dim1,dim2 integer :: varid1, varid2,time_varid,varid3,varid4,varid5,varid6,varid7,varid8 integer :: dimid1,dimid2,dimid3,dimid4, time_dimid real*4,dimension(dim1,dim2,dim3,dim4) :: var5 real,dimension(dim1) :: lonmod real,dimension(dim2) :: latmod real,dimension(dim3) :: tempmod tempmod(:)=0. do i=1,30 tempmod(i)=i enddo ! CREATION DU FICHIER nf90_create call check(nf90_create(fname, NF90_CLOBBER, & nc_id)) call check(nf90_put_att(nc_id, NF90_GLOBAL, 'Description', & 'LMDZ CALIOP Orbit files')) call check(nf90_put_att(nc_id, NF90_GLOBAL, 'Version','caliporbit2')) call check(nf90_put_att(nc_id, NF90_GLOBAL, 'Author', & 'Gregory CESANA; Cesana et al., 2015, Evaluation of 3 CALIPSO & Cloud Phase Products Using In-Situ Airborne Measurements')) ! DEFINITION DES DIMENSIONS call check(nf90_def_dim(nc_id,ndim1, dim1, dimid1)) call check(nf90_def_dim(nc_id,ndim2, dim2, dimid2 )) call check(nf90_def_dim(nc_id,ndim3, NF90_UNLIMITED, dimid3 )) ! call check(nf90_def_dim(nc_id, 'time', NF90_UNLIMITED, time_dimid)) ! DIMENSION DES MATRICES dim = (/dimid1,dimid2, dimid3/) ! VARIABLE 4D +time call check(nf90_def_var(nc_id, trim(vname)//'tot', NF90_FLOAT, dim, varid5)) call check(nf90_put_att(nc_id, varid5, '_FillValue',nan)) call check(nf90_def_var(nc_id, trim(vname)//'high', NF90_FLOAT, dim, varid6)) call check(nf90_put_att(nc_id, varid6, '_FillValue',nan)) call check(nf90_def_var(nc_id, trim(vname)//'mid', NF90_FLOAT, dim, varid7)) call check(nf90_put_att(nc_id, varid7, '_FillValue',nan)) call check(nf90_def_var(nc_id, trim(vname)//'low', NF90_FLOAT, dim, varid8)) call check(nf90_put_att(nc_id, varid8, '_FillValue',nan)) call check(nf90_def_var(nc_id, ndim1, NF90_FLOAT, dimid1, varid1)) call check(nf90_def_var(nc_id, ndim2, NF90_FLOAT, dimid2, varid2)) call check(nf90_def_var(nc_id, ndim3, NF90_FLOAT, dimid3, varid3)) ! FIN DU MODE DEFINITION call check(nf90_enddef(nc_id)) ! ENREGISTREMENT DES VARIABLES/VECTEURS DIMENSIONS call check(nf90_put_var(nc_id, varid1, lonmod)) call check(nf90_put_var(nc_id, varid2, latmod)) call check(nf90_put_var(nc_id, varid3, tempmod)) call check(nf90_put_var(nc_id, varid5, var5(:,:,:,1))) call check(nf90_put_var(nc_id, varid6, var5(:,:,:,2))) call check(nf90_put_var(nc_id, varid7, var5(:,:,:,3))) call check(nf90_put_var(nc_id, varid8, var5(:,:,:,4))) call check(nf90_close(nc_id)) end subroutine create_2d end program