MODFLOW 6  version 6.7.0.dev0
USGS Modular Hydrologic Model
obsutilitymodule Module Reference

This module contains the ObsUtilityModule module. More...

Functions/Subroutines

subroutine, public write_fmtd_obs (fmtc, obsrv, obsOutputList, value)
 @ brief Write formatted observation More...
 
subroutine, public write_unfmtd_obs (obsrv, iprec, obsOutputList, value)
 @ brief Write unformatted observation More...
 

Detailed Description

This module contains subroutines for writing simulated values stored in objects of ObserveType to output files. The subroutines handle continuous observations, and can write values to either formatted or unformatted files.

Function/Subroutine Documentation

◆ write_fmtd_obs()

subroutine, public obsutilitymodule::write_fmtd_obs ( character(len=*), intent(in)  fmtc,
type(observetype), intent(inout)  obsrv,
type(obsoutputlisttype), intent(inout), pointer  obsOutputList,
real(dp), intent(in)  value 
)

Subroutine to write observation data for the end of a time step to a formatted file. If the simulation time has not been written to for the current time step, totim is written. The simulated value is written in the format specified in the fmtc argument.

Parameters
[in]fmtcobservation format
[in,out]obsrvobservation type
[in,out]obsoutputlistobservation list
[in]valueobservation

Definition at line 33 of file ObsUtility.f90.

34  ! -- dummy
35  character(len=*), intent(in) :: fmtc !< observation format
36  type(ObserveType), intent(inout) :: obsrv !< observation type
37  type(ObsOutputListType), pointer, intent(inout) :: obsOutputList !< observation list
38  real(DP), intent(in) :: value !< observation
39  ! -- local
40  integer(I4B) :: indx
41  integer(I4B) :: nunit
42  character(len=20) :: ctotim
43  character(len=50) :: cval
44  type(ObsOutputType), pointer :: ObsOutput => null()
45  ! -- output unit
46  nunit = obsrv%UnitNumber
47  !
48  indx = obsrv%indxObsOutput
49  obsoutput => obsoutputlist%Get(indx)
50  if (obsoutput%empty_line) then
51  obsoutput%empty_line = .false.
52  write (ctotim, '(G20.13)') totim
53  else
54  ctotim = ''
55  end if
56  ! -- append value to output line
57  write (cval, fmtc) value
58  write (nunit, '(3a)', advance='NO') &
59  trim(adjustl(ctotim)), ',', trim(adjustl(cval))
60  !
61  ! -- flush the file
62  ! Added flush after each non-advancing write to resolve
63  ! issue with ifort (IFORT) 19.1.0.166 20191121 for Linux
64  ! that occurred on some Linux systems.
65  flush (nunit)
Here is the caller graph for this function:

◆ write_unfmtd_obs()

subroutine, public obsutilitymodule::write_unfmtd_obs ( type(observetype), intent(inout)  obsrv,
integer(i4b), intent(in)  iprec,
type(obsoutputlisttype), intent(inout), pointer  obsOutputList,
real(dp), intent(in)  value 
)

Subroutine to write observation data for the end of a time step to a unformatted file. If the simulation time has not been written for the current time step, totim is written. The simulated value is written using the precision specified in the iprec argument.

iprec = 1: real32 specifies 32-bit real = 4 bytes = single precision. iprec = 2: real64 specifies 64-bit real = 8 bytes = double precision.

Parameters
[in,out]obsrvobservation type
[in]iprecobservation precision
[in,out]obsoutputlistobservation list
[in]valueobservation

Definition at line 79 of file ObsUtility.f90.

80  use iso_fortran_env, only: real32, real64
81  ! -- dummy
82  type(ObserveType), intent(inout) :: obsrv !< observation type
83  integer(I4B), intent(in) :: iprec !< observation precision
84  type(ObsOutputListType), pointer, intent(inout) :: obsOutputList !< observation list
85  real(DP), intent(in) :: value !< observation
86  ! -- local
87  integer(I4B) :: indx, nunit
88  real(real32) :: totimsngl, valsngl
89  real(real64) :: totimdbl, valdbl
90  type(ObsOutputType), pointer :: obsOutput => null()
91  !
92  ! -- output unit
93  nunit = obsrv%UnitNumber
94  ! -- continuous observation
95  indx = obsrv%indxObsOutput
96  obsoutput => obsoutputlist%Get(indx)
97  if (obsoutput%empty_line) then
98  obsoutput%empty_line = .false.
99  if (iprec == 1) then
100  totimsngl = real(totim, real32)
101  write (nunit) totimsngl
102  elseif (iprec == 2) then
103  totimdbl = totim
104  write (nunit) totimdbl
105  end if
106  end if
107  ! -- write value to unformatted output
108  if (iprec == 1) then
109  valsngl = real(value, real32)
110  write (nunit) valsngl
111  elseif (iprec == 2) then
112  valdbl = value
113  write (nunit) valdbl
114  end if
Here is the caller graph for this function: