MODULE xios_data !----------------------------------------------------------------------- ! NAME ! xios_data ! ! DESCRIPTION ! Read XIOS output data and process it for PEM initialization. ! ! AUTHORS & DATE ! JB Clement, 2025 ! ! NOTES ! !----------------------------------------------------------------------- ! DEPENDENCIES ! ------------ use netcdf, only: nf90_open, nf90_close, nf90_inquire_dimension, nf90_inq_dimid, nf90_noerr, nf90_nowrite, nf90_get_var, nf90_inq_varid ! DECLARATION ! ----------- implicit none ! MODULE VARIABLES ! ---------------- character(19), parameter :: file1_daily = "Xoutdaily4pem_Y1.nc" ! XIOS daily output file, year 1 character(19), parameter :: file2_daily = "Xoutdaily4pem_Y2.nc" ! XIOS daily output file, year 2 character(20), parameter :: file1_yearly = "Xoutyearly4pem_Y1.nc" ! XIOS yearly output file, year 1 character(20), parameter :: file2_yearly = "Xoutyearly4pem_Y2.nc" ! XIOS yearly output file, year 2 character(256) :: msg ! Message for reading errors integer :: fID, vID ! File and variable IDs for reading ! INTERFACES ! ---------- interface get_var module procedure get_var_1d, get_var_2d, get_var_3d, get_var_4d end interface get_var contains !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ !======================================================================= SUBROUTINE load_xios_data(ngrid,nslope,nsoil_PCM,nsol,h2ofrost_PCM,co2frost_PCM,ps_avg,tsurf_avg,tsurf_avg_y1,tsoil_avg,tsoil_ts,watersurf_density_avg,d_h2oice,d_co2ice, & ps_ts,q_h2o_ts,q_co2_ts,watersoil_density_ts,min_h2oice,min_co2ice) !----------------------------------------------------------------------- ! NAME ! load_xios_data ! ! DESCRIPTION ! Reads yearly and daily XIOS data, computes frost and ice tendencies. ! ! AUTHORS & DATE ! JB Clement, 2025 ! ! NOTES ! !----------------------------------------------------------------------- ! DEPENDENCIES ! ------------ use grid_conversion, only: lonlat2vect use soil, only: do_soil use tendencies, only: compute_tend use metamorphism, only: compute_frost ! DECLARATION ! ----------- implicit none ! ARGUMENTS ! --------- integer, intent(in) :: ngrid, nslope, nsoil_PCM, nsol ! Grid dimensions real, dimension(ngrid,nslope), intent(in) :: h2ofrost_PCM, co2frost_PCM ! PCM frost fields real, dimension(ngrid), intent(out) :: ps_avg ! Average surface pressure real, dimension(ngrid,nslope), intent(out) :: tsurf_avg, tsurf_avg_y1 ! Surface temperature real, dimension(ngrid,nslope), intent(out) :: watersurf_density_avg ! Water density real, dimension(ngrid,nslope), intent(out) :: d_h2oice, d_co2ice ! Ice tendencies real, dimension(ngrid,nslope), intent(out) :: min_h2oice, min_co2ice ! Ice minima real, dimension(ngrid,nsoil_PCM,nslope), intent(out) :: tsoil_avg ! Soil temperature real, dimension(ngrid,nsol), intent(out) :: ps_ts, q_h2o_ts, q_co2_ts ! Time series real, dimension(ngrid,nsoil_PCM,nslope,nsol), intent(out) :: tsoil_ts, watersoil_density_ts ! Soil time series ! LOCAL VARIABLES ! --------------- integer :: islope, isoil, isol, nlon, nlat real, dimension(:,:), allocatable :: var_read_2d real, dimension(:,:,:), allocatable :: var_read_3d real, dimension(:,:,:,:), allocatable :: var_read_4d character(:), allocatable :: num ! For reading slope variables real, dimension(ngrid,nslope,2) :: min_h2operice, min_co2perice, min_h2ofrost, min_co2frost ! CODE ! ---- ! Initialization min_h2operice = 0. min_co2perice = 0. min_h2ofrost = 0. min_co2frost = 0. if (nslope == 1) then ! No slope allocate(character(0) :: num) else ! Using slopes allocate(character(8) :: num) endif !~~~~~~~~~~~~~~~~~~~~~~~~ Year 1 - Yearly data ~~~~~~~~~~~~~~~~~~~~~~~~~ ! Open the NetCDF file of XIOS outputs write(*,*) "Opening "//file1_yearly//"..." call error_msg(nf90_open(file1_yearly,nf90_nowrite,fID),'open',file1_yearly) ! Get the dimensions call error_msg(nf90_inq_dimid(fID,'lon',vID),'inq',file1_yearly) call error_msg(nf90_inquire_dimension(fID,vID,len = nlon),'inq',file1_yearly) call error_msg(nf90_inq_dimid(fID,'lat',vID),'inq',file1_yearly) call error_msg(nf90_inquire_dimension(fID,vID,len = nlat),'inq',file1_yearly) ! Allocate and read the variables allocate(var_read_2d(nlon,nlat),var_read_3d(nlon,nlat,nsoil_PCM)) do islope = 1,nslope if (nslope /= 1) then num=' ' write(num,'(i2.2)') islope num = '_slope'//num endif call get_var('co2ice'//num,var_read_2d) ; call lonlat2vect(nlon,nlat,ngrid,var_read_2d,min_co2frost(:,islope,1)) call get_var('h2o_ice_s'//num,var_read_2d) ; call lonlat2vect(nlon,nlat,ngrid,var_read_2d,min_h2ofrost(:,islope,1)) call get_var('watercap'//num,var_read_2d) ; call lonlat2vect(nlon,nlat,ngrid,var_read_2d,min_h2operice(:,islope,1)) call get_var('perennial_co2ice'//num,var_read_2d); call lonlat2vect(nlon,nlat,ngrid,var_read_2d,min_co2perice(:,islope,1)) call get_var('tsurf'//num,var_read_2d) ; call lonlat2vect(nlon,nlat,ngrid,var_read_2d,tsurf_avg_y1(:,islope)) enddo ! Close the NetCDF file of XIOS outputs call error_msg(nf90_close(fID),"close",file1_yearly) write(*,*) '> Data from "'//file1_yearly//'" downloaded!' !~~~~~~~~~~~~~~~~~~~~~~~~ Year 2 - Yearly data ~~~~~~~~~~~~~~~~~~~~~~~~~ ! Open the NetCDF file of XIOS outputs write(*,*) "Opening "//file2_yearly//"..." call error_msg(nf90_open(file2_yearly,nf90_nowrite,fID),'open',file2_yearly) ! Allocate and read the variables call get_var('ps',var_read_2d); call lonlat2vect(nlon,nlat,ngrid,var_read_2d,ps_avg) do islope = 1,nslope if (nslope /= 1) then num=' ' write(num,'(i2.2)') islope num = '_slope'//num endif call get_var('tsurf'//num,var_read_2d) ; call lonlat2vect(nlon,nlat,ngrid,var_read_2d,tsurf_avg(:,islope)) call get_var('co2ice'//num,var_read_2d) ; call lonlat2vect(nlon,nlat,ngrid,var_read_2d,min_co2frost(:,islope,2)) call get_var('h2o_ice_s'//num,var_read_2d) ; call lonlat2vect(nlon,nlat,ngrid,var_read_2d,min_h2ofrost(:,islope,2)) call get_var('watercap'//num,var_read_2d) ; call lonlat2vect(nlon,nlat,ngrid,var_read_2d,min_h2operice(:,islope,2)) call get_var('perennial_co2ice'//num,var_read_2d); call lonlat2vect(nlon,nlat,ngrid,var_read_2d,min_co2perice(:,islope,2)) if (do_soil) then call get_var('soiltemp'//num,var_read_3d) do isoil = 1,nsoil_PCM call lonlat2vect(nlon,nlat,ngrid,var_read_3d(:,:,isoil),tsoil_avg(:,isoil,islope)) enddo call get_var('waterdensity_surface'//num,var_read_2d); call lonlat2vect(nlon,nlat,ngrid,var_read_2d,watersurf_density_avg(:,islope)) endif enddo ! Close the NetCDF file of XIOS outputs call error_msg(nf90_close(fID),"close",file1_yearly) write(*,*) '> Data from "'//file1_yearly//'" downloaded!' !~~~~~~~~~~~~~~~~~~~~~~~~~ Year 2 - Daily data ~~~~~~~~~~~~~~~~~~~~~~~~~ ! Open the NetCDF file of XIOS outputs write(*,*) "Opening "//file2_daily//"..." call error_msg(nf90_open(file2_daily,nf90_nowrite,fID),'open',file2_daily) ! Allocate and read the variables deallocate(var_read_2d,var_read_3d) allocate(var_read_3d(nlon,nlat,nsol),var_read_4d(nlon,nlat,nsoil_PCM,nsol)) call get_var('ps',var_read_3d) do isol = 1,nsoil_PCM call lonlat2vect(nlon,nlat,ngrid,var_read_3d(:,:,isol),ps_ts(:,isol)) enddo call get_var('h2o_layer1',var_read_3d) do isol = 1,nsoil_PCM call lonlat2vect(nlon,nlat,ngrid,var_read_3d(:,:,isol),q_h2o_ts(:,isol)) enddo call get_var('co2_layer1',var_read_3d) do isol = 1,nsoil_PCM call lonlat2vect(nlon,nlat,ngrid,var_read_3d(:,:,isol),q_co2_ts(:,isol)) enddo if (do_soil) then do islope = 1,nslope if (nslope /= 1) then num=' ' write(num,'(i2.2)') islope num = '_slope'//num endif call get_var('soiltemp'//num,var_read_4d) do isol = 1,nsol do isoil = 1,nsoil_PCM call lonlat2vect(nlon,nlat,ngrid,var_read_4d(:,:,isoil,isol),tsoil_ts(:,isoil,islope,isol)) enddo enddo call get_var('waterdensity_soil'//num,var_read_4d) do isol = 1,nsol do isoil = 1,nsoil_PCM call lonlat2vect(nlon,nlat,ngrid,var_read_4d(:,:,isoil,isol),watersoil_density_ts(:,isoil,islope,isol)) enddo enddo enddo endif deallocate(var_read_3d,var_read_4d,num) ! Close the NetCDF file of XIOS outputs call error_msg(nf90_close(fID),"close",file1_daily) write(*,*) '> Data from "'//file2_daily//'" downloaded!' !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ ! Compute frost from yearly minima call compute_frost(ngrid,nslope,h2ofrost_PCM,min_h2ofrost,co2frost_PCM,min_co2frost) !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ ! Compute ice tendencies from yearly minima write(*,*) '> Computing surface ice tendencies' call compute_tend(ngrid,nslope,min_h2operice + min_h2ofrost,d_h2oice) write(*,*) 'H2O ice tendencies (min/max):', minval(d_h2oice), maxval(d_h2oice) call compute_tend(ngrid,nslope,min_co2perice + min_co2frost,d_co2ice) write(*,*) 'CO2 ice tendencies (min/max):', minval(d_co2ice), maxval(d_co2ice) END SUBROUTINE load_xios_data !======================================================================= !======================================================================= SUBROUTINE get_timelen(filename,timelen) !----------------------------------------------------------------------- ! NAME ! get_timelen ! ! DESCRIPTION ! Get time dimension length from a NetCDF file. ! ! AUTHORS & DATE ! JB Clement, 2025 ! ! NOTES ! !----------------------------------------------------------------------- ! DEPENDENCIES ! ------------ use netcdf ! DECLARATION ! ----------- implicit none ! ARGUMENTS ! --------- character(*), intent(in) :: filename ! NetCDF filename integer, intent(out) :: timelen ! Length of time dimension ! LOCAL VARIABLES ! --------------- integer :: ncid ! File ID integer :: dimid ! Dimension ID integer :: ierr ! Return code ! CODE ! ---- ! Open the NetCDF file ierr = nf90_open(trim(filename),NF90_NOWRITE,ncid) if (ierr /= nf90_noerr) then write(*,*) "Error opening file:", trim(nf90_strerror(ierr)) error stop endif ! Get the dimension ID for 'time_counter' ierr = nf90_inq_dimid(ncid,"time_counter",dimid) if (ierr /= nf90_noerr) then write(*,*) "Error getting dimid 'time_counter':", trim(nf90_strerror(ierr)) error stop endif ! Get the size of the dimension 'time_counter' ierr = nf90_inquire_dimension(ncid,dimid,len = timelen) if (ierr /= nf90_noerr) then write(*,*) "Error getting dimension length:", trim(nf90_strerror(ierr)) error stop endif ! Close the file ierr = nf90_close(ncid) if (ierr /= nf90_noerr) then write(*,*) "Error closing file:", trim(nf90_strerror(ierr)) error stop endif END SUBROUTINE get_timelen !======================================================================= !======================================================================= SUBROUTINE error_msg(ierr,typ,nam) !----------------------------------------------------------------------- ! NAME ! error_msg ! ! DESCRIPTION ! Handle and report NetCDF errors. ! ! AUTHORS & DATE ! JB Clement, 2025 ! ! NOTES ! !----------------------------------------------------------------------- ! DECLARATION ! ----------- implicit none ! ARGUMENTS ! --------- integer, intent(in) :: ierr ! NetCDF error code character(*), intent(in) :: typ ! Type of operation (inq, get, put, open, close) character(*), intent(in) :: nam ! Field/file name ! CODE ! ---- if (ierr == nf90_noerr) return select case(typ) case('inq') ; msg = "Dim/Field <"//trim(nam)//"> is missing" case('get') ; msg = "Reading failed for <"//trim(nam)//">" case('put') ; msg = "Writing failed for <"//trim(nam)//">" case('open') ; msg = "File opening failed for <"//trim(nam)//">" case('close'); msg = "File closing failed for <"//trim(nam)//">" case default write(*,*) 'There is no message for this error.' error stop end select call abort_gcm(__FILE__,trim(msg),ierr) END SUBROUTINE error_msg !======================================================================= !======================================================================= SUBROUTINE get_var_1d(var,v) !----------------------------------------------------------------------- ! NAME ! get_var_1d ! ! DESCRIPTION ! Read a 1D variable from open NetCDF file. ! ! AUTHORS & DATE ! JB Clement, 2025 ! ! NOTES ! !----------------------------------------------------------------------- ! DECLARATION ! ----------- implicit none ! ARGUMENTS ! --------- character(*), intent(in) :: var ! Variable name real, dimension(:), intent(out) :: v ! Output array ! CODE ! ---- call error_msg(nf90_inq_varid(fID,var,vID),"inq",var) call error_msg(nf90_get_var(fID,vID,v),"get",var) END SUBROUTINE get_var_1d !======================================================================= !======================================================================= SUBROUTINE get_var_2d(var,v) !----------------------------------------------------------------------- ! NAME ! get_var_2d ! ! DESCRIPTION ! Read a 2D variable from open NetCDF file. ! ! AUTHORS & DATE ! JB Clement, 2025 ! ! NOTES ! !----------------------------------------------------------------------- ! DECLARATION ! ----------- implicit none ! ARGUMENTS ! --------- character(*), intent(in) :: var ! Variable name real, dimension(:,:), intent(out) :: v ! Output array ! CODE ! ---- call error_msg(nf90_inq_varid(fID,var,vID),"inq",var) call error_msg(nf90_get_var(fID,vID,v),"get",var) END SUBROUTINE get_var_2d !======================================================================= !======================================================================= SUBROUTINE get_var_3d(var,v) !----------------------------------------------------------------------- ! NAME ! get_var_3d ! ! DESCRIPTION ! Read a 3D variable from open NetCDF file. ! ! AUTHORS & DATE ! JB Clement, 2025 ! ! NOTES ! !----------------------------------------------------------------------- ! DECLARATION ! ----------- implicit none ! ARGUMENTS ! --------- character(*), intent(in) :: var ! Variable name real, dimension(:,:,:), intent(out) :: v ! Output array ! CODE ! ---- call error_msg(nf90_inq_varid(fID,var,vID),"inq",var) call error_msg(nf90_get_var(fID,vID,v),"get",var) END SUBROUTINE get_var_3d !======================================================================= !======================================================================= SUBROUTINE get_var_4d(var,v) !----------------------------------------------------------------------- ! NAME ! get_var_4d ! ! DESCRIPTION ! Read a 4D variable from open NetCDF file. ! ! AUTHORS & DATE ! JB Clement, 2025 ! ! NOTES ! !----------------------------------------------------------------------- ! DECLARATION ! ----------- implicit none ! ARGUMENTS ! --------- character(*), intent(in) :: var ! Variable name real, dimension(:,:,:,:), intent(out) :: v ! Output array ! CODE ! ---- call error_msg(nf90_inq_varid(fID,var,vID),"inq",var) call error_msg(nf90_get_var(fID,vID,v),"get",var) END SUBROUTINE get_var_4d !======================================================================= END MODULE xios_data