program reshape_XIOS_output !======================================================================= ! ! Purpose: Read XIOS files, and convert them into the correct GCM grid ! XIOS longitudes start at -180 but stop before -180 (not duplicated) ! We basically add the last point, and complete the XIOS file. Looped ! over the two GCM runs ! ! Authors: RV & LL !======================================================================= use netcdf implicit none integer :: status, ncid, ncid1, ncid2 integer :: nDims, nVars, nGlobalAtts, unlimDimID integer i,j integer :: include_parents integer, dimension(:),allocatable :: dimids integer, dimension(:),allocatable :: varids integer, dimension(:),allocatable :: dimids_2 integer, dimension(:),allocatable :: varids_2 integer, dimension(:),allocatable :: dimid_var real, dimension(:), allocatable :: tempvalues_1d real, dimension(:), allocatable :: values_1d real, dimension(:,:), allocatable :: tempvalues_2d real, dimension(:,:), allocatable :: values_2d real, dimension(:,:,:), allocatable :: tempvalues_3d real, dimension(:,:,:), allocatable :: values_3d real, dimension(:,:,:,:), allocatable :: tempvalues_4d real, dimension(:,:,:,:), allocatable :: values_4d character*1 str2 character*30 :: name_ character*30 :: namevar integer :: xtype_var integer :: len_ integer :: len_1,len_2 integer :: len_lat, len_lon, len_time, len_soil integer :: dimid_lon, dimid_lat, dimid_time, dimid_soil integer :: dimid_2 integer :: numdims integer :: numatts integer :: numyear DO numyear=1, 2 write(*,*) 'numyear',numyear write(str2(1:1),'(i1.1)') numyear !nf90_open ! open existing netCDF dataset !integer :: ncid, status !... status = nf90_open(path = "Xdiurnalave"//str2//".nc", mode = nf90_nowrite, ncid = ncid1) if(status /= nf90_noerr) call handle_err(status) status = nf90_create(path = "Xdiurnalave"//str2//".nc_new", cmode=or(nf90_noclobber,nf90_64bit_offset), ncid = ncid2) if(status /= nf90_noerr) call handle_err(status) status = nf90_inquire(ncid1, ndims, nvars, nglobalatts, unlimdimid) if(status /= nf90_noerr) call handle_err(status) allocate(dimids(ndims)) allocate(varids(nvars)) allocate(dimids_2(ndims)) allocate(varids_2(nvars)) status = nf90_inq_dimids(ncid1, ndims, dimids, include_parents) if(status /= nf90_noerr) call handle_err(status) status = nf90_inq_varids(ncid1, nvars, varids) if(status /= nf90_noerr) call handle_err(status) do i=1,ndims status = nf90_inquire_dimension(ncid1, dimids(i), name_, len_) if(status /= nf90_noerr) call handle_err(status) if(name_.eq."lon") then dimid_lon=dimids(i) len_lon=len_ len_=len_+1 elseif(name_.eq."lat") then dimid_lat=dimids(i) len_lat=len_ elseif(name_.eq."time_counter") then dimid_time=dimids(i) len_time=len_ elseif(name_.eq."soil_layers") then dimid_soil=dimids(i) len_soil=len_ endif status = nf90_def_dim(ncid2, name_, len_, dimid_2) if(status /= nf90_noerr) call handle_err(status) dimids_2(i)=dimid_2 enddo allocate(tempvalues_3d(len_lon,len_lat,len_time)) allocate(values_3d(len_lon+1,len_lat,len_time)) allocate(tempvalues_4d(len_lon,len_lat,len_soil,len_time)) allocate(values_4d(len_lon+1,len_lat,len_soil,len_time)) do i=1,nvars status = nf90_inquire_variable(ncid1, varids(i), name=namevar, xtype=xtype_var, ndims = numdims,natts = numatts) print *, "namevar00= ", namevar if(status /= nf90_noerr) call handle_err(status) allocate(dimid_var(numdims)) status = nf90_inquire_variable(ncid1, varids(i), name=namevar, xtype=xtype_var, ndims = numdims, dimids=dimid_var, natts = numatts) if(status /= nf90_noerr) call handle_err(status) if(numdims.eq.1) then if(namevar.eq."lon") then allocate(tempvalues_1d(len_lon)) allocate(values_1d(len_lon+1)) status = nf90_get_var(ncid1, varids(i), tempvalues_1d) if(status /= nf90_noerr) call handle_err(status) status = nf90_def_var(ncid2, namevar, xtype_var, dimid_var, varids_2(i)) if(status /= nf90_noerr) call handle_err(status) values_1d(1:len_lon)=tempvalues_1d(:) values_1d(len_lon+1)=values_1d(1) status = nf90_enddef(ncid2) if(status /= nf90_noerr) call handle_err(status) status = nf90_put_var(ncid2, varids_2(i), values_1d) if(status /= nf90_noerr) call handle_err(status) status = nf90_redef(ncid2) if(status /= nf90_noerr) call handle_err(status) deallocate(tempvalues_1d) deallocate(values_1d) else status = nf90_inquire_dimension(ncid1, dimid_var(1), name_, len_) if(status /= nf90_noerr) call handle_err(status) allocate(tempvalues_1d(len_)) status = nf90_get_var(ncid1, varids(i), tempvalues_1d) if(status /= nf90_noerr) call handle_err(status) status = nf90_def_var(ncid2, namevar, xtype_var, dimid_var, varids_2(i)) if(status /= nf90_noerr) call handle_err(status) status = nf90_enddef(ncid2) if(status /= nf90_noerr) call handle_err(status) status = nf90_put_var(ncid2, varids_2(i), tempvalues_1d) if(status /= nf90_noerr) call handle_err(status) status = nf90_redef(ncid2) if(status /= nf90_noerr) call handle_err(status) deallocate(tempvalues_1d) endif elseif(numdims.eq.2) then if(namevar.eq."area") then allocate(tempvalues_2d(len_lon,len_lat)) allocate(values_2d(len_lon+1,len_lat)) status = nf90_get_var(ncid1, varids(i), tempvalues_2d) if(status /= nf90_noerr) call handle_err(status) status = nf90_def_var(ncid2, namevar, xtype_var, dimid_var, varids_2(i)) if(status /= nf90_noerr) call handle_err(status) values_2d(1:len_lon,:)=tempvalues_2d(:,:) values_2d(len_lon+1,:)=values_2d(1,:) status = nf90_enddef(ncid2) if(status /= nf90_noerr) call handle_err(status) status = nf90_put_var(ncid2, varids_2(i), values_2d) if(status /= nf90_noerr) call handle_err(status) status = nf90_redef(ncid2) if(status /= nf90_noerr) call handle_err(status) deallocate(tempvalues_2d) deallocate(values_2d) else status = nf90_inquire_dimension(ncid1, dimid_var(1), name_, len_1) if(status /= nf90_noerr) call handle_err(status) status = nf90_inquire_dimension(ncid1, dimid_var(2), name_, len_2) if(status /= nf90_noerr) call handle_err(status) allocate(tempvalues_2d(len_1,len_2)) status = nf90_get_var(ncid1, varids(i), tempvalues_2d) if(status /= nf90_noerr) call handle_err(status) status = nf90_def_var(ncid2, namevar, xtype_var, dimid_var, varids_2(i)) if(status /= nf90_noerr) call handle_err(status) status = nf90_enddef(ncid2) if(status /= nf90_noerr) call handle_err(status) status = nf90_put_var(ncid2, varids_2(i), tempvalues_2d) if(status /= nf90_noerr) call handle_err(status) status = nf90_redef(ncid2) if(status /= nf90_noerr) call handle_err(status) deallocate(tempvalues_2d) endif elseif(numdims.eq.3) then status = nf90_get_var(ncid1, varids(i), tempvalues_3d) if(status /= nf90_noerr) call handle_err(status) status = nf90_def_var(ncid2, namevar, xtype_var, dimid_var, varids_2(i)) if(status /= nf90_noerr) call handle_err(status) values_3d(1:len_lon,:,:)=tempvalues_3d(:,:,:) values_3d(len_lon+1,:,:)=values_3d(1,:,:) status = nf90_enddef(ncid2) if(status /= nf90_noerr) call handle_err(status) status = nf90_put_var(ncid2, varids_2(i), values_3d) if(status /= nf90_noerr) call handle_err(status) status = nf90_redef(ncid2) if(status /= nf90_noerr) call handle_err(status) elseif(numdims.eq.4) then status = nf90_get_var(ncid1, varids(i), tempvalues_4d) if(status /= nf90_noerr) call handle_err(status) status = nf90_def_var(ncid2, namevar, xtype_var, dimid_var, varids_2(i)) if(status /= nf90_noerr) call handle_err(status) status = nf90_enddef(ncid2) values_4d(1:len_lon,:,:,:)=tempvalues_4d(:,:,:,:) values_4d(len_lon+1,:,:,:)=values_4d(1,:,:,:) if(status /= nf90_noerr) call handle_err(status) status = nf90_put_var(ncid2, varids_2(i), values_4d) if(status /= nf90_noerr) call handle_err(status) status = nf90_redef(ncid2) if(status /= nf90_noerr) call handle_err(status) endif deallocate(dimid_var) enddo status = nf90_enddef(ncid2) if(status /= nf90_noerr) call handle_err(status) status = nf90_close(ncid1) if(status /= nf90_noerr) call handle_err(status) status = nf90_close(ncid2) if(status /= nf90_noerr) call handle_err(status) deallocate(dimids) deallocate(varids) deallocate(dimids_2) deallocate(varids_2) deallocate(tempvalues_3d) deallocate(values_3d) deallocate(tempvalues_4d) deallocate(values_4d) enddo end program reshape_XIOS_output