PROGRAM reshape_XIOS_output !======================================================================= ! ! Purpose: Read XIOS files, and convert them into the correct PCM 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 PCM runs ! ! Authors: RV & LL !======================================================================= use netcdf implicit none integer :: state, ncid1, ncid2, nDims, nVars, nGlobalAtts, unlimDimID integer :: i, include_parents, cstat integer, dimension(:), allocatable :: dimids, varids, dimids_2, varids_2, dimid_var real, dimension(:), allocatable :: tempvalues_1d, values_1d real, dimension(:,:), allocatable :: tempvalues_2d, values_2d real, dimension(:,:,:), allocatable :: tempvalues_3d, values_3d real, dimension(:,:,:,:), allocatable :: tempvalues_4d, values_4d character(1) :: str character(30) :: name_, namevar integer :: xtype_var, len_, len_1, len_2, len_lat, len_lon, len_time, len_soil integer :: dimid_lon, dimid_lat, dimid_time, dimid_soil, dimid_2, numdims, numatts, numyear logical :: yes do numyear = 1,2 write(*,*) 'numyear',numyear write(str(1:1),'(i1.1)') numyear state = nf90_open(path = "data2reshape_Y"//str//".nc", mode = nf90_nowrite, ncid = ncid1) if (state /= nf90_noerr) call handle_err(state) inquire(file = 'data_PCM_Y'//str//'.nc', exist = yes) if (yes) then call execute_command_line('rm data_PCM_Y'//str//'.nc',cmdstat = cstat) if (cstat > 0) then error stop 'Command exection failed!' else if (cstat < 0) then error stop 'Command execution not supported!' endif endif state = nf90_create(path = "data_PCM_Y"//str//".nc", cmode=or(nf90_noclobber,nf90_64bit_offset), ncid = ncid2) if (state /= nf90_noerr) call handle_err(state) state = nf90_inquire(ncid1, ndims, nvars, nglobalatts, unlimdimid) if (state /= nf90_noerr) call handle_err(state) allocate(dimids(ndims)) allocate(varids(nvars)) allocate(dimids_2(ndims)) allocate(varids_2(nvars)) state = nf90_inq_dimids(ncid1, ndims, dimids, include_parents) if (state /= nf90_noerr) call handle_err(state) state = nf90_inq_varids(ncid1, nvars, varids) if (state /= nf90_noerr) call handle_err(state) do i = 1,ndims state = nf90_inquire_dimension(ncid1, dimids(i), name_, len_) if (state /= nf90_noerr) call handle_err(state) if (name_ == "lon" .or. name_ == "longitude") then dimid_lon = dimids(i) len_lon = len_ len_ = len_ + 1 else if (name_ == "lat".or. name_ == "latitude") then dimid_lat=dimids(i) len_lat=len_ else if (name_ == "time_counter".or. name_ == "Time") then dimid_time=dimids(i) len_time=len_ else if (name_ == "soil_layers".or. name_ == "subsurface_layers") then dimid_soil=dimids(i) len_soil = len_ endif state = nf90_def_dim(ncid2, name_,len_,dimid_2) if (state /= nf90_noerr) call handle_err(state) dimids_2(i) = dimid_2 enddo do i = 1,nvars state = nf90_inquire_variable(ncid1, varids(i),name = namevar,xtype = xtype_var,ndims = numdims,natts = numatts) write(*,*) "namevar00= ", namevar if (state /= nf90_noerr) call handle_err(state) allocate(dimid_var(numdims)) state = nf90_inquire_variable(ncid1,varids(i),name = namevar,xtype = xtype_var,ndims = numdims,dimids = dimid_var,natts = numatts) if (state /= nf90_noerr) call handle_err(state) if (numdims == 1) then if (namevar == "lon") then allocate(tempvalues_1d(len_lon)) allocate(values_1d(len_lon + 1)) state = nf90_get_var(ncid1,varids(i),tempvalues_1d) if (state /= nf90_noerr) call handle_err(state) state = nf90_def_var(ncid2,namevar,xtype_var, dimid_var, varids_2(i)) if (state /= nf90_noerr) call handle_err(state) values_1d(1:len_lon) = tempvalues_1d(:) values_1d(len_lon + 1) = values_1d(1) state = nf90_enddef(ncid2) if (state /= nf90_noerr) call handle_err(state) state = nf90_put_var(ncid2, varids_2(i), values_1d) if (state /= nf90_noerr) call handle_err(state) state = nf90_redef(ncid2) if (state /= nf90_noerr) call handle_err(state) deallocate(tempvalues_1d) deallocate(values_1d) else state = nf90_inquire_dimension(ncid1,dimid_var(1),name_,len_) if (state /= nf90_noerr) call handle_err(state) allocate(tempvalues_1d(len_)) state = nf90_get_var(ncid1,varids(i),tempvalues_1d) if (state /= nf90_noerr) call handle_err(state) state = nf90_def_var(ncid2,namevar,xtype_var,dimid_var,varids_2(i)) if (state /= nf90_noerr) call handle_err(state) state = nf90_enddef(ncid2) if (state /= nf90_noerr) call handle_err(state) state = nf90_put_var(ncid2, varids_2(i), tempvalues_1d) if (state /= nf90_noerr) call handle_err(state) state = nf90_redef(ncid2) if (state /= nf90_noerr) call handle_err(state) deallocate(tempvalues_1d) endif else if (numdims == 2) then if (namevar == "area") then allocate(tempvalues_2d(len_lon,len_lat)) allocate(values_2d(len_lon + 1,len_lat)) state = nf90_get_var(ncid1,varids(i),tempvalues_2d) if (state /= nf90_noerr) call handle_err(state) state = nf90_def_var(ncid2,namevar,xtype_var,dimid_var,varids_2(i)) if (state /= nf90_noerr) call handle_err(state) values_2d(1:len_lon,:) = tempvalues_2d(:,:) values_2d(len_lon+1,:) = values_2d(1,:) state = nf90_enddef(ncid2) if (state /= nf90_noerr) call handle_err(state) state = nf90_put_var(ncid2,varids_2(i),values_2d) if (state /= nf90_noerr) call handle_err(state) state = nf90_redef(ncid2) if (state /= nf90_noerr) call handle_err(state) deallocate(tempvalues_2d) deallocate(values_2d) else state = nf90_inquire_dimension(ncid1,dimid_var(1),name_,len_1) if (state /= nf90_noerr) call handle_err(state) state = nf90_inquire_dimension(ncid1,dimid_var(2),name_,len_2) if (state /= nf90_noerr) call handle_err(state) allocate(tempvalues_2d(len_1,len_2)) state = nf90_get_var(ncid1, varids(i), tempvalues_2d) if (state /= nf90_noerr) call handle_err(state) state = nf90_def_var(ncid2,namevar,xtype_var,dimid_var,varids_2(i)) if (state /= nf90_noerr) call handle_err(state) state = nf90_enddef(ncid2) if (state /= nf90_noerr) call handle_err(state) state = nf90_put_var(ncid2, varids_2(i), tempvalues_2d) if (state /= nf90_noerr) call handle_err(state) state = nf90_redef(ncid2) if (state /= nf90_noerr) call handle_err(state) deallocate(tempvalues_2d) endif else if (numdims == 3) then allocate(tempvalues_3d(len_lon,len_lat,len_time)) allocate(values_3d(len_lon + 1,len_lat,len_time)) state = nf90_get_var(ncid1,varids(i),tempvalues_3d) if (state /= nf90_noerr) call handle_err(state) state = nf90_def_var(ncid2,namevar,xtype_var,dimid_var,varids_2(i)) if (state /= nf90_noerr) call handle_err(state) values_3d(1:len_lon,:,:) = tempvalues_3d(:,:,:) values_3d(len_lon+1,:,:) = values_3d(1,:,:) state = nf90_enddef(ncid2) if (state /= nf90_noerr) call handle_err(state) state = nf90_put_var(ncid2, varids_2(i), values_3d) if (state /= nf90_noerr) call handle_err(state) state = nf90_redef(ncid2) if (state /= nf90_noerr) call handle_err(state) deallocate(tempvalues_3d) deallocate(values_3d) else if (numdims == 4) then allocate(tempvalues_4d(len_lon,len_lat,len_soil,len_time)) allocate(values_4d(len_lon+1,len_lat,len_soil,len_time)) state = nf90_get_var(ncid1, varids(i), tempvalues_4d) if (state /= nf90_noerr) call handle_err(state) state = nf90_def_var(ncid2, namevar, xtype_var, dimid_var, varids_2(i)) if (state /= nf90_noerr) call handle_err(state) state = nf90_enddef(ncid2) values_4d(1:len_lon,:,:,:) = tempvalues_4d(:,:,:,:) values_4d(len_lon+1,:,:,:) = values_4d(1,:,:,:) if (state /= nf90_noerr) call handle_err(state) state = nf90_put_var(ncid2, varids_2(i), values_4d) if (state /= nf90_noerr) call handle_err(state) state = nf90_redef(ncid2) if (state /= nf90_noerr) call handle_err(state) deallocate(tempvalues_4d) deallocate(values_4d) endif deallocate(dimid_var) enddo state = nf90_enddef(ncid2) if (state /= nf90_noerr) call handle_err(state) state = nf90_close(ncid1) if (state /= nf90_noerr) call handle_err(state) state = nf90_close(ncid2) if (state /= nf90_noerr) call handle_err(state) deallocate(dimids,varids,dimids_2,varids_2) enddo END PROGRAM reshape_XIOS_output