To all users with # define PSOURCE_NCFILE, there is a small typo that makes your river constant (and not variable)
you want to change in analytical.F:
Qbar(is)=rampQbar0(is)
by
Qbar(is)=rampQbar(is)
Hi all,
No there is not a bug.
This is in the analytical.F routine, where there is no temporal variation. This is used in with the PSOURCE and ANA_SOURCE cpp keys. You use this variable also in case of PSOURCE_NCFILE only IF there is not data and you need a constant analytical value.
Variable discharge and/or concentration are done using PSOURCE_NCFILE and PSOURCE_NCFILE_TS and are related to the routines get_psource.F and get_psource_ts.F.
In crocotools, the variable (time and/or concentration) runoff netcdf files can be created using the make_runoff.m routines.
Cheers,
Gildas
Hi all,
No there is not a bug.
This is in the analytical.F routine, where there is no temporal variation. This is used in with the PSOURCE and ANA_SOURCE cpp keys. You use this variable also in case of PSOURCE_NCFILE only IF there is not data and you need a constant analytical value.
Variable discharge and/or concentration are done using PSOURCE_NCFILE and PSOURCE_NCFILE_TS and are related to the routines get_psource.F and get_psource_ts.F.
In crocotools, the variable (time and/or concentration) runoff netcdf files can be created using the make_runoff.m routines.
Cheers,
Gildas
Hi Gildas,
Thank you for the answer,
If I am not mistaken, ANA_PSOURCE will always be defined if you use either PSOURCE or PSOURCE_NCFILE (from cppdefs_dev.h) because it needs to compute the new Qbar using Qshape.
At each time step, Qsrc(is,k)=Qbar(is)*Qshape(is,k) will be read from analytical.F
The issue is that the line above this is:
Qbar(is)=ramp*Qbar0(is) (which is Qbar0(is) if ramp=1 by default) which will always be a constant (set above when icc.eq.start).
By replacing Qbar(is)=ramp*Qbar(is) (instead of Qbar0(is), my river discharge variability works perfectly)
Vincent
This is a complementary answer to Vincent and Gildas. Vincent is right that Qsrc still needs to be defined with the correct vertical shape (Qshape) from Qbar in ana_psource. There is thus a problem in ana_psource with Qbar(is)=rampQbar0(is) because Qbar0 is only given at the start of the simulation. This was done for the constant forcing (analytical) case and was not adapted to the case where a variable Qbar is read in a file (PSOURCE_NCFILE). However, setting Qbar(is)=rampQbar(is) is not the correct fix if ramp is not 1 because Qbar is then wrongly calculated as an iterative procedure (and underestimated). In a development branch that will soon be released, we recently fixed this issue by setting outside of the (iic.eq.ntstart) condition (and before the use of ramp):
!
! If Qbar is read in nc file, update Qbar0 at every time step
!
ifdef PSOURCE_NCFILE
do is=1,Nsrc
Qbar0(is)=Qbar(is)
enddo
endif
Patrick