diff --git a/prep_syndata/get_flux.f90 b/prep_syndata/get_flux.f90 index bbd99509b1009b069fbb1b550516d508bb11ad23..bb88d6c89682131d48bdc57e98e661ab514a7dbe 100644 --- a/prep_syndata/get_flux.f90 +++ b/prep_syndata/get_flux.f90 @@ -46,7 +46,7 @@ subroutine get_flux(config, files, flx, state) type (config_t), intent(in) :: config type (files_t), intent(in) :: files - real, dimension(nxgrid,nygrid,nt_flx), intent(in out) :: flx + real, dimension(nxgrid,nygrid,nflx), intent(in out) :: flx real, dimension(nvar), intent(in out) :: state character(len=max_path_len) :: filename @@ -62,7 +62,7 @@ subroutine get_flux(config, files, flx, state) real, dimension(:,:,:), allocatable :: tmp real, parameter :: smallnum = 1.e-15 - allocate( flx_time(nt_flx) ) + allocate( flx_time(nflx) ) flx_time(:) = 0. ! read prior fluxes @@ -101,7 +101,7 @@ subroutine get_flux(config, files, flx, state) if ( config%spec.eq.'ghg' ) then ! GHG species - do while ( (jd.lt.juldatef).and.(n.lt.nt_flx) ) + do while ( (jd.lt.juldatef).and.(n.lt.nflx) ) n = n + 1 call caldate(jd, jjjjmmdd, hhmiss) write(adate,'(I4.4)') jjjjmmdd/10000 @@ -143,7 +143,7 @@ subroutine get_flux(config, files, flx, state) else if ( config%spec.eq.'co2' ) then ! CO2 species - do while ( (jd.le.juldatef).and.(n.lt.nt_flx) ) + do while ( (jd.le.juldatef).and.(n.lt.nflx) ) n = n + 1 call caldate(jd, jjjjmmdd, hhmiss) write(adate,'(I4.4)') jjjjmmdd/10000 diff --git a/prep_syndata/initialize.f90 b/prep_syndata/initialize.f90 index 42ced70264fb999fc5af182da9d06e0317223395..18ab77fc5eef53bcf165875d44ecae9b97dff4ac 100644 --- a/prep_syndata/initialize.f90 +++ b/prep_syndata/initialize.f90 @@ -132,15 +132,18 @@ subroutine initialize(files, config) flxtres = aint(365./real(config%nt_flx)) ntstate = nint(real(nday)/statres) nt_flx = nday*nd_nee + nflx = nt_flx else ! statres is time resolution of state variables ! flxtres is temporal resolution of fluxes (same as state variables for non-CO2 tracers) + ! nflx is the number of input flux timesteps falling within the inversion interval ! ntstate is also number of state time invervals statres = config%statres - flxtres = config%statres + flxtres = aint(365./real(config%nt_flx)) ntstate = nint(real(nday)/statres) ndt = 1 nt_flx = config%nt_flx + nflx = max(1, nday/int(365./real(config%nt_flx))) endif if ( ntstate.lt.1 ) then write(logid,*) 'ERROR initialize: number flux time steps < 1 -> increase time interval' @@ -155,6 +158,7 @@ subroutine initialize(files, config) print*, 'flxtres: ',flxtres print*, 'ntstate: ',ntstate print*, 'nt_flx: ',nt_flx + print*, 'nflx: ',nflx print*, 'nday: ',nday print*, 'nyr: ',nyr diff --git a/prep_syndata/mod_save.f90 b/prep_syndata/mod_save.f90 index 84d34072acf1f2b00a4b8a66f2376aa2f5b7dd72..d63a2f20da88537d030c5330e4831802faead529 100644 --- a/prep_syndata/mod_save.f90 +++ b/prep_syndata/mod_save.f90 @@ -60,7 +60,7 @@ module mod_save type (config_t), intent (in) :: config type (files_t), intent (in) :: files - real, dimension(nxgrid,nygrid,nt_flx), intent (in out) :: flx + real, dimension(nxgrid,nygrid,nflx), intent (in out) :: flx real, dimension(nvar), intent (in out) :: state character(len=max_name_len) :: filename, lonname, latname, timename, varname @@ -75,8 +75,8 @@ module mod_save integer :: d, n, nd, nn, nt, ii, i1, i2 integer :: ix, jy, ix1, ix2, jy1, jy2 real :: frac, hh, hutc - real, dimension(nxregrid,nyregrid,nt_flx) :: flx_out - real, dimension(nt_flx) :: time + real, dimension(nxregrid,nyregrid,nflx) :: flx_out + real, dimension(nflx) :: time real, dimension(:,:,:), allocatable :: tmp real(kind=8) :: jd integer :: yyyymmdd, hhmiss @@ -149,9 +149,10 @@ module mod_save if ( config%spec.eq.'ghg' ) then ! GHG species - jd = juldatei - n = 1 - do while ( (jd.lt.juldatef).and.(n.lt.nt_flx) ) +! jd = juldatei +! n = 1 +! do while ( (jd.lt.juldatef).and.(n.lt.nflx) ) + do n = 1, nflx nt = minloc(stat_time,dim=1,mask=abs(stat_time-flx_time(n)).eq.minval(abs(stat_time-flx_time(n)))) nt = min(nt,ntstate) print*, 'stat_time(nt), flx_time(n) = ',stat_time(nt), flx_time(n) @@ -208,14 +209,14 @@ module mod_save endif end do end do - jd = jd + statres - n = n + 1 end do endif ! GHG - print*, 'mod_save: n = ',n ! replace fluxes in inversion domain with perturbed fluxes - flx(ix1:ix2,jy1:jy2,1:n) = flx_out(:,:,1:n) + flx(ix1:ix2,jy1:jy2,:) = flx_out(:,:,:) + +!*** FOR SF6 + where ( flx(:,:,:).lt.0. ) flx(:,:,:) = 0. ! undo numerical scaling flx = flx/numscale diff --git a/prep_syndata/mod_var.f90 b/prep_syndata/mod_var.f90 index 95c68f81949e3d0f536b79983bc4ff74e257158b..8a5f3d0cf56d8772d4c36c961f6fd48a0050c3bd 100644 --- a/prep_syndata/mod_var.f90 +++ b/prep_syndata/mod_var.f90 @@ -31,7 +31,7 @@ module mod_var integer, parameter :: max_path_len=200 ! max character length of paths and files integer, parameter :: max_name_len=50 ! max character length of variable names integer, parameter :: recname_len=3 ! length of receptor names - real, parameter :: numscale=1.e20 ! numeric scaling factor + real, parameter :: numscale=1.e12 ! numeric scaling factor integer, parameter :: maxpoint=1000 ! max number of releases per receptor per month integer, parameter :: maxspec=10 ! max number of species in a flexpart run integer, parameter :: maxlev=50 ! max number of vertical levels in flexpart output @@ -67,6 +67,7 @@ module mod_var integer :: nyr ! number of years in inversion interval integer :: nd_nee ! number of NEE values per day (in prior NEE file) integer :: nt_flx ! number of flux values per year in prior flux files + integer :: nflx ! number of input flux time steps falling within inversion interval integer :: ntstate ! number of time steps in state vector real :: statres ! time resolution of the state vector (for CO2 1 day) (days) real :: statres_hr ! sub-daily time resolution of the state vector (hours) diff --git a/prep_syndata/prep_syndata.f90 b/prep_syndata/prep_syndata.f90 index ece75880352a873fdfab42f1a5bd05a93b59ef80..3e30d7d125183a3d5396645d45fa923eb21d187a 100644 --- a/prep_syndata/prep_syndata.f90 +++ b/prep_syndata/prep_syndata.f90 @@ -60,8 +60,8 @@ program prep_syndata stop 'ERROR: need to specify SETTINGS_config file' endif - write(logid,*) 'Using: '//trim(settings_files) - write(logid,*) 'Using: '//trim(settings_config) + write(*,*) 'Using: '//trim(settings_files) + write(*,*) 'Using: '//trim(settings_config) ! initialization ! -------------- @@ -79,7 +79,7 @@ program prep_syndata call get_error_cov(files, config) ! get prior fluxes - allocate( flx(nxgrid,nygrid,nt_flx)) + allocate( flx(nxgrid,nygrid,nflx)) allocate( state(nvar) ) flx(:,:,:) = 0. state(:) = 0.