Line data Source code
1 : !> \file mo_edk_modules.f90
2 : !> \brief Modules with global variables for EDK
3 : !> \details This file contains the MAIN SHARED variables
4 : !! - `mainvar`: main setup variables
5 : !! - `runcontrol`: run configuration
6 : !! - `kriging`: kriging setup variables
7 : !! - `varfit`: variogram fitting variables
8 : !! - `netcdfvar`: NetCDF I/O definition variables
9 : !!
10 : !> \author Luis Samaniego
11 : !> \date 22.03.2006
12 : !> \date 24.03.2006
13 : !> \author Sebastian Mueller
14 : !> \date June 2022
15 : !! - refactored
16 :
17 : !> \brief main variables
18 : !> \copyright Copyright 2005-\today, the CHS Developers, Sabine Attinger: All rights reserved.
19 : !! EDK is released under the LGPLv3+ license \license_note
20 : module mainVar
21 : use mo_kind, only: i4, sp, dp
22 : implicit none
23 : ! parameters
24 : integer(i4) :: yStart !< starting year
25 : integer(i4) :: mStart !< starting month
26 : integer(i4) :: dStart !< starting day
27 : integer(i4) :: yEnd !< ending year
28 : integer(i4) :: mEnd !< ending month
29 : integer(i4) :: dEnd !< ending day
30 : integer(i4) :: jStart !< julian day start
31 : integer(i4) :: jEnd !< julian day end
32 : integer(i4) :: tBuffer !< number of days for time buffering
33 : integer(i4) :: tDays !< number of days
34 : integer(i4) :: nSta !< number of stations for a block
35 : integer(i4) :: nCell !< number of cells to estimate z
36 : integer(i4) :: cellFactor !< > 1 , size grid metereological data
37 : integer(i4) :: DEMNcFlag !< flag for DEM format 0 = text file, 1 = netCDF
38 : real(dp) :: DataConvertFactor !< precipitation & temperature(in 1/10 mm) **** only in NECKAR BASIN *****
39 : real(dp) :: OffSet !< constant to be added (Ex: add 273 to convert tavg from C to K )
40 : real(dp) :: noDataValue !< no data value of station input data
41 : real(dp) :: thresholdDist !< treshold cellsize distance
42 : ! constants
43 : real(dp), parameter :: DayHours = 24.0_dp !< hours per day
44 : real(dp), parameter :: YearDays = 365.0_dp !< days in a year
45 : real(dp), parameter :: DaySecs = 86400.0_dp !< sec in a day
46 : real(sp), public, parameter :: nodata_sp = -9999._sp !< [-] global no data value (single precision)
47 :
48 : !> \class cellfiner
49 : !> \brief cell elevation
50 : type CellFiner
51 : real(dp) :: h !< elevation (sinks removed) [m]
52 : end type CellFiner
53 : type(CellFiner), dimension(:,:), allocatable :: G !< Cell characteristics
54 :
55 : !> \class meteostation
56 : !> \brief Meteo Stations type
57 : type MeteoStation
58 : integer(i4) :: Id !< Id number
59 : real(dp) :: x !< x coordinate
60 : real(dp) :: y !< y coordinate
61 : real(dp) :: h !< elevation
62 : real(dp), dimension(:), allocatable :: z !< observed daily value (prec. temp, etc)
63 : end type MeteoStation
64 : type(MeteoStation), dimension(:), allocatable :: MetSta !< Meteo Stations
65 :
66 : !> \class gridgeoref
67 : !> \brief GRID description
68 : type gridGeoRef
69 : integer(i4) :: ncols !< number of columns
70 : integer(i4) :: nrows !< number of rows
71 : real(dp) :: xllcorner !< x coordinate of the lowerleft corner
72 : real(dp) :: yllcorner !< y coordinate of the lowerleft corner
73 : integer(i4) :: cellsize !< cellsize x = cellsize y
74 : integer(i4) :: nodata_value !< code to define the mask
75 : real(dp), dimension(:,:), allocatable :: easting !< irregular grid easting
76 : real(dp), dimension(:,:), allocatable :: northing !< irregular grid northing
77 : real(dp), dimension(:), allocatable :: latitude !< latitude for the output
78 : real(dp), dimension(:), allocatable :: longitude !< longitude for the output
79 : end type gridGeoRef
80 : type (gridGeoRef) :: grid !< grid
81 : type (gridGeoRef) :: gridMeteo !< reference of the metereological variables
82 :
83 : !> \class period
84 : !> \brief PERIOD description
85 : type period
86 : integer(i4) :: dStart !< first day
87 : integer(i4) :: mStart !< first month
88 : integer(i4) :: yStart !< first year
89 : integer(i4) :: dEnd !< last day
90 : integer(i4) :: mEnd !< last month
91 : integer(i4) :: yEnd !< last year
92 : integer(i4) :: julStart !< first julian day
93 : integer(i4) :: julEnd !< last julian day
94 : integer(i4) :: nObs !< total number of observations
95 : CONTAINS
96 : !> \copydoc mainvar::init
97 : procedure :: init !< \see mainvar::init
98 : end type period
99 :
100 : contains
101 :
102 : !> \brief initialize a period
103 0 : subroutine init(self, dStart, mStart, yStart, dEnd, mEnd, yEnd)
104 :
105 : use mo_julian, only: julday
106 :
107 : implicit none
108 :
109 : class(period), intent(inout) :: self
110 : integer(i4), intent(in) :: dStart, mStart, yStart, dEnd, mEnd, yEnd
111 :
112 0 : self%dStart = dStart
113 0 : self%mStart = mStart
114 0 : self%yStart = yStart
115 0 : self%dEnd = dEnd
116 0 : self%mEnd = mEnd
117 0 : self%yEnd = yEnd
118 :
119 0 : self%julStart = julday(dd = dStart, mm = mStart, yy = yStart)
120 0 : self%julEnd = julday(dd = dEnd, mm = mEnd, yy = yEnd)
121 0 : self%nObs = self%julEnd - self%julStart + 1_i4
122 :
123 0 : end subroutine init
124 :
125 0 : end module mainVar
126 :
127 : !> \brief RUN Control
128 : !> \copyright Copyright 2005-\today, the CHS Developers, Sabine Attinger: All rights reserved.
129 : !! EDK is released under the LGPLv3+ license \license_note
130 : module runControl
131 : use mo_kind, only: i4, sp
132 : ! timer
133 : character(256) :: DataPathIn
134 : character(256) :: DataPathOut
135 : character(256) :: DataPathDEM
136 : character(256) :: fNameDEM !< DEM file name
137 : character(256) :: fNameSTA !< Station id and coordinates
138 : character(256) :: fNameVARIO !< file name of variogram parameters
139 : character(256) :: prefix !< prefix data file
140 : integer(i4) :: nBlocks !< number of interpolation blocks
141 : logical :: flagEDK !< estimate EDK (T/F)
142 : logical :: flagVario !< estimate variogran (T/F)
143 : logical :: correctNeg !< correct negative interpolated values
144 : logical :: distZero !< values further away than dist threshold are set to zero
145 : integer(i4) :: interMth !< interp. method
146 : end module runControl
147 :
148 : !> \brief kriging variables
149 : !> \copyright Copyright 2005-\today, the CHS Developers, Sabine Attinger: All rights reserved.
150 : !! EDK is released under the LGPLv3+ license \license_note
151 : module kriging
152 : use mo_kind, only: i4, sp, dp
153 : use mo_edk_types, only: dist_t
154 : real(dp) :: maxDist !< max distance [m] search stations
155 :
156 : !> \class cellcoarser
157 : !> \brief cell coarser type
158 : type CellCoarser
159 : integer(i4) :: nNS !< No. Nearest Stations (NS) d <= maxDist
160 : integer(i4), dimension(:), allocatable :: listNS !< list of NS
161 : real(dp) :: x !< x- coordinate
162 : real(dp) :: y !< y- coordinate
163 : real(dp) :: h !< (estimated) elevation [m] (from the nearest cells DEM)
164 : real(sp), allocatable :: z(:) !< z values to be interpolated (OUTPUT)
165 : end type CellCoarser
166 : type(CellCoarser), dimension(:), allocatable :: cell !< EDK output
167 :
168 : type(dist_t) :: edk_dist !< distance calculations for EDK
169 : real(dp) :: xl, xr, yd, yu !< coordinates of the interpolation block
170 0 : end module kriging
171 :
172 : !> \brief variogram fitting variables
173 : !> \copyright Copyright 2005-\today, the CHS Developers, Sabine Attinger: All rights reserved.
174 : !! EDK is released under the LGPLv3+ license \license_note
175 : module VarFit
176 : use mo_kind, only: i4, dp
177 : implicit none
178 : ! parameters
179 : integer(i4) :: vType !< variogram type
180 : integer(i4) :: nParam !< number of parameters
181 : integer(i4), dimension(:), allocatable :: Nh !< number of pairs per bin
182 : integer(i4) :: nbins
183 : integer(i4) :: nobs !< number of observations
184 : real(dp) :: dh !< bin size
185 : real(dp) :: hMax !< max distance h variogram
186 : real(dp), dimension(:,:), allocatable :: gamma !< variograms
187 : real(dp) :: gmax(2)
188 : real(dp) :: m0
189 : real(dp) :: v0
190 : real(dp), dimension(:), allocatable :: beta !< parameters of the variogram
191 : real(dp), dimension(8) :: E !< efficiency measures
192 : real(dp), parameter :: sRadius = 1.d4 !< searching distance limits
193 : real(dp), parameter :: gradE = 2.d-1 !< gradient limit = delta E / distance (searching)
194 :
195 : end module VarFit
196 :
197 : !> \brief NetCDF IO specifications variables
198 : !> \copyright Copyright 2005-\today, the CHS Developers, Sabine Attinger: All rights reserved.
199 : !! EDK is released under the LGPLv3+ license \license_note
200 : module NetCDFVar
201 : use mo_kind, only : dp
202 : implicit none
203 : ! Variables for NetCDF writing
204 : character(256) :: fileOut !< File Name out
205 : real(dp), dimension(:), allocatable, target :: yCoor !< GK4 (DHDN3-zone 4) easting
206 : real(dp), dimension(:), allocatable, target :: xCoor !< GK4 (DHDN3-zone 4) northing
207 : real(dp), dimension(:,:), allocatable, target :: lons !< WGS84 lons
208 : real(dp), dimension(:,:), allocatable, target :: lats !< WGS84 lats
209 : !
210 : character(256) :: variable_name !< name of netcdf variable
211 : character(256) :: variable_unit !< unit of netcdf variable
212 : character(256) :: variable_long_name !< long name of netcdf variable
213 : character(256) :: originator !< originator of netcdf file
214 : character(256) :: crs !< coordinate reference system (EPSG:XXXX)
215 : character(256) :: title !< title of netcdf file
216 : character(256) :: source !< source and methods of original data
217 : character(256) :: contact !< contact (email address)
218 : character(256) :: institution !< institution
219 : character(256) :: variable_standard_name !< standard name of netcdf variable
220 : character(256) :: variable_calendar_type !< calendar type (time variable attribute)
221 : logical :: invert_y
222 : !
223 : ! netcdf input specifications
224 : character(256) :: ncIn_variable_name !< netcdf input
225 : character(256) :: ncIn_dem_variable_name !< netcdf input
226 : character(256) :: ncIn_yCoord_name !< netcdf input
227 : character(256) :: ncIn_xCoord_name !< netcdf input
228 : character(256) :: ncOut_dem_variable_name !< netcdf input
229 : character(256) :: ncOut_dem_yCoord_name !< netcdf input
230 : character(256) :: ncOut_dem_xCoord_name !< netcdf input
231 : character(256) :: ncOut_dem_Latitude !< netcdf input
232 : character(256) :: ncOut_dem_Longitude !< netcdf input
233 :
234 : end module NetCDFVar
|