Skip to content
Permalink
Branch: master
Find file Copy path
Fetching contributors…
Cannot retrieve contributors at this time
5294 lines (3991 sloc) 240 KB
!-------------------------------------------------------------------------
! IST/MARETEC, Marine Modelling Group, Mohid2000 modelling system
!-------------------------------------------------------------------------
!BOI
! !TITLE: Mohid2000 Statistics module
! !AUTHORS: Paulo Chambel, Frank Braunschweig
! !AFFILIATION: IST/MARETEC, Marine Modelling Group
! !DATE: May2002
! !INTRODUCTION: Module which makes basic statistic operations over matrixes
!
!EOI
!------------------------------------------------------------------------------
!
!This program is free software; you can redistribute it and/or
!modify it under the terms of the GNU General Public License
!version 2, as published by the Free Software Foundation.
!
!This program is distributed in the hope that it will be useful,
!but WITHOUT ANY WARRANTY; without even the implied warranty of
!MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
!GNU General Public License for more details.
!
!You should have received a copy of the GNU General Public License
!along with this program; if not, write to the Free Software
!Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
!
!------------------------------------------------------------------------------
Module ModuleStatistic
use ModuleGlobalData
use ModuleTime
use ModuleStopWatch, only: StartWatch, StopWatch
use ModuleEnterData
use ModuleFunctions
use ModuleHDF5
implicit none
private
!Subroutines---------------------------------------------------------------
!Constructor
public :: ConstructStatistic
private :: AllocateInstance
private :: ReadDataFile
private :: AllocateStatisticMatrixes
!Modifier
public :: ModifyStatistic
private :: StatisticAnalysis3D
private :: ModifyGlobalStatistic
private :: ModifyDailyStatistic
private :: ModifyMonthlyStatistic
private :: ModifySpecificHourStatistic
private :: ModifyClassification
private :: StatisticAnalysis2D
private :: ModifyGlobalStatistic2D
private :: ModifyDailyStatistic2D
private :: ModifyMonthlyStatistic2D
private :: ModifySpecificHourStatistic2D
private :: ModifyClassification2D
private :: WriteValuesToFileHDF5
public :: AddStatisticLayers
private :: AverageValueBetweenDepths
private :: AverageValueBetweenLayers
!Selector
public :: GetStatisticMethod
public :: GetStatisticLayerDef
public :: GetStatisticParameters
public :: GetStatisticLayersNumber
public :: GetStatisticClasses
public :: GetStatisticClassesNumber
public :: UnGetStatistic
!Destructor
public :: KillStatistic
private :: DeallocateInstance
private :: DeallocateStatisticMatrixes
!Interfaces----------------------------------------------------------------
private :: UngetStatistic2D
interface UngetStatistic
module procedure UngetStatistic2D
end interface UngetStatistic
!Parameter-----------------------------------------------------------------
integer, parameter :: Value3DStat3D_ = 1, Value3DStatLayers_ = 2, Value2DStat2D_ = 3
integer, parameter :: Depth_ = 1, Layer_ = 2
type T_Classification
logical :: On = .false.
real :: Percentil = null_real
integer :: nClasses = null_int
real, dimension(:, :), pointer :: Classes => null()
real, dimension(:, :, :, :), pointer :: Frequency => null()
real, dimension(:, : , :), pointer :: Frequency2D => null()
real :: RunPeriod = null_real
type (T_Time) :: LastCalculation
end type T_Classification
type T_SimpleStatistic
logical :: On = .false.
real, dimension(:, :, :), pointer :: Minimum => null() !inicialization: Carina
real, dimension(:, :, :), pointer :: Maximum => null() !inicialization: Carina
real, dimension(:, :, :), pointer :: Average => null() !inicialization: Carina
real, dimension(:, :, :), pointer :: SquareAverage => null() !inicialization: Carina
real, dimension(:, :, :), pointer :: StandardDeviation => null() !inicialization: Carina
real, dimension(:, :, :), pointer :: GeomAverage => null() !inicialization: Carina
real, dimension(:, :, :), pointer :: SquareGeomAverage => null() !inicialization: Carina
real, dimension(:, :, :), pointer :: GeomStandardDeviation => null() !inicialization: Carina
real, dimension(:, :, :), pointer :: Accumulated => null() !inicialization: Carina
!guillaume juan
real, dimension(:, :, :), pointer :: PercentBelowCriticalValue => null() !inicialization: Carina
real, dimension(:, : ), pointer :: Minimum2D => null() !inicialization: Carina
real, dimension(:, : ), pointer :: Maximum2D => null() !inicialization: Carina
real, dimension(:, : ), pointer :: Average2D => null() !inicialization: Carina
real, dimension(:, : ), pointer :: SquareAverage2D => null() !inicialization: Carina
real, dimension(:, : ), pointer :: StandardDeviation2D => null() !inicialization: Carina
real, dimension(:, : ), pointer :: GeomAverage2D => null() !inicialization: Carina
real, dimension(:, : ), pointer :: SquareGeomAverage2D => null() !inicialization: Carina
real, dimension(:, : ), pointer :: GeomStandardDeviation2D => null() !inicialization: Carina
real, dimension(:, : ), pointer :: Accumulated2D => null() !inicialization: Carina
!guillaume juan
real, dimension(:, : ), pointer :: PercentBelowCriticalValue2D => null() !inicialization: Carina
real :: RunPeriod = null_real !inicialization: Carina
integer :: OutputNumber = null_int !inicialization: Carina
type (T_Time) :: LastCalculation
type (T_Time) :: NextOutputTime
end type T_SimpleStatistic
type T_Layers
integer :: Number = null_int !inicialization: Carina
integer, dimension(:), pointer :: Definition => null() !inicialization: Carina !1 - Matrix of layer numbers
!2 - Matrix of depths
! These two matrixes are depths and are allocated only when Methodology = 2 (Values 3D and Statistic 2D)
real, dimension(:,:,: ), pointer :: UpperDepth => null() !inicialization: Carina
real, dimension(:,:,: ), pointer :: LowerDepth => null() !inicialization: Carina
!These two matrixes are layers number and are allocated only when Methodology = 2 (Values 3D and Statistic 2D)
real, dimension(:,:,: ), pointer :: UpperLayer => null() !inicialization: Carina
real, dimension(:,:,: ), pointer :: LowerLayer => null() !inicialization: Carina
!The average value in the layer
real, dimension(:,:,: ), pointer :: Value => null() !inicialization: Carina
!if exist one water point between the layer limits than WaterPoints2D (i, j) = 1
integer, dimension(:,:,: ), pointer :: WaterPoints => null() !inicialization: Carina
integer, dimension(:), pointer :: MinLayer => null() !inicialization: Carina
integer, dimension(:), pointer :: MaxLayer => null() !inicialization: Carina
real, dimension(:), pointer :: MinDepth => null() !inicialization: Carina
real, dimension(:), pointer :: MaxDepth => null() !inicialization: Carina
end type T_Layers
type T_ExternalVar
type (T_Time ) :: Now
type (T_Time) :: BeginTime, EndTime
type (T_Size3D ) :: Size
type (T_Size3D ) :: WorkSize
real, dimension(:, :, :), pointer :: Value3D => null() !inicialization: Carina
real, dimension(:, : ), pointer :: Value2D => null() !inicialization: Carina
end type T_ExternalVar
type T_Statistic
integer :: InstanceID = null_int !inicialization: Carina
type (T_ExternalVar ) :: ExternalVar
integer :: Methodology = null_int !inicialization: Carina
!1 - Values 3D and Statistic 3D
!2 - Values 3D and Statistic 2D (2D layer)
!3 - Values 2D and Statistic 2D
type (T_Layers ) :: Layers
type (T_SimpleStatistic) :: Global
type (T_SimpleStatistic) :: Daily
type (T_SimpleStatistic) :: Monthly
type (T_SimpleStatistic) :: SpecificHour
integer :: SpecificHourValue = null_int !inicialization: Carina
type (T_Classification ) :: Classification
integer :: ObjTime = 0
integer :: ObjHDF5 = 0
integer :: ObjEnterData = 0
character(len=StringLength) :: Name = null_str !inicialization: Carina
character(len=StringLength) :: GroupName = null_str !inicialization: Carina
logical :: Accumulated = .false. !Accumulated Values
logical :: GeomMean = .false. !Geometric statistics
!guillaume juan
logical :: Critical = .false.
real :: CriticalValue = null_real !inicialization: Carina
logical :: NormalizeFreq = .false. !inicialization: Carina
type (T_Statistic ), pointer :: Next => null() !inicialization: Carina
end type T_Statistic
type (T_Statistic), pointer :: FirstStatistic => null()
type (T_Statistic), pointer :: Me => null()
contains
!++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
!++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
!CONSTRUCTOR CONSTRUCTOR CONSTRUCTOR CONSTRUCTOR CONSTRUCTOR CONSTRUCTOR
!++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
!++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
!--------------------------------------------------------------------------
subroutine ConstructStatistic (StatisticID, ObjTime, ObjHDF5, &
Size, WorkSize, DataFile, Name, GroupName, Rank, &
STAT)
!Arguments-------------------------------------------------------------
integer :: StatisticID
integer :: ObjTime
integer :: ObjHDF5
type (T_Size3D) :: Size, WorkSize
character(len=*) :: DataFile
character(len=*) :: Name
character(len=*), optional :: GroupName
integer, optional :: Rank
integer, optional :: STAT
!Local-----------------------------------------------------------------
integer :: STAT_CALL
integer :: STAT_, ready_
type (T_Time) :: AuxTime
real :: Year,Month,Day
!----------------------------------------------------------------------
STAT_ = UNKNOWN_
!Assures nullification of the global variable
if (.not. ModuleIsRegistered(mStatistic_)) then
nullify (FirstStatistic)
call RegisterModule (mStatistic_)
endif
write(*,*) 'constructing...'
call Ready (StatisticID, ready_)
write(*,*) 'constructed...'
if (ready_ .EQ. OFF_ERR_) then
call AllocateInstance
!Name
Me%Name = Name
if (present(GroupName)) then
Me%GroupName = "/Statistics/"//GroupName
else
Me%GroupName = "/Statistics/"
endif
!Associates External Instances
Me%ObjTime = AssociateInstance (mTIME_, ObjTime )
Me%ObjHDF5 = AssociateInstance (mHDF5_, ObjHDF5 )
! Actualized the time
call GetComputeCurrentTime(Me%ObjTime, &
Me%ExternalVar%Now, STAT = STAT_CALL)
if (STAT_CALL /= SUCCESS_) stop 'ConstructStatistic - ModuleStatistic - ERR03'
call GetComputeTimeLimits(TimeID = Me%ObjTime, &
EndTime = Me%ExternalVar%EndTime, &
BeginTime = Me%ExternalVar%BeginTime, &
STAT = STAT_CALL)
if (STAT_CALL /= SUCCESS_) stop 'ConstructStatistic - ModuleStatistic - ERR05'
!Sets Size
Me%ExternalVar%Size = Size
Me%ExternalVar%WorkSize = WorkSize
!Reads Data File
call ReadDataFile (DataFile)
if (present(Rank)) then
if (Rank == 2) then
Me%Methodology = Value2DStat2D_
elseif (Rank == 3) then
if (Me%Methodology /= Value3DStatLayers_) then
Me%Methodology = Value3DStat3D_
endif
endif
endif
!João Sobrinho
if (Me%Methodology == Value3DStatLayers_) then
call AllocateLayerMatrixes
endif
!Allocates Variables 3D
if (Me%Global%On) &
call AllocateStatisticMatrixes (Me%Global)
if (Me%Daily%On) then
call AllocateStatisticMatrixes (Me%Daily)
call ExtractDate (Time1=Me%ExternalVar%Now, Year=Year, Month=Month, Day=Day)
call SetDate(Time1=AuxTime, Year=Year, Month=Month, Day=Day, &
Hour=0.0, Minute=0.0, Second=0.0)
Me%Daily%NextOutputTime = AuxTime + 24.*3600.
endif
if (Me%Monthly%On) then
call AllocateStatisticMatrixes (Me%Monthly)
endif
if (Me%SpecificHour%On) &
call AllocateStatisticMatrixes (Me%SpecificHour)
if (Me%Classification%On) &
Call AllocateFrequencyMatrixes
!Returns Statistic InstanceID
StatisticID = Me%InstanceID
else
write(*,*) Name
stop 'ModuleStatistic - ConstructStatistic - ERR99'
endif
STAT_ = SUCCESS_
if (present(STAT)) STAT = STAT_
end subroutine ConstructStatistic
!--------------------------------------------------------------------------
subroutine AllocateInstance
!Arguments-------------------------------------------------------------
!Local-----------------------------------------------------------------
type (T_Statistic), pointer :: NewStatistic
type (T_Statistic), pointer :: PreviousStatistic
!Allocates new instance
allocate (NewStatistic)
nullify (NewStatistic%Next)
!Insert New Instance into list and makes Current point to it
if (.not. associated(FirstStatistic)) then
FirstStatistic => NewStatistic
Me => NewStatistic
else
PreviousStatistic => FirstStatistic
Me => FirstStatistic%Next
do while (associated(Me))
PreviousStatistic => Me
Me => Me%Next
enddo
Me => NewStatistic
PreviousStatistic%Next => NewStatistic
endif
Me%InstanceID = RegisterNewInstance (mSTATISTIC_)
end subroutine AllocateInstance
!--------------------------------------------------------------------------
subroutine ReadDataFile (DataFile)
!Arguments-------------------------------------------------------------
character (len=*) :: DataFile
!Local-----------------------------------------------------------------
integer :: FromFile, iflag
integer :: STAT_CALL
integer :: ClientNumber
logical :: BlockFound
integer :: FirstLine, LastLine, iLine
integer :: nClasses
integer :: ILB, IUB
integer :: JLB, JUB
integer :: KLB, KUB
!Shorten
ILB = Me%ExternalVar%Size%ILB
IUB = Me%ExternalVar%Size%IUB
JLB = Me%ExternalVar%Size%JLB
JUB = Me%ExternalVar%Size%JUB
KLB = Me%ExternalVar%Size%KLB
KUB = Me%ExternalVar%Size%KUB
call GetExtractType (FromFile = FromFile)
!Constructs the data file
call ConstructEnterData (Me%ObjEnterData, DataFile, STAT = STAT_CALL)
if (STAT_CALL /= SUCCESS_) &
call SetError (FATAL_, KEYWORD_, 'ModuleStatistic - ReadDataFile - ERR10')
!1 - Values 3D and Statistic 3D
!2 - Values 3D and Statistic 2D (2D layer)
!3 - Values 2D and Statistic 2D
call GetData(Me%Methodology, &
Me%ObjEnterData, &
iflag, &
SearchType = FromFile, &
keyword = 'METHOD_STATISTIC', &
default = Value3DStat3D_, &
ClientModule = 'ModuleStatistic', &
STAT = STAT_CALL)
if (STAT_CALL /= SUCCESS_) &
call SetError (FATAL_, KEYWORD_, 'ModuleStatistic - ReadDataFile - ERR20')
!Global Statistic
call GetData(Me%Global%On, &
Me%ObjEnterData, &
iflag, &
SearchType = FromFile, &
keyword = 'GLOBAL_STATISTIC', &
default = .false., &
ClientModule = 'ModuleStatistic', &
STAT = STAT_CALL)
if (STAT_CALL /= SUCCESS_) &
call SetError (FATAL_, KEYWORD_, 'ModuleStatistic - ReadDataFile - ERR30')
!Daily Statistic
call GetData(Me%Daily%On, &
Me%ObjEnterData, &
iflag, &
SearchType = FromFile, &
keyword = 'DAILY_STATISTIC', &
default = .false., &
ClientModule = 'ModuleStatistic', &
STAT = STAT_CALL)
if (STAT_CALL /= SUCCESS_) &
call SetError (FATAL_, KEYWORD_, 'ModuleStatistic - ReadDataFile - ERR40')
!Monthly Statistic
call GetData(Me%Monthly%On, &
Me%ObjEnterData, &
iflag, &
SearchType = FromFile, &
keyword = 'MONTHLY_STATISTIC', &
default = .false., &
ClientModule = 'ModuleStatistic', &
STAT = STAT_CALL)
if (STAT_CALL /= SUCCESS_) &
call SetError (FATAL_, KEYWORD_, 'ModuleStatistic - ReadDataFile - ERR50')
!Statistics for a specific hour of day
call GetData(Me%SpecificHour%On, &
Me%ObjEnterData, &
iflag, &
SearchType = FromFile, &
keyword = 'SPECIFIC_HOUR_STATISTIC', &
default = .false., &
ClientModule = 'ModuleStatistic', &
STAT = STAT_CALL)
if (STAT_CALL /= SUCCESS_) &
call SetError (FATAL_, KEYWORD_, 'ModuleStatistic - ReadDataFile - ERR60')
if (Me%SpecificHour%On) then
!Statistics for a specific hour of day
call GetData(Me%SpecificHourValue, &
Me%ObjEnterData, &
iflag, &
SearchType = FromFile, &
keyword = 'SPECIFIC_HOUR', &
default = 12, &
ClientModule = 'ModuleStatistic', &
STAT = STAT_CALL)
if (STAT_CALL /= SUCCESS_) &
call SetError (FATAL_, KEYWORD_, 'ModuleStatistic - ReadDataFile - ERR70')
endif
!Accumulated values
call GetData(Me%Accumulated, &
Me%ObjEnterData, &
iflag, &
SearchType = FromFile, &
keyword = 'ACCUMULATED', &
default = .false., &
ClientModule = 'ModuleStatistic', &
STAT = STAT_CALL)
if (STAT_CALL /= SUCCESS_) &
call SetError (FATAL_, KEYWORD_, 'ModuleStatistic - ReadDataFile - ERR80')
!Geometric Average and Standard Deviation
call GetData(Me%GeomMean, &
Me%ObjEnterData, &
iflag, &
SearchType = FromFile, &
keyword = 'GEOMETRIC_MEAN', &
default = .false., &
ClientModule = 'ModuleStatistic', &
STAT = STAT_CALL)
if (STAT_CALL /= SUCCESS_) &
call SetError (FATAL_, KEYWORD_, 'ModuleStatistic - ReadDataFile - ERR90')
!guillaume juan
call GetData(Me%Critical, &
Me%ObjEnterData, &
iflag, &
SearchType = FromFile, &
keyword = 'CRITICAL', &
default = .false., &
ClientModule = 'ModuleStatistic', &
STAT = STAT_CALL)
if (STAT_CALL /= SUCCESS_) &
call SetError (FATAL_, KEYWORD_, 'ModuleStatistic - ReadDataFile - ERR100')
!guillaume juan
if (Me%Critical) then
!Statistics for a specific hour of day
call GetData(Me%CriticalValue, &
Me%ObjEnterData, &
iflag, &
SearchType = FromFile, &
keyword = 'CRITICAL_VALUE', &
default = 0.02, &
ClientModule = 'ModuleStatistic', &
STAT = STAT_CALL)
if (STAT_CALL /= SUCCESS_) &
call SetError (FATAL_, KEYWORD_, 'ModuleStatistic - ReadDataFile - ERR110')
endif
call GetData(Me%NormalizeFreq, &
Me%ObjEnterData, &
iflag, &
SearchType = FromFile, &
keyword = 'NORMALIZE_FREQ', &
default = .false., &
ClientModule = 'ModuleStatistic', &
STAT = STAT_CALL)
if (STAT_CALL /= SUCCESS_) &
call SetError (FATAL_, KEYWORD_, 'ModuleStatistic - ReadDataFile - ERR120')
Me%Classification%On = .false.
nullify (Me%Classification%Classes)
nullify (Me%Classification%Frequency)
nullify (Me%Classification%Frequency2D)
!Classification
!The infinit loop is necessary to scan the file again. When
!the block ("<BeginClass>", "<EndClass>") is not found (BlockFound = .false.)
!all the variables are initialize and the Object EnterData is read to be scan again
!(ex: scan the block "<beginlayer>, <endlayer>).
do1 : do
call ExtractBlockFromBuffer(Me%ObjEnterData, ClientNumber, &
"<BeginClass>", "<EndClass>", &
BlockFound, &
FirstLine = FirstLine, &
LastLine = LastLine, &
STAT = STAT_CALL)
if (STAT_CALL /= SUCCESS_) &
call SetError (FATAL_, KEYWORD_, 'ModuleStatistic - ReadDataFile - ERR07')
!Reads Classification Data
cd1: if (BlockFound) then
Me%Classification%On = .true.
nClasses = (LastLine-1) - (FirstLine+1) + 1
Me%Classification%nClasses = nClasses
allocate (Me%Classification%Classes(nClasses, 2))
nClasses = 0
do iLine = FirstLine+1, LastLine-1
nClasses = nClasses + 1
call GetData(Me%Classification%Classes(nClasses, :), &
Me%ObjEnterData, iflag, &
Buffer_Line = iLine, &
STAT = STAT_CALL)
if (iflag /= 2) stop 'ModuleStatistic - ReadDataFile - ERR08'
if (STAT_CALL /= SUCCESS_) stop 'ModuleStatistic - ReadDataFile - ERR09'
enddo
else cd1
call Block_Unlock(Me%ObjEnterData, ClientNumber, STAT = STAT_CALL)
if (STAT_CALL /= SUCCESS_) &
call SetError (FATAL_, KEYWORD_, 'ModuleStatistic - ReadDataFile - ERR10')
exit do1 !No more blocks
endif cd1
enddo do1
if (Me%Classification%On) then
!Percentil
call GetData(Me%Classification%Percentil, &
Me%ObjEnterData, &
iflag, &
SearchType = FromFile, &
keyword = 'PERCENTILE', &
default = 90., &
ClientModule = 'ModuleStatistic', &
STAT = STAT_CALL)
if (STAT_CALL /= SUCCESS_) &
call SetError (FATAL_, KEYWORD_, 'ModuleStatistic - ReadDataFile - ERR11')
endif
if (Me%Methodology == Value3DStatLayers_) &
call Construct_Layers
call KillEnterData (Me%ObjEnterData, STAT = STAT_CALL)
if (STAT_CALL /= SUCCESS_) &
call SetError (FATAL_, KEYWORD_, 'ModuleStatistic - ReadDataFile - ERR12')
end subroutine ReadDataFile
!--------------------------------------------------------------------------
subroutine AllocateLayerMatrixes
!Arguments-------------------------------------------------------------
!Local-----------------------------------------------------------------
integer :: ILB, IUB, i
integer :: JLB, JUB
!Shorten
ILB = Me%ExternalVar%Size%ILB
IUB = Me%ExternalVar%Size%IUB
JLB = Me%ExternalVar%Size%JLB
JUB = Me%ExternalVar%Size%JUB
allocate (Me%Layers%Value (ILB:IUB, JLB:JUB, 1:Me%Layers%Number))
Me%Layers%Value(:,:,:) = FillValueReal
allocate (Me%Layers%WaterPoints(ILB:IUB, JLB:JUB, 1:Me%Layers%Number))
Me%Layers%WaterPoints (:,:,:) = FillValueInt
allocate (Me%Layers%UpperDepth(ILB:IUB, JLB:JUB, 1:Me%Layers%Number))
allocate (Me%Layers%LowerDepth(ILB:IUB, JLB:JUB, 1:Me%Layers%Number))
allocate (Me%Layers%UpperLayer(ILB:IUB, JLB:JUB, 1:Me%Layers%Number))
allocate (Me%Layers%LowerLayer(ILB:IUB, JLB:JUB, 1:Me%Layers%Number))
do i=1, Me%Layers%Number
Me%Layers%UpperDepth(:,:,i) = Me%Layers%MinDepth(i)
Me%Layers%LowerDepth(:,:,i) = Me%Layers%MaxDepth(i)
Me%Layers%UpperLayer(:,:,i) = Me%Layers%MaxLayer(i)
Me%Layers%LowerLayer(:,:,i) = Me%Layers%MinLayer(i)
enddo
end subroutine AllocateLayerMatrixes
!--------------------------------------------------------------------------
subroutine AllocateStatisticMatrixes (Statistic)
!Arguments-------------------------------------------------------------
type (T_SimpleStatistic) :: Statistic
!Local-----------------------------------------------------------------
integer :: ILB, IUB
integer :: JLB, JUB
integer :: KLB, KUB
!Shorten
ILB = Me%ExternalVar%Size%ILB
IUB = Me%ExternalVar%Size%IUB
JLB = Me%ExternalVar%Size%JLB
JUB = Me%ExternalVar%Size%JUB
KLB = Me%ExternalVar%Size%KLB
KUB = Me%ExternalVar%Size%KUB
Statistic%RunPeriod = 0.
Statistic%OutputNumber = 1.
Statistic%LastCalculation = Me%ExternalVar%Now
if (Me%Methodology == Value3DStat3D_) then
allocate (Statistic%Minimum (ILB:IUB, JLB:JUB, KLB:KUB))
allocate (Statistic%Maximum (ILB:IUB, JLB:JUB, KLB:KUB))
allocate (Statistic%Average (ILB:IUB, JLB:JUB, KLB:KUB))
allocate (Statistic%SquareAverage (ILB:IUB, JLB:JUB, KLB:KUB))
allocate (Statistic%StandardDeviation(ILB:IUB, JLB:JUB, KLB:KUB))
allocate (Statistic%Minimum2D (ILB:IUB, JLB:JUB))
allocate (Statistic%Maximum2D (ILB:IUB, JLB:JUB))
Statistic%Minimum = - FillValueReal
Statistic%Maximum = + FillValueReal
Statistic%Average = + FillValueReal
Statistic%SquareAverage = + FillValueReal
Statistic%StandardDeviation = + FillValueReal
Statistic%Minimum2D = - FillValueReal
Statistic%Maximum2D = + FillValueReal
if (Me%Accumulated) then
allocate (Statistic%Accumulated (ILB:IUB, JLB:JUB, KLB:KUB))
Statistic%Accumulated = 0.0
endif
if (Me%GeomMean) then !Geometric Average is to be calculated
allocate (Statistic%GeomAverage (ILB:IUB, JLB:JUB, KLB:KUB))
allocate (Statistic%SquareGeomAverage (ILB:IUB, JLB:JUB, KLB:KUB))
allocate (Statistic%GeomStandardDeviation(ILB:IUB, JLB:JUB, KLB:KUB))
Statistic%GeomAverage = - FillValueReal !must be positive because log
Statistic%SquareGeomAverage = + FillValueReal
Statistic%GeomStandardDeviation = + FillValueReal
endif
!guillaume juan
if (Me%Critical) then
allocate (Statistic%PercentBelowCriticalValue (ILB:IUB, JLB:JUB, KLB:KUB))
Statistic%PercentBelowCriticalValue = + FillValueReal
endif
else if (Me%Methodology == Value3DStatLayers_) then
allocate (Statistic%Minimum (ILB:IUB, JLB:JUB, 1:Me%Layers%Number))
allocate (Statistic%Maximum (ILB:IUB, JLB:JUB, 1:Me%Layers%Number))
allocate (Statistic%Average (ILB:IUB, JLB:JUB, 1:Me%Layers%Number))
allocate (Statistic%SquareAverage (ILB:IUB, JLB:JUB, 1:Me%Layers%Number))
allocate (Statistic%StandardDeviation(ILB:IUB, JLB:JUB, 1:Me%Layers%Number))
allocate (Statistic%Minimum2D (ILB:IUB, JLB:JUB))
allocate (Statistic%Maximum2D (ILB:IUB, JLB:JUB))
Statistic%Minimum = - FillValueReal
Statistic%Maximum = + FillValueReal
Statistic%Average = + FillValueReal
Statistic%SquareAverage = + FillValueReal
Statistic%StandardDeviation = + FillValueReal
Statistic%Minimum2D = - FillValueReal
Statistic%Maximum2D = + FillValueReal
if (Me%Accumulated) then
allocate (Statistic%Accumulated (ILB:IUB, JLB:JUB, 1:Me%Layers%Number))
Statistic%Accumulated = 0.0
endif
if (Me%GeomMean) then !Geometric Average is to be calculated
allocate (Statistic%GeomAverage (ILB:IUB, JLB:JUB, 1:Me%Layers%Number))
allocate (Statistic%SquareGeomAverage (ILB:IUB, JLB:JUB, 1:Me%Layers%Number))
allocate (Statistic%GeomStandardDeviation(ILB:IUB, JLB:JUB, 1:Me%Layers%Number))
Statistic%GeomAverage = - FillValueReal !must be positive because log
Statistic%SquareGeomAverage = + FillValueReal
Statistic%GeomStandardDeviation = + FillValueReal
endif
!guillaume juan
if (Me%Critical) then
allocate (Statistic%PercentBelowCriticalValue (ILB:IUB, JLB:JUB, 1:Me%Layers%Number))
Statistic%PercentBelowCriticalValue = + FillValueReal
endif
else if (Me%Methodology == Value2DStat2D_) then
allocate (Statistic%Minimum2D (ILB:IUB, JLB:JUB))
allocate (Statistic%Maximum2D (ILB:IUB, JLB:JUB))
allocate (Statistic%Average2D (ILB:IUB, JLB:JUB))
allocate (Statistic%SquareAverage2D (ILB:IUB, JLB:JUB))
allocate (Statistic%StandardDeviation2D(ILB:IUB, JLB:JUB))
Statistic%Minimum2D = - FillValueReal
Statistic%Maximum2D = + FillValueReal
Statistic%Average2D = + FillValueReal
Statistic%SquareAverage2D = + FillValueReal
Statistic%StandardDeviation2D = + FillValueReal
if (Me%Accumulated) then
allocate (Statistic%Accumulated2D (ILB:IUB, JLB:JUB))
Statistic%Accumulated2D = 0.0
endif
if (Me%GeomMean) then !Geometric Average is to be calculated
allocate (Statistic%GeomAverage2D (ILB:IUB, JLB:JUB))
allocate (Statistic%SquareGeomAverage2D (ILB:IUB, JLB:JUB))
allocate (Statistic%GeomStandardDeviation2D(ILB:IUB, JLB:JUB))
Statistic%GeomAverage2D = - FillValueReal !must be positive because log
Statistic%SquareGeomAverage2D = + FillValueReal
Statistic%GeomStandardDeviation2D = + FillValueReal
endif
!guillaume juan
if (Me%Critical) then
allocate (Statistic%PercentBelowCriticalValue2D (ILB:IUB, JLB:JUB))
Statistic%PercentBelowCriticalValue2D = + FillValueReal
endif
endif
end subroutine AllocateStatisticMatrixes
!--------------------------------------------------------------------------
subroutine AllocateFrequencyMatrixes
!Arguments-------------------------------------------------------------
!Local-----------------------------------------------------------------
integer :: ILB, IUB
integer :: JLB, JUB
integer :: KLB, KUB, nClasses
!Shorten
ILB = Me%ExternalVar%Size%ILB
IUB = Me%ExternalVar%Size%IUB
JLB = Me%ExternalVar%Size%JLB
JUB = Me%ExternalVar%Size%JUB
KLB = Me%ExternalVar%Size%KLB
KUB = Me%ExternalVar%Size%KUB
nClasses = Me%Classification%nClasses
Me%Classification%LastCalculation = Me%ExternalVar%Now
Me%Classification%RunPeriod = 0.
if (Me%Methodology == Value3DStat3D_) then
allocate (Me%Classification%Frequency (ILB:IUB, JLB:JUB, KLB:KUB, 1:nClasses))
Me%Classification%Frequency(:,:,:,:) = 0.
else if (Me%Methodology == Value3DStatLayers_) then
allocate (Me%Classification%Frequency (ILB:IUB, JLB:JUB, 1:Me%Layers%Number, 1:nClasses))
Me%Classification%Frequency(:,:,:,:) = 0.
else if (Me%Methodology == Value2DStat2D_) then
allocate (Me%Classification%Frequency2D(ILB:IUB, JLB:JUB, 1:nClasses))
Me%Classification%Frequency2D(:,:,:) = 0.
endif
end subroutine AllocateFrequencyMatrixes
!--------------------------------------------------------------------------
subroutine Construct_Layers
!Arguments---------------------------------------------------------------
!Local-------------------------------------------------------------------
character(LEN = StringLength), parameter :: block_begin = '<beginlayer>'
character(LEN = StringLength), parameter :: block_end = '<endlayer>'
real, dimension(:), pointer :: Aux1, Aux2, Aux3, Aux4, Aux5
integer :: STAT_CALL
integer :: ClientNumber
logical :: BlockFound
integer :: LayerNumber
!Begin-------------------------------------------------------------------
LayerNumber = 0
do1 : do
call ExtractBlockFromBuffer(Me%ObjEnterData, &
ClientNumber = ClientNumber, &
block_begin = block_begin, &
block_end = block_end, &
BlockFound = BlockFound, &
STAT = STAT_CALL)
cd1 : if (STAT_CALL .EQ. SUCCESS_ ) then
cd2 : if (BlockFound) then
LayerNumber = LayerNumber + 1
if (associated (Me%Layers%MinDepth)) &
deallocate (Me%Layers%MinDepth)
if (associated (Me%Layers%MaxDepth)) &
deallocate (Me%Layers%MaxDepth)
if (associated (Me%Layers%MaxLayer)) &
deallocate (Me%Layers%MaxLayer)
if (associated (Me%Layers%MinLayer)) &
deallocate(Me%Layers%MinLayer)
if (associated (Me%Layers%Definition)) &
deallocate(Me%Layers%Definition)
allocate(Me%Layers%MinDepth (LayerNumber))
allocate(Me%Layers%MaxDepth (LayerNumber))
allocate(Me%Layers%MaxLayer (LayerNumber))
allocate(Me%Layers%MinLayer (LayerNumber))
allocate(Me%Layers%Definition(LayerNumber))
if (LayerNumber>1) then
Me%Layers%MinDepth (1:LayerNumber-1) = Aux1(1:LayerNumber-1)
Me%Layers%MaxDepth (1:LayerNumber-1) = Aux2(1:LayerNumber-1)
Me%Layers%MinLayer (1:LayerNumber-1) = Aux3(1:LayerNumber-1)
Me%Layers%MaxLayer (1:LayerNumber-1) = Aux4(1:LayerNumber-1)
Me%Layers%Definition(1:LayerNumber-1) = Aux5(1:LayerNumber-1)
deallocate(Aux1)
deallocate(Aux2)
deallocate(Aux3)
deallocate(Aux4)
deallocate(Aux5)
endif
! Read the characteristics of a New Layer
Call Read_Layer(LayerNumber)
allocate(Aux1(LayerNumber))
allocate(Aux2(LayerNumber))
allocate(Aux3(LayerNumber))
allocate(Aux4(LayerNumber))
allocate(Aux5(LayerNumber))
Aux1(:) = Me%Layers%MinDepth (:)
Aux2(:) = Me%Layers%MaxDepth (:)
Aux3(:) = Me%Layers%MinLayer (:)
Aux4(:) = Me%Layers%MaxLayer (:)
Aux5(:) = Me%Layers%Definition(:)
else cd2
if (associated (Aux1)) deallocate(Aux1)
if (associated (Aux2)) deallocate(Aux2)
if (associated (Aux3)) deallocate(Aux3)
if (associated (Aux4)) deallocate(Aux4)
if (associated (Aux5)) deallocate(Aux5)
Me%Layers%Number = LayerNumber
call Block_Unlock(Me%ObjEnterData, ClientNumber, STAT = STAT_CALL)
if (STAT_CALL .NE. SUCCESS_) &
stop 'Subroutine Construct_Layers; ModuleStatistic. ERR01.'
exit do1 !No more blocks
end if cd2
else if (STAT_CALL .EQ. BLOCK_END_ERR_) then cd1
write(*,*)
write(*,*) 'Error calling ExtractBlockFromBuffer. '
stop 'Construct_Layers; ModuleStatistic. ERR02'
else cd1
stop 'Construct_Layers; ModuleStatistic. ERR03'
end if cd1
end do do1
!------------------------------------------------------------------------
end subroutine Construct_Layers
!--------------------------------------------------------------------------
subroutine Read_Layer(LayerNumber)
!Arguments---------------------------------------------------------------
integer :: LayerNumber
!Local----------------------------------------------------------------
integer :: FromBlock, iflag
integer :: STAT_CALL
!Begin-------------------------------------------------------------------
call GetExtractType (FromBlock = FromBlock)
!2 - Matrix of layer numbers
!1 - depth
call GetData(Me%Layers%Definition(LayerNumber), &
Me%ObjEnterData, &
iflag, &
SearchType = FromBlock, &
keyword = 'LAYER_DEFINITION', &
default = Depth_, &
ClientModule = 'ModuleStatistic', &
STAT = STAT_CALL)
if (STAT_CALL /= SUCCESS_) &
call SetError (FATAL_, KEYWORD_, 'ModuleStatistic - Read_Layer - ERR01')
call GetData(Me%Layers%MinLayer(LayerNumber), &
Me%ObjEnterData, &
iflag, &
SearchType = FromBlock, &
keyword = 'MIN_LAYER', &
default = Me%ExternalVar%WorkSize%KLB, &
ClientModule = 'ModuleStatistic', &
STAT = STAT_CALL)
if (STAT_CALL /= SUCCESS_) &
call SetError (FATAL_, KEYWORD_, 'ModuleStatistic - Read_Layer - ERR02')
call GetData(Me%Layers%MaxLayer(LayerNumber), &
Me%ObjEnterData, &
iflag, &
SearchType = FromBlock, &
keyword = 'MAX_LAYER', &
default = Me%ExternalVar%WorkSize%KUB, &
ClientModule = 'ModuleStatistic', &
STAT = STAT_CALL)
if (STAT_CALL /= SUCCESS_) &
call SetError (FATAL_, KEYWORD_, 'ModuleStatistic - Read_Layer - ERR03')
call GetData(Me%Layers%MinDepth(LayerNumber), &
Me%ObjEnterData, &
iflag, &
SearchType = FromBlock, &
keyword = 'MIN_DEPTH', &
default = 0., &
ClientModule = 'ModuleStatistic', &
STAT = STAT_CALL)
if (STAT_CALL /= SUCCESS_) &
call SetError (FATAL_, KEYWORD_, 'ModuleStatistic - Read_Layer - ERR04')
call GetData(Me%Layers%MaxDepth(LayerNumber), &
Me%ObjEnterData, &
iflag, &
SearchType = FromBlock, &
keyword = 'MAX_DEPTH', &
default = 10000., &
ClientModule = 'ModuleStatistic', &
STAT = STAT_CALL)
if (STAT_CALL /= SUCCESS_) &
call SetError (FATAL_, KEYWORD_, 'ModuleStatistic - Read_Layer - ERR05')
end subroutine Read_Layer
!++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
!++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
!MODIFIER MODIFIER MODIFIER MODIFIER MODIFIER MODIFIER MODIFIER MODIFIER MODI
!++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
!++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
subroutine ModifyStatistic (StatisticID, Value2D, Value3D, WaterPoints2D, &
WaterPoints3D, Now, STAT)
!Arguments---------------------------------------------------------------
integer :: StatisticID
real, dimension(:, : ), pointer, optional :: Value2D
real, dimension(:, :, :), pointer, optional :: Value3D
!real(4), dimension(:, :, :), pointer, optional :: Value3D_R4
integer, dimension(:, : ), pointer, optional :: WaterPoints2D
integer, dimension(:, :, :), pointer, optional :: WaterPoints3D
type (T_Time ), optional :: Now
integer, optional :: STAT
!Local-------------------------------------------------------------------
integer :: STAT_, ready_
integer :: STAT_CALL
!Monitores Performance
if (MonitorPerformance) call StartWatch ("ModuleStatistic", "ModifyStatistic")
STAT_ = UNKNOWN_
call Ready (StatisticID, ready_)
if (ready_ .EQ. IDLE_ERR_) then
if (present(Now)) then
Me%ExternalVar%Now = Now
else
call GetComputeCurrentTime(Me%ObjTime, Me%ExternalVar%Now, STAT = STAT_CALL)
if (STAT_CALL /= SUCCESS_) stop 'ModifyStatistic - ModuleStatistic - ERR01'
endif
if (Me%Methodology == Value2DStat2D_) then
if (present(Value2D)) then
call StatisticAnalysis2D (Value2D, WaterPoints2D)
STAT_ = SUCCESS_
else
STAT_ = UNKNOWN_
endif
endif
if (Me%Methodology == Value3DStat3D_) then
if (present(Value3D)) then
Me%ExternalVar%Value3D => Value3D
call StatisticAnalysis3D (WaterPoints3D)
nullify(Me%ExternalVar%Value3D)
STAT_ = SUCCESS_
else
STAT_ = UNKNOWN_
endif
endif
if (Me%Methodology == Value3DStatLayers_) then
!João Sobrinho------------------------------------------------------------
if (present(Value3D)) then
call StatisticAnalysis3D (WaterPoints3D)
STAT_ = SUCCESS_
else
STAT_ = UNKNOWN_
endif
endif
else
STAT_ = ready_
endif
if (MonitorPerformance) call StopWatch ("ModuleStatistic", "ModifyStatistic")
if (present(STAT)) STAT = STAT_
end subroutine ModifyStatistic
!--------------------------------------------------------------------------
subroutine AddStatisticLayers (StatisticID, Value3D, WaterPoints3D, DZ3D, &
LayerNumber, UpperDepth, LowerDepth, UpperLayer, &
LowerLayer, Now, STAT)
!Arguments---------------------------------------------------------------
integer :: StatisticID
real, dimension(:, :, :), pointer, optional :: Value3D
integer, dimension(:, :, :), pointer :: WaterPoints3D
real, dimension(:, :, :), pointer :: DZ3D
integer :: LayerNumber
real, dimension(:, : ), pointer, optional :: UpperDepth, LowerDepth
integer, dimension(:, : ), pointer, optional :: UpperLayer, LowerLayer
type (T_Time ), optional :: Now
integer, optional :: STAT
!Local-------------------------------------------------------------------
integer :: STAT_, ready_
integer :: STAT_CALL
!Monitores Performance
if (MonitorPerformance) call StartWatch ("ModuleStatistic", "AddStatisticLayers")
STAT_ = UNKNOWN_
call Ready (StatisticID, ready_)
if (ready_ .EQ. IDLE_ERR_) then
if (present(Now)) then
Me%ExternalVar%Now = Now
else
call GetComputeCurrentTime(Me%ObjTime, Me%ExternalVar%Now, STAT = STAT_CALL)
if (STAT_CALL /= SUCCESS_) stop 'AddStatisticLayers - ModuleStatistic - ERR01'
endif
if (Me%Methodology == Value3DStatLayers_) then
if (Me%Layers%Definition(LayerNumber) == Depth_) then
if (present(UpperDepth)) &
Me%Layers%UpperDepth(:,:, LayerNumber) = UpperDepth(:,:)
if (present(LowerDepth)) &
Me%Layers%LowerDepth(:,:, LayerNumber) = LowerDepth(:,:)
call AverageValueBetweenDepths(Value3D = Value3D, WaterPoints3D = WaterPoints3D, &
DZ3D = DZ3D, LayerNumber = LayerNumber)
else if (Me%Layers%Definition(LayerNumber) == Layer_) then
if (present(UpperLayer)) &
Me%Layers%UpperLayer(:,:, LayerNumber) = UpperLayer(:,:)
if (present(LowerLayer)) &
Me%Layers%LowerLayer(:,:, LayerNumber) = LowerLayer(:,:)
call AverageValueBetweenLayers(Value3D = Value3D, WaterPoints3D = WaterPoints3D, &
DZ3D = DZ3D, LayerNumber = LayerNumber)
endif
STAT_ = SUCCESS_
else
STAT_ = UNKNOWN_
endif
else
STAT_ = ready_
endif
if (MonitorPerformance) call StopWatch ("ModuleStatistic", "AddStatisticLayers")
if (present(STAT)) STAT = STAT_
end subroutine AddStatisticLayers
!--------------------------------------------------------------------------
subroutine AverageValueBetweenDepths(Value3D, WaterPoints3D, DZ3D, LayerNumber)
!Arguments-------------------------------------------------------------
real, dimension(:, :, :), pointer :: Value3D
integer, dimension(:, :, :), pointer :: WaterPoints3D
real, dimension(:, :, :), pointer :: DZ3D
integer :: LayerNumber
!Local-----------------------------------------------------------------
real :: DepthTopLayer, DepthBottomLayer, &
DZ_Total, DZ, DepthMin, DepthMax
integer :: ILB, IUB, i
integer :: JLB, JUB, j
integer :: KLB, KUB, k
!Begin-----------------------------------------------------------------
!Shorten
ILB = Me%ExternalVar%WorkSize%ILB
IUB = Me%ExternalVar%WorkSize%IUB
JLB = Me%ExternalVar%WorkSize%JLB
JUB = Me%ExternalVar%WorkSize%JUB
KLB = Me%ExternalVar%WorkSize%KLB
KUB = Me%ExternalVar%WorkSize%KUB
do j = JLB, JUB
do i = ILB, IUB
DepthTopLayer = 0.
DepthBottomLayer = 0.
DZ_Total = 0.
Me%Layers%WaterPoints (i, j, LayerNumber) = 0
Me%Layers%Value (i, j, LayerNumber) = 0
do k = KUB, KLB, -1
DepthMin = Me%Layers%UpperDepth(i, j, LayerNumber)
DepthMax = Me%Layers%LowerDepth(i, j, LayerNumber)
if (WaterPoints3D(i,j,k) == WaterPoint) then
DepthTopLayer = DepthBottomLayer
DepthBottomLayer = DepthBottomLayer + DZ3D(i, j, k)
!If the layer intercept the vertical interval predefined than an average value can be compute
if (DepthTopLayer <= DepthMax .and. DepthBottomLayer >= DepthMin) then
Me%Layers%WaterPoints(i, j, LayerNumber) = 1
if (DepthTopLayer >= DepthMin .and. DepthBottomLayer <= DepthMax) DZ = DZ3D(i, j, k)
if (DepthTopLayer < DepthMin .and. DepthBottomLayer <= DepthMax) DZ = DepthBottomLayer - DepthMin
if (DepthTopLayer >= DepthMin .and. DepthBottomLayer > DepthMax) DZ = DepthMax - DepthTopLayer
if ((DZ_Total + DZ)>0 ) then
Me%Layers%Value(i, j, LayerNumber) = &
(Me%Layers%Value(i, j, LayerNumber) * DZ_Total + &
Value3D(i, j, k) * DZ) / (DZ_Total + DZ)
else
Me%Layers%Value(i, j, LayerNumber) = Value3D(i, j, k)
endif
endif
DZ_Total = DZ_Total + DZ
endif
enddo
enddo
enddo
end subroutine AverageValueBetweenDepths
!--------------------------------------------------------------------------
subroutine AverageValueBetweenLayers(Value3D, WaterPoints3D, DZ3D, LayerNumber)
!Arguments-------------------------------------------------------------
real, dimension(:, :, :), pointer :: Value3D
integer, dimension(:, :, :), pointer :: WaterPoints3D
real, dimension(:, :, :), pointer :: DZ3D
integer :: LayerNumber
!Local-----------------------------------------------------------------
real :: DepthTopLayer, DepthBottomLayer, &
DZ_Total
integer :: ILB, IUB, i
integer :: JLB, JUB, j
integer :: KUB, k, kup, klo
!Begin-----------------------------------------------------------------
!Shorten
ILB = Me%ExternalVar%WorkSize%ILB
IUB = Me%ExternalVar%WorkSize%IUB
JLB = Me%ExternalVar%WorkSize%JLB
JUB = Me%ExternalVar%WorkSize%JUB
KUB = Me%ExternalVar%WorkSize%KUB
do j = JLB, JUB
do i = ILB, IUB
if (WaterPoints3D(i,j,KUB) == WaterPoint) then
DepthTopLayer = 0.
DepthBottomLayer = 0.
DZ_Total = 0.
Me%Layers%WaterPoints (i, j, LayerNumber) = 0
Me%Layers%Value (i, j, LayerNumber) = 0
kup = Me%Layers%UpperLayer(i, j, LayerNumber)
klo = Me%Layers%LowerLayer(i, j, LayerNumber)
do k = kup, klo, -1
if (WaterPoints3D(i,j,k) == WaterPoint) then
Me%Layers%WaterPoints(i, j, LayerNumber) = 1
Me%Layers%Value(i, j, LayerNumber) = &
(Me%Layers%Value(i, j, LayerNumber) * DZ_Total + &
Value3D(i, j, k) * DZ3D(i, j, k)) / (DZ_Total + DZ3D(i, j, k))
DZ_Total = DZ_Total + DZ3D(i, j, k)
endif
enddo
endif
enddo
enddo
end subroutine AverageValueBetweenLayers
!--------------------------------------------------------------------------
subroutine StatisticAnalysis3D (WaterPoints3D)
!Arguments-------------------------------------------------------------
integer, dimension(:, :, :), pointer :: WaterPoints3D
!Local-----------------------------------------------------------------
real, dimension(:, :, :), pointer :: Value3D
integer :: KLB, KUB
integer, dimension(:, :, :), pointer :: lWaterPoints3D
!Begin-----------------------------------------------------------------
if (Me%Methodology == Value3DStat3D_) then
Value3D => Me%ExternalVar%Value3D
lWaterPoints3D => WaterPoints3D
KLB = Me%ExternalVar%WorkSize%KLB
KUB = Me%ExternalVar%WorkSize%KUB
else if (Me%Methodology == Value3DStatLayers_) then
Value3D => Me%Layers%Value
lWaterPoints3D => Me%Layers%WaterPoints
KLB = 1
KUB = Me%Layers%Number
endif
!Global
if (Me%Global%On) then
call ModifyGlobalStatistic (Value3D, lWaterPoints3D, KLB, KUB)
endif
!Daily
if (Me%Daily%On) then
call ModifyDailyStatistic (Value3D, lWaterPoints3D, KLB, KUB)
endif
!Monthly
if (Me%Monthly%On) then
call ModifyMonthlyStatistic (Value3D, lWaterPoints3D, KLB, KUB)
endif
!SpecificHour
if (Me%SpecificHour%On) then
call ModifySpecificHourStatistic (Value3D, lWaterPoints3D, KLB, KUB)
endif
!Classification
if (Me%Classification%On) then
call ModifyClassification (Value3D, lWaterPoints3D, KLB, KUB)
endif
nullify(Value3D)
nullify(lWaterPoints3D)
end subroutine StatisticAnalysis3D
!------------------------------------------------------------------------------------------------------
!subroutine StatisticAnalysis3D_R4 (WaterPoints3D)
!
! !Arguments-------------------------------------------------------------
! integer, dimension(:, :, :), pointer :: WaterPoints3D
!
! !Local-----------------------------------------------------------------
! real(4), dimension(:, :, :), pointer :: Value3D_R4
! integer :: KLB, KUB
! integer, dimension(:, :, :), pointer :: lWaterPoints3D
!
! !Begin-----------------------------------------------------------------
!
! if (Me%Methodology == Value3DStat3D_) then
!
! Value3D_R4 => Me%ExternalVar%Value3D_R4
! lWaterPoints3D => WaterPoints3D
! KLB = Me%ExternalVar%WorkSize%KLB
! KUB = Me%ExternalVar%WorkSize%KUB
!
! else if (Me%Methodology == Value3DStatLayers_) then
!
! Value3D_R4 => Me%Layers%Value_R4
! lWaterPoints3D => Me%Layers%WaterPoints
! KLB = 1
! KUB = Me%Layers%Number
!
! endif
!
! !Global
! if (Me%Global%On) then
! call ModifyGlobalStatistic_R4 (Value3D_R4, lWaterPoints3D, KLB, KUB)
! endif
!
! !Daily
! if (Me%Daily%On) then
! call ModifyDailyStatistic_R4 (Value3D_R4, lWaterPoints3D, KLB, KUB)
! endif
!
! !Monthly
! if (Me%Monthly%On) then
! call ModifyMonthlyStatistic_R4 (Value3D_R4, lWaterPoints3D, KLB, KUB)
! endif
!
! !SpecificHour
! if (Me%SpecificHour%On) then
! call ModifySpecificHourStatistic_R4 (Value3D_R4, lWaterPoints3D, KLB, KUB)
! endif
!
! !Classification
! if (Me%Classification%On) then
! call ModifyClassification_R4 (Value3D_R4, lWaterPoints3D, KLB, KUB)
! endif
!
! nullify(Value3D_R4)
! nullify(lWaterPoints3D)
!
!end subroutine StatisticAnalysis3D_R4
!--------------------------------------------------------------------------
subroutine StatisticAnalysis2D (Value2D, WaterPoints2D)
!Arguments-------------------------------------------------------------
real, dimension(:, :), pointer :: Value2D
integer, dimension(:, :), pointer :: WaterPoints2D
!Local-----------------------------------------------------------------
!Global
if (Me%Global%On) then
call ModifyGlobalStatistic2D (Value2D, WaterPoints2D)
endif
!Daily
if (Me%Daily%On) then
call ModifyDailyStatistic2D (Value2D, WaterPoints2D)
endif
!Monthly
if (Me%Monthly%On) then
call ModifyMonthlyStatistic2D (Value2D, WaterPoints2D)
endif
!SpecificHour
if (Me%SpecificHour%On) then
call ModifySpecificHourStatistic2D (Value2D, WaterPoints2D)
endif
!Classification
if (Me%Classification%On) then
call ModifyClassification2D (Value2D, WaterPoints2D)
endif
end subroutine StatisticAnalysis2D
!--------------------------------------------------------------------------
subroutine ModifyGlobalStatistic (Value, WaterPoints3D, KLB, KUB)
!Arguments-------------------------------------------------------------
real, dimension(:, :, :), pointer :: Value
integer, dimension(:, :, :), pointer :: WaterPoints3D
integer :: KLB, KUB
!Local-----------------------------------------------------------------
integer :: ILB, IUB, i
integer :: JLB, JUB, j
integer :: k
real :: DT, DX, AuxValue
!Shorten
ILB = Me%ExternalVar%WorkSize%ILB
IUB = Me%ExternalVar%WorkSize%IUB
JLB = Me%ExternalVar%WorkSize%JLB
JUB = Me%ExternalVar%WorkSize%JUB
!Time since last calculation
DT = Me%ExternalVar%Now - Me%Global%LastCalculation
cd1: if (DT>0) then
!Loops
do k = KLB, KUB
do j = JLB, JUB
do i = ILB, IUB
if (WaterPoints3D(i, j, k) == WaterPoint) then
!Minimum Value
if (Value (i, j, k) < Me%Global%Minimum (i, j, k)) &
Me%Global%Minimum (i, j, k) = Value (i, j, k)
!Maximum Value
if (Value (i, j, k) > Me%Global%Maximum (i, j, k)) &
Me%Global%Maximum (i, j, k) = Value (i, j, k)
!Minimum Value
if (Value (i, j, k) < Me%Global%Minimum2D (i, j)) &
Me%Global%Minimum2D (i, j) = Value (i, j, k)
!Maximum Value
if (Value (i, j, k) > Me%Global%Maximum2D (i, j)) &
Me%Global%Maximum2D (i, j) = Value (i, j, k)
!Average
Me%Global%Average (i, j, k) = &
(Me%Global%Average (i, j, k) * &
Me%Global%RunPeriod + &
Value (i, j, k) * DT) / (Me%Global%RunPeriod + DT)
!Square Average
Me%Global%SquareAverage (i, j, k) = &
(Me%Global%SquareAverage (i, j, k) * &
Me%Global%RunPeriod + &
Value (i, j, k)**2 * DT) / (Me%Global%RunPeriod + DT)
!Standard deviation
DX = Me%Global%SquareAverage (i, j, k) - &
Me%Global%Average (i, j, k) ** 2
DX = abs(DX)
Me%Global%StandardDeviation(i, j, k) = sqrt(DX)
!Accumulated Values
if (Me%Accumulated) &
Me%Global%Accumulated (i,j,k) = Me%Global%Accumulated (i,j,k) &
+ Value (i, j, k)
if (Me%GeomMean) then !Geometric Average to calculate
!Geometric Average
AuxValue = Value (i, j, k)
if (AuxValue == 0.0) then
AuxValue = 1.0
elseif (AuxValue .lt. 0.0) then
write(*,*) 'Negative valued property.'
write(*,*) 'Geometric Average cannot be calculated.'
stop 'ModifyGlobalStatistic - ModuleStatistic - ERR01'
endif
Me%Global%GeomAverage (i, j, k) = 10** &
((LOG10(Me%Global%GeomAverage (i, j, k)) * &
Me%Global%RunPeriod + &
LOG10(AuxValue) * DT) / (Me%Global%RunPeriod + DT))
!Squared Geometric Average
Me%Global%SquareGeomAverage (i, j, k) = &
(Me%Global%SquareGeomAverage (i, j, k) * &
Me%Global%RunPeriod + &
AuxValue**2 * DT) / (Me%Global%RunPeriod + DT)
!Geometric Standard Deviation
DX = Me%Global%SquareGeomAverage (i, j, k) - &
Me%Global%GeomAverage (i, j, k) ** 2
DX = abs(DX)
Me%Global%GeomStandardDeviation(i, j, k) = sqrt(DX)
endif
!guillaume juan
if (Me%Critical .and. (Value (i, j, k) > FillValueReal / 2.)) then
if (Value (i, j, k) < Me%CriticalValue) then
Me%Global%PercentBelowCriticalValue (i, j, k) = &
(Me%Global%PercentBelowCriticalValue (i, j, k) * &
Me%Global%RunPeriod + DT) / (Me%Global%RunPeriod + DT)
else
Me%Global%PercentBelowCriticalValue (i, j, k) = &
(Me%Global%PercentBelowCriticalValue (i, j, k) * &
Me%Global%RunPeriod + 0.) / (Me%Global%RunPeriod + DT)
endif
endif
endif
enddo
enddo
enddo
!Updates Time
Me%Global%RunPeriod = Me%Global%RunPeriod + DT
Me%Global%LastCalculation = Me%ExternalVar%Now
endif cd1
end subroutine ModifyGlobalStatistic
!--------------------------------------------------------------------------
!--------------------------------------------------------------------------
! subroutine ModifyGlobalStatistic_R4 (Value_R4, WaterPoints3D, KLB, KUB)
!
! !Arguments-------------------------------------------------------------
! real(4), dimension(:, :, :), pointer :: Value_R4
! integer, dimension(:, :, :), pointer :: WaterPoints3D
! integer :: KLB, KUB
!
! !Local-----------------------------------------------------------------
! integer :: ILB, IUB, i
! integer :: JLB, JUB, j
! integer :: k
! real :: DT
! real(4) :: AuxValue, DX
!
! !Shorten
! ILB = Me%ExternalVar%WorkSize%ILB
! IUB = Me%ExternalVar%WorkSize%IUB
! JLB = Me%ExternalVar%WorkSize%JLB
! JUB = Me%ExternalVar%WorkSize%JUB
!
!
! !Time since last calculation
! DT = Me%ExternalVar%Now - Me%Global%LastCalculation
!
!cd1: if (DT>0) then
!
! !Loops
! do k = KLB, KUB
! do j = JLB, JUB
! do i = ILB, IUB
! if (WaterPoints3D(i, j, k) == WaterPoint) then
!
! !Minimum Value
! if (Value_R4 (i, j, k) < Me%Global%Minimum_R4 (i, j, k)) &
! Me%Global%Minimum_R4 (i, j, k) = Value_R4 (i, j, k)
!
! !Maximum Value
! if (Value_R4 (i, j, k) > Me%Global%Maximum_R4 (i, j, k)) &
! Me%Global%Maximum_R4 (i, j, k) = Value_R4 (i, j, k)
!
! !Average
! Me%Global%Average_R4 (i, j, k) = &
! (Me%Global%Average_R4 (i, j, k) * &
! Me%Global%RunPeriod + &
! Value_R4 (i, j, k) * DT) / (Me%Global%RunPeriod + DT)
!
! !Square Average
! Me%Global%SquareAverage_R4 (i, j, k) = &
! (Me%Global%SquareAverage_R4 (i, j, k) * &
! Me%Global%RunPeriod + &
! Value_R4 (i, j, k)**2 * DT) / (Me%Global%RunPeriod + DT)
!
! !Standard deviation
! DX = Me%Global%SquareAverage_R4 (i, j, k) - &
! Me%Global%Average_R4 (i, j, k) ** 2
!
! DX = abs(DX)
!
! Me%Global%StandardDeviation_R4(i, j, k) = sqrt(DX)
!
! !Accumulated Values
! if (Me%Accumulated) &
! Me%Global%Accumulated_R4 (i,j,k) = Me%Global%Accumulated_R4 (i,j,k) &
! + Value_R4 (i, j, k)
!
! if (Me%GeomMean) then !Geometric Average to calculate
! !Geometric Average
! AuxValue = Value_R4 (i, j, k)
! if (AuxValue == 0.0) then
! AuxValue = 1.0
! elseif (AuxValue .lt. 0.0) then
! write(*,*) 'Negative valued property.'
! write(*,*) 'Geometric Average cannot be calculated.'
! stop 'ModifyGlobalStatistic_R4 - ModuleStatistic - ERR01'
! endif
!
! Me%Global%GeomAverage_R4 (i, j, k) = 10** &
! ((LOG10(Me%Global%GeomAverage_R4 (i, j, k)) * &
! Me%Global%RunPeriod + &
! LOG10(AuxValue) * DT) / (Me%Global%RunPeriod + DT))
!
! !Squared Geometric Average
! Me%Global%SquareGeomAverage_R4 (i, j, k) = &
! (Me%Global%SquareGeomAverage_R4 (i, j, k) * &
! Me%Global%RunPeriod + &
! AuxValue**2 * DT) / (Me%Global%RunPeriod + DT)
!
! !Geometric Standard Deviation
! DX = Me%Global%SquareGeomAverage_R4 (i, j, k) - &
! Me%Global%GeomAverage_R4 (i, j, k) ** 2
!
! DX = abs(DX)
!
! Me%Global%GeomStandardDeviation_R4(i, j, k) = sqrt(DX)
!
! endif
!
! !guillaume juan
! if (Me%Critical .and. (Value_R4 (i, j, k) > FillValueReal / 2.)) then
! if (Value_R4 (i, j, k) < Me%CriticalValue) then
! Me%Global%PercentBelowCriticalValue_R4 (i, j, k) = &
! (Me%Global%PercentBelowCriticalValue_R4 (i, j, k) * &
! Me%Global%RunPeriod + DT) / (Me%Global%RunPeriod + DT)
! else
! Me%Global%PercentBelowCriticalValue_R4 (i, j, k) = &
! (Me%Global%PercentBelowCriticalValue_R4 (i, j, k) * &
! Me%Global%RunPeriod + 0.) / (Me%Global%RunPeriod + DT)
! endif
! endif
!
! endif
! enddo
! enddo
! enddo
!
! !Updates Time
! Me%Global%RunPeriod = Me%Global%RunPeriod + DT
! Me%Global%LastCalculation = Me%ExternalVar%Now
!
! endif cd1
!
! end subroutine ModifyGlobalStatistic_R4
!-------------------------------------------------------------------------------
subroutine ModifyDailyStatistic (Value, WaterPoints3D, KLB, KUB)
!Arguments-------------------------------------------------------------
real, dimension(:, :, :), pointer :: Value
integer, dimension(:, :, :), pointer :: WaterPoints3D
integer :: KLB, KUB
!Local-----------------------------------------------------------------
integer :: ILB, IUB, i
integer :: JLB, JUB, j
integer :: k
real :: DT, DX, AuxValue
real :: OldDay, PresentDay
!Shorten
ILB = Me%ExternalVar%WorkSize%ILB
IUB = Me%ExternalVar%WorkSize%IUB
JLB = Me%ExternalVar%WorkSize%JLB
JUB = Me%ExternalVar%WorkSize%JUB
!Time since last calculation
DT = Me%ExternalVar%Now - Me%Daily%LastCalculation
cd1: if (DT>0) then
!Loops
do k = KLB, KUB
do j = JLB, JUB
do i = ILB, IUB
if (WaterPoints3D(i, j, k) == WaterPoint) then
!Minimum Value
if (Value (i, j, k) < Me%Daily%Minimum (i, j, k)) &
Me%Daily%Minimum (i, j, k) = Value (i, j, k)
!Maximum Value
if (Value (i, j, k) > Me%Daily%Maximum (i, j, k)) &
Me%Daily%Maximum (i, j, k) = Value (i, j, k)
!Average
Me%Daily%Average (i, j, k) = &
(Me%Daily%Average (i, j, k) * &
Me%Daily%RunPeriod + &
Value (i, j, k) * DT) / (Me%Daily%RunPeriod + DT)
!Square Average
Me%Daily%SquareAverage (i, j, k) = &
(Me%Daily%SquareAverage (i, j, k) * &
Me%Daily%RunPeriod + &
Value (i, j, k)**2 * DT) / (Me%Daily%RunPeriod + DT)
!Standard deviation
DX = Me%Daily%SquareAverage (i, j, k) - &
Me%Daily%Average (i, j, k) ** 2
DX = abs(DX)
Me%Daily%StandardDeviation(i, j, k) = sqrt(DX)
!Accumulated Values
if (Me%Accumulated) &
Me%Daily%Accumulated (i,j,k) = Me%Daily%Accumulated (i,j,k) &
+ Value (i, j, k)
if (Me%GeomMean) then !Geometric Average to calculate
!Geometric Average
AuxValue = Value (i, j, k)
if (AuxValue == 0.0) then
AuxValue = 1.0
elseif (AuxValue .lt. 0.0) then
write(*,*) 'Negative valued property.'
write(*,*) 'Geometric Average cannot be calculated.'
stop 'ModifyDailyStatistic - ModuleStatistic - ERR01'
endif
if (Me%Daily%GeomAverage (i, j, k) == 0.0) then
Me%Daily%GeomAverage (i, j, k) = 1.0
endif
Me%Daily%GeomAverage (i, j, k) = 10** &
((LOG10(Me%Daily%GeomAverage (i, j, k)) * &
Me%Daily%RunPeriod + &
LOG10(AuxValue) * DT) / (Me%Daily%RunPeriod + DT))
!Squared Geometric Average
Me%Daily%SquareGeomAverage (i, j, k) = &
(Me%Daily%SquareGeomAverage (i, j, k) * &
Me%Daily%RunPeriod + &
AuxValue**2 * DT) / (Me%Daily%RunPeriod + DT)
!Geometric Standard Deviation
DX = Me%Daily%SquareGeomAverage (i, j, k) - &
Me%Daily%GeomAverage (i, j, k) ** 2
DX = abs(DX)
Me%Daily%GeomStandardDeviation(i, j, k) = sqrt(DX)
endif
!guillaume juan
if (Me%Critical .and. (Value (i, j, k) > FillValueReal / 2.)) then
if (Value (i, j, k) < Me%CriticalValue) then
Me%Daily%PercentBelowCriticalValue (i, j, k) = &
(Me%Daily%PercentBelowCriticalValue (i, j, k) * &
Me%Daily%RunPeriod + DT) / (Me%Daily%RunPeriod + DT)
else
Me%Daily%PercentBelowCriticalValue (i, j, k) = &
(Me%Daily%PercentBelowCriticalValue (i, j, k) * &
Me%Daily%RunPeriod + 0.) / (Me%Daily%RunPeriod + DT)
endif
endif
endif
enddo
enddo
enddo
!Verifies if the present time is a new output
call ExtractDate (Time1=Me%ExternalVar%Now, Day = PresentDay)
call ExtractDate (Time1=Me%Daily%LastCalculation, Day = OldDay)
if (int(PresentDay) /= int(OldDay)) then
call WriteValuesToFileHDF5 (.false., .true., .false., .false., .false.)
Me%Daily%Minimum = Value
Me%Daily%Maximum = Value
Me%Daily%Average = Value
Me%Daily%SquareAverage = Value
Me%Daily%StandardDeviation = Value
if (Me%Accumulated) Me%Daily%Accumulated = 0.0
if (Me%GeomMean) then !Geometric Average to calculate
Me%Daily%GeomAverage = Value
Me%Daily%SquareGeomAverage = Value
Me%Daily%GeomStandardDeviation = Value
endif
!guillaume juan
if (Me%Critical) then
do k = KLB, KUB
do j = JLB, JUB
do i = ILB, IUB
if ( (WaterPoints3D(i, j, k) == WaterPoint) .and. (Value (i, j, k) > FillValueReal / 2.)) then
if(Value(i, j, k) < Me%CriticalValue) then
Me%Daily%PercentBelowCriticalValue (i, j, k) = 1.
else
Me%Daily%PercentBelowCriticalValue (i, j, k) = 0.
endif
endif
enddo
enddo
enddo
endif
Me%Daily%RunPeriod = 0.
else
Me%Daily%RunPeriod = Me%Daily%RunPeriod + DT
endif
!Updates Time
Me%Daily%LastCalculation = Me%ExternalVar%Now
endif cd1
end subroutine ModifyDailyStatistic
!--------------------------------------------------------------------------
! subroutine ModifyDailyStatistic_R4 (Value_R4, WaterPoints3D, KLB, KUB)
!
! !Arguments-------------------------------------------------------------
! real(4), dimension(:, :, :), pointer :: Value_R4
! integer, dimension(:, :, :), pointer :: WaterPoints3D
! integer :: KLB, KUB
!
! !Local-----------------------------------------------------------------
! integer :: ILB, IUB, i
! integer :: JLB, JUB, j
! integer :: k
! real :: DT
! real(4) :: DX, AuxValue
! real :: OldDay, PresentDay
!
! !Shorten
! ILB = Me%ExternalVar%WorkSize%ILB
! IUB = Me%ExternalVar%WorkSize%IUB
! JLB = Me%ExternalVar%WorkSize%JLB
! JUB = Me%ExternalVar%WorkSize%JUB
!
!
! !Time since last calculation
! DT = Me%ExternalVar%Now - Me%Daily%LastCalculation
!
!cd1: if (DT>0) then
!
! !Loops
! do k = KLB, KUB
! do j = JLB, JUB
! do i = ILB, IUB
! if (WaterPoints3D(i, j, k) == WaterPoint) then
!
! !Minimum Value
! if (Value_R4 (i, j, k) < Me%Daily%Minimum_R4 (i, j, k)) &
! Me%Daily%Minimum_R4 (i, j, k) = Value_R4 (i, j, k)
!
! !Maximum Value_R4
! if (Value_R4 (i, j, k) > Me%Daily%Maximum_R4 (i, j, k)) &
! Me%Daily%Maximum_R4 (i, j, k) = Value_R4 (i, j, k)
!
! !Average
! Me%Daily%Average_R4 (i, j, k) = &
! (Me%Daily%Average_R4 (i, j, k) * &
! Me%Daily%RunPeriod + &
! Value_R4 (i, j, k) * DT) / (Me%Daily%RunPeriod + DT)
!
! !Square Average
! Me%Daily%SquareAverage_R4 (i, j, k) = &
! (Me%Daily%SquareAverage_R4 (i, j, k) * &
! Me%Daily%RunPeriod + &
! Value_R4 (i, j, k)**2 * DT) / (Me%Daily%RunPeriod + DT)
!
! !Standard deviation
! DX = Me%Daily%SquareAverage_R4 (i, j, k) - &
! Me%Daily%Average_R4 (i, j, k) ** 2
!
! DX = abs(DX)
!
! Me%Daily%StandardDeviation_R4(i, j, k) = sqrt(DX)
!
! !Accumulated Values
! if (Me%Accumulated) &
! Me%Daily%Accumulated_R4 (i,j,k) = Me%Daily%Accumulated_R4 (i,j,k) &
! + Value_R4 (i, j, k)
!
! if (Me%GeomMean) then !Geometric Average to calculate
! !Geometric Average
! AuxValue = Value_R4 (i, j, k)
! if (AuxValue == 0.0) then
! AuxValue = 1.0
! elseif (AuxValue .lt. 0.0) then
! write(*,*) 'Negative valued property.'
! write(*,*) 'Geometric Average cannot be calculated.'
! stop 'ModifyDailyStatistic_R4 - ModuleStatistic - ERR01'
! endif
!
! if (Me%Daily%GeomAverage_R4 (i, j, k) == 0.0) then
! Me%Daily%GeomAverage_R4 (i, j, k) = 1.0
! endif
!
! Me%Daily%GeomAverage_R4 (i, j, k) = 10** &
! ((LOG10(Me%Daily%GeomAverage_R4 (i, j, k)) * &
! Me%Daily%RunPeriod + &
! LOG10(AuxValue) * DT) / (Me%Daily%RunPeriod + DT))
!
! !Squared Geometric Average
! Me%Daily%SquareGeomAverage_R4 (i, j, k) = &
! (Me%Daily%SquareGeomAverage_R4 (i, j, k) * &
! Me%Daily%RunPeriod + &
! AuxValue**2 * DT) / (Me%Daily%RunPeriod + DT)
!
! !Geometric Standard Deviation
! DX = Me%Daily%SquareGeomAverage_R4 (i, j, k) - &
! Me%Daily%GeomAverage_R4 (i, j, k) ** 2
!
! DX = abs(DX)
!
! Me%Daily%GeomStandardDeviation_R4(i, j, k) = sqrt(DX)
!
! endif
!
! !guillaume juan
! if (Me%Critical .and. (Value_R4 (i, j, k) > FillValueReal / 2.)) then
! if (Value_R4 (i, j, k) < Me%CriticalValue) then
! Me%Daily%PercentBelowCriticalValue_R4 (i, j, k) = &
! (Me%Daily%PercentBelowCriticalValue_R4 (i, j, k) * &
! Me%Daily%RunPeriod + DT) / (Me%Daily%RunPeriod + DT)
! else
! Me%Daily%PercentBelowCriticalValue_R4 (i, j, k) = &
! (Me%Daily%PercentBelowCriticalValue_R4 (i, j, k) * &
! Me%Daily%RunPeriod + 0.) / (Me%Daily%RunPeriod + DT)
! endif
! endif
!
! endif
! enddo
! enddo
! enddo
!
! !Verifies if the present time is a new output
! call ExtractDate (Time1=Me%ExternalVar%Now, Day = PresentDay)
! call ExtractDate (Time1=Me%Daily%LastCalculation, Day = OldDay)
! if (int(PresentDay) /= int(OldDay)) then
! call WriteValuesToFileHDF5_R4 (.false., .true., .false., .false., .false.)
! Me%Daily%Minimum_R4 = Value_R4
! Me%Daily%Maximum_R4 = Value_R4
! Me%Daily%Average_R4 = Value_R4
! Me%Daily%SquareAverage_R4 = Value_R4
! Me%Daily%StandardDeviation_R4 = Value_R4
!
! if (Me%Accumulated) Me%Daily%Accumulated_R4 = 0.0
!
! if (Me%GeomMean) then !Geometric Average to calculate
! Me%Daily%GeomAverage_R4 = Value_R4
! Me%Daily%SquareGeomAverage_R4 = Value_R4
! Me%Daily%GeomStandardDeviation_R4 = Value_R4
! endif
!
! !guillaume juan
! if (Me%Critical) then
! do k = KLB, KUB
! do j = JLB, JUB
! do i = ILB, IUB
! if ( (WaterPoints3D(i, j, k) == WaterPoint) .and. (Value_R4 (i, j, k) > FillValueReal / 2.)) then
! if(Value_R4(i, j, k) < Me%CriticalValue) then
! Me%Daily%PercentBelowCriticalValue_R4 (i, j, k) = 1.
! else
! Me%Daily%PercentBelowCriticalValue_R4 (i, j, k) = 0.
! endif
! endif
! enddo
! enddo
! enddo
! endif
!
! Me%Daily%RunPeriod = 0.
! else
! Me%Daily%RunPeriod = Me%Daily%RunPeriod + DT
! endif
!
! !Updates Time
! Me%Daily%LastCalculation = Me%ExternalVar%Now
!
!
! endif cd1
!
!
! end subroutine ModifyDailyStatistic_R4
!--------------------------------------------------------------------------
subroutine ModifyMonthlyStatistic (Value, WaterPoints3D, KLB, KUB)
!Arguments-------------------------------------------------------------
real, dimension(:, :, :), pointer :: Value
integer, dimension(:, :, :), pointer :: WaterPoints3D
integer :: KLB, KUB
!Local-----------------------------------------------------------------
integer :: ILB, IUB, i
integer :: JLB, JUB, j
integer :: k
real :: DT, DX, AuxValue
real :: OldMonth, PresentMonth
!Shorten
ILB = Me%ExternalVar%WorkSize%ILB
IUB = Me%ExternalVar%WorkSize%IUB
JLB = Me%ExternalVar%WorkSize%JLB
JUB = Me%ExternalVar%WorkSize%JUB
!Time since last calculation
DT = Me%ExternalVar%Now - Me%Monthly%LastCalculation
cd1: if (DT>0) then
!Loops
do k = KLB, KUB
do j = JLB, JUB
do i = ILB, IUB
if (WaterPoints3D(i, j, k) == WaterPoint) then
!Minimum Value
if (Value (i, j, k) < Me%Monthly%Minimum (i, j, k)) &
Me%Monthly%Minimum (i, j, k) = Value (i, j, k)
!Maximum Value
if (Value (i, j, k) > Me%Monthly%Maximum (i, j, k)) &
Me%Monthly%Maximum (i, j, k) = Value (i, j, k)
!Average
Me%Monthly%Average (i, j, k) = &
(Me%Monthly%Average (i, j, k) * &
Me%Monthly%RunPeriod + &
Value (i, j, k) * DT) / (Me%Monthly%RunPeriod + DT)
!Square Average
Me%Monthly%SquareAverage (i, j, k) = &
(Me%Monthly%SquareAverage (i, j, k) * &
Me%Monthly%RunPeriod + &
Value (i, j, k)**2 * DT) / (Me%Monthly%RunPeriod + DT)
!Standard deviation
DX = Me%Monthly%SquareAverage (i, j, k) - &
Me%Monthly%Average (i, j, k) ** 2
DX = abs(DX)
Me%Monthly%StandardDeviation(i, j, k) = sqrt(DX)
!Accumulated Values
if (Me%Accumulated) &
Me%Monthly%Accumulated (i,j,k) = Me%Monthly%Accumulated (i,j,k) &
+ Value (i, j, k)
if (Me%GeomMean) then !Geometric Average to calculate
!Geometric Average
AuxValue = Value (i, j, k)
if (AuxValue == 0.0) then
AuxValue = 1.0
elseif (AuxValue .lt. 0.0) then
write(*,*) 'Negative valued property.'
write(*,*) 'Geometric Average cannot be calculated.'
stop 'ModifyMonthlyStatistic - ModuleStatistic - ERR01'
endif
if (Me%Monthly%GeomAverage (i, j, k) == 0.0) then
Me%Monthly%GeomAverage (i, j, k) = 1.0
endif
Me%Monthly%GeomAverage (i, j, k) = 10** &
((LOG10(Me%Monthly%GeomAverage (i, j, k)) * &
Me%Monthly%RunPeriod + &
LOG10(AuxValue) * DT) / (Me%Monthly%RunPeriod + DT))
!Squared Geometric Average
Me%Monthly%SquareGeomAverage (i, j, k) = &
(Me%Monthly%SquareGeomAverage (i, j, k) * &
Me%Monthly%RunPeriod + &
AuxValue**2 * DT) / (Me%Monthly%RunPeriod + DT)
!Geometric Standard Deviation
DX = Me%Monthly%SquareGeomAverage (i, j, k) - &
Me%Monthly%GeomAverage (i, j, k) ** 2
DX = abs(DX)
Me%Monthly%GeomStandardDeviation(i, j, k) = sqrt(DX)
endif
!guillaume juan
if (Me%Critical .and. (Value (i, j, k) > FillValueReal / 2.)) then
if (Value (i, j, k) < Me%CriticalValue) then
Me%Monthly%PercentBelowCriticalValue (i, j, k) = &
(Me%Monthly%PercentBelowCriticalValue (i, j, k) * &
Me%Monthly%RunPeriod + DT) / (Me%Monthly%RunPeriod + DT)
else
Me%Monthly%PercentBelowCriticalValue (i, j, k) = &
(Me%Monthly%PercentBelowCriticalValue (i, j, k) * &
Me%Monthly%RunPeriod + 0.) / (Me%Monthly%RunPeriod + DT)
endif
endif
endif
enddo
enddo
enddo
!Verifies if the present time is a new output
call ExtractDate (Time1=Me%ExternalVar%Now, Month = PresentMonth)
call ExtractDate (Time1=Me%Monthly%LastCalculation, Month = OldMonth)
if (int(PresentMonth) /= int(OldMonth)) then
call WriteValuesToFileHDF5 (.false., .false., .true., .false., .false.)
Me%Monthly%Minimum = Value
Me%Monthly%Maximum = Value
Me%Monthly%Average = Value
Me%Monthly%SquareAverage = Value
Me%Monthly%StandardDeviation = Value
if (Me%Accumulated) Me%Monthly%Accumulated = 0.0
if (Me%GeomMean) then !Geometric Average to calculate
Me%Monthly%GeomAverage = Value
Me%Monthly%SquareGeomAverage = Value
Me%Monthly%GeomStandardDeviation = Value
endif
!guillaume juan
if (Me%Critical) then
do k = KLB, KUB
do j = JLB, JUB
do i = ILB, IUB
if ((WaterPoints3D(i, j, k) == WaterPoint) .and. (Value (i, j, k) > FillValueReal / 2.)) then
if(Value(i, j, k) < Me%CriticalValue) then
Me%Monthly%PercentBelowCriticalValue (i, j, k) = 1.
else
Me%Monthly%PercentBelowCriticalValue (i, j, k) = 0.
endif
endif
enddo
enddo
enddo
endif
Me%Monthly%RunPeriod = 0.
else
Me%Monthly%RunPeriod = Me%Monthly%RunPeriod + DT
endif
!Updates Time
Me%Monthly%LastCalculation = Me%ExternalVar%Now
endif cd1
end subroutine ModifyMonthlyStatistic
!--------------------------------------------------------------------------
! subroutine ModifyMonthlyStatistic_R4 (Value_R4, WaterPoints3D, KLB, KUB)
!
! !Arguments-------------------------------------------------------------
! real(4), dimension(:, :, :), pointer :: Value_R4
! integer, dimension(:, :, :), pointer :: WaterPoints3D
! integer :: KLB, KUB
!
! !Local-----------------------------------------------------------------
! integer :: ILB, IUB, i
! integer :: JLB, JUB, j
! integer :: k
! real :: DT
! real(4) :: DX, AuxValue
! real :: OldMonth, PresentMonth
!
! !Shorten
! ILB = Me%ExternalVar%WorkSize%ILB
! IUB = Me%ExternalVar%WorkSize%IUB
! JLB = Me%ExternalVar%WorkSize%JLB
! JUB = Me%ExternalVar%WorkSize%JUB
!
!
! !Time since last calculation
! DT = Me%ExternalVar%Now - Me%Monthly%LastCalculation
!
!cd1: if (DT>0) then
!
! !Loops
! do k = KLB, KUB
! do j = JLB, JUB
! do i = ILB, IUB
! if (WaterPoints3D(i, j, k) == WaterPoint) then
!
! !Minimum Value
! if (Value_R4 (i, j, k) < Me%Monthly%Minimum_R4 (i, j, k)) &
! Me%Monthly%Minimum_R4 (i, j, k) = Value_R4 (i, j, k)
!
! !Maximum Value
! if (Value_R4 (i, j, k) > Me%Monthly%Maximum_R4 (i, j, k)) &
! Me%Monthly%Maximum_R4 (i, j, k) = Value_R4 (i, j, k)
!
! !Average
! Me%Monthly%Average_R4 (i, j, k) = &
! (Me%Monthly%Average_R4 (i, j, k) * &
! Me%Monthly%RunPeriod + &
! Value_R4 (i, j, k) * DT) / (Me%Monthly%RunPeriod + DT)
!
! !Square Average
! Me%Monthly%SquareAverage_R4 (i, j, k) = &
! (Me%Monthly%SquareAverage_R4 (i, j, k) * &
! Me%Monthly%RunPeriod + &
! Value_R4 (i, j, k)**2 * DT) / (Me%Monthly%RunPeriod + DT)
!
! !Standard deviation
! DX = Me%Monthly%SquareAverage_R4 (i, j, k) - &
! Me%Monthly%Average_R4 (i, j, k) ** 2
!
! DX = abs(DX)
!
!
! Me%Monthly%StandardDeviation_R4(i, j, k) = sqrt(DX)
!
!
! !Accumulated Values
! if (Me%Accumulated) &
! Me%Monthly%Accumulated_R4 (i,j,k) = Me%Monthly%Accumulated_R4 (i,j,k) &
! + Value_R4 (i, j, k)
!
!
! if (Me%GeomMean) then !Geometric Average to calculate
! !Geometric Average
! AuxValue = Value_R4 (i, j, k)
! if (AuxValue == 0.0) then
! AuxValue = 1.0
! elseif (AuxValue .lt. 0.0) then
! write(*,*) 'Negative valued property.'
! write(*,*) 'Geometric Average cannot be calculated.'
! stop 'ModifyMonthlyStatistic_R4 - ModuleStatistic - ERR01'
! endif
!
! if (Me%Monthly%GeomAverage_R4 (i, j, k) == 0.0) then
! Me%Monthly%GeomAverage_R4 (i, j, k) = 1.0
! endif
!
! Me%Monthly%GeomAverage_R4 (i, j, k) = 10** &
! ((LOG10(Me%Monthly%GeomAverage_R4 (i, j, k)) * &
! Me%Monthly%RunPeriod + &
! LOG10(AuxValue) * DT) / (Me%Monthly%RunPeriod + DT))
!
! !Squared Geometric Average
! Me%Monthly%SquareGeomAverage_R4 (i, j, k) = &
! (Me%Monthly%SquareGeomAverage_R4 (i, j, k) * &
! Me%Monthly%RunPeriod + &
! AuxValue**2 * DT) / (Me%Monthly%RunPeriod + DT)
!
! !Geometric Standard Deviation
! DX = Me%Monthly%SquareGeomAverage_R4 (i, j, k) - &
! Me%Monthly%GeomAverage_R4 (i, j, k) ** 2
!
! DX = abs(DX)
!
! Me%Monthly%GeomStandardDeviation_R4(i, j, k) = sqrt(DX)
!
! endif
!
! !guillaume juan
! if (Me%Critical .and. (Value_R4 (i, j, k) > FillValueReal / 2.)) then
! if (Value_R4 (i, j, k) < Me%CriticalValue) then
! Me%Monthly%PercentBelowCriticalValue_R4 (i, j, k) = &
! (Me%Monthly%PercentBelowCriticalValue_R4 (i, j, k) * &
! Me%Monthly%RunPeriod + DT) / (Me%Monthly%RunPeriod + DT)
! else
! Me%Monthly%PercentBelowCriticalValue_R4 (i, j, k) = &
! (Me%Monthly%PercentBelowCriticalValue_R4 (i, j, k) * &
! Me%Monthly%RunPeriod + 0.) / (Me%Monthly%RunPeriod + DT)
! endif
! endif
!
! endif
! enddo
! enddo
! enddo
!
! !Verifies if the present time is a new output
! call ExtractDate (Time1=Me%ExternalVar%Now, Month = PresentMonth)
! call ExtractDate (Time1=Me%Monthly%LastCalculation, Month = OldMonth)
! if (int(PresentMonth) /= int(OldMonth)) then
! call WriteValuesToFileHDF5_R4 (.false., .false., .true., .false., .false.)
! Me%Monthly%Minimum_R4 = Value_R4
! Me%Monthly%Maximum_R4 = Value_R4
! Me%Monthly%Average_R4 = Value_R4
! Me%Monthly%SquareAverage_R4 = Value_R4
! Me%Monthly%StandardDeviation_R4 = Value_R4
!
! if (Me%Accumulated) Me%Monthly%Accumulated_R4 = 0.0
!
! if (Me%GeomMean) then !Geometric Average to calculate
! Me%Monthly%GeomAverage_R4 = Value_R4
! Me%Monthly%SquareGeomAverage_R4 = Value_R4
! Me%Monthly%GeomStandardDeviation_R4 = Value_R4
! endif
!
! !guillaume juan
! if (Me%Critical) then
! do k = KLB, KUB
! do j = JLB, JUB
! do i = ILB, IUB
! if ((WaterPoints3D(i, j, k) == WaterPoint) .and. (Value_R4 (i, j, k) > FillValueReal / 2.)) then
! if(Value_R4(i, j, k) < Me%CriticalValue) then
! Me%Monthly%PercentBelowCriticalValue_R4 (i, j, k) = 1.
! else
! Me%Monthly%PercentBelowCriticalValue_R4 (i, j, k) = 0.
! endif
! endif
! enddo
! enddo
! enddo
! endif
!
! Me%Monthly%RunPeriod = 0.
! else
! Me%Monthly%RunPeriod = Me%Monthly%RunPeriod + DT
! endif
!
! !Updates Time
! Me%Monthly%LastCalculation = Me%ExternalVar%Now
!
! endif cd1
!
!
! end subroutine ModifyMonthlyStatistic_R4
!--------------------------------------------------------------------------
subroutine ModifySpecificHourStatistic (Value, WaterPoints3D, KLB, KUB)
!Arguments-------------------------------------------------------------
real, dimension(:, :, :), pointer :: Value
integer, dimension(:, :, :), pointer :: WaterPoints3D
integer :: KLB, KUB
!Local-----------------------------------------------------------------
integer :: ILB, IUB, i
integer :: JLB, JUB, j
integer :: k
real :: DT, DX, AuxValue
real :: PresentHour
!Shorten
ILB = Me%ExternalVar%WorkSize%ILB
IUB = Me%ExternalVar%WorkSize%IUB
JLB = Me%ExternalVar%WorkSize%JLB
JUB = Me%ExternalVar%WorkSize%JUB
call ExtractDate (Time1=Me%ExternalVar%Now, Hour = PresentHour)
if1: if (int(PresentHour) == int(Me%SpecificHourValue)) then
!Time since last calculation
DT = Me%ExternalVar%Now - Me%SpecificHour%LastCalculation
cd1: if (DT>0) then
!Loops
do k = KLB, KUB
do j = JLB, JUB
do i = ILB, IUB
if (WaterPoints3D(i, j, k) == WaterPoint) then
!Minimum Value
if (Value (i, j, k) < Me%SpecificHour%Minimum (i, j, k)) &
Me%SpecificHour%Minimum (i, j, k) = Value (i, j, k)
!Maximum Value
if (Value (i, j, k) > Me%SpecificHour%Maximum (i, j, k)) &
Me%SpecificHour%Maximum (i, j, k) = Value (i, j, k)
!Average
Me%SpecificHour%Average (i, j, k) = &
(Me%SpecificHour%Average (i, j, k) * &
Me%SpecificHour%RunPeriod + &
Value (i, j, k) * DT) / (Me%SpecificHour%RunPeriod + DT)
!Square Average
Me%SpecificHour%SquareAverage (i, j, k) = &
(Me%SpecificHour%SquareAverage (i, j, k) * &
Me%SpecificHour%RunPeriod + &
Value (i, j, k)**2 * DT) / (Me%SpecificHour%RunPeriod + DT)
!Standard deviation
DX = Me%SpecificHour%SquareAverage (i, j, k) - &
Me%SpecificHour%Average (i, j, k) ** 2
DX = abs(DX)
Me%SpecificHour%StandardDeviation(i, j, k) = sqrt(DX)
!Accumulated Values
if (Me%Accumulated) &
Me%SpecificHour%Accumulated (i,j,k) = Me%SpecificHour%Accumulated (i,j,k) &
+ Value (i, j, k)
if (Me%GeomMean) then !Geometric Average to calculate
!Geometric Average
AuxValue = Value (i, j, k)
if (AuxValue == 0.0) then
AuxValue = 1.0
elseif (AuxValue .lt. 0.0) then
write(*,*) 'Negative valued property.'
write(*,*) 'Geometric Average cannot be calculated.'
stop 'ModifySpecificHourStatistic - ModuleStatistic - ERR01'
endif
if (Me%SpecificHour%GeomAverage (i, j, k) == 0.0) then
Me%SpecificHour%GeomAverage (i, j, k) = 1.0
endif
Me%SpecificHour%GeomAverage (i, j, k) = 10** &
((LOG10(Me%SpecificHour%GeomAverage (i, j, k)) * &
Me%SpecificHour%RunPeriod + &
LOG10(AuxValue) * DT) / (Me%SpecificHour%RunPeriod + DT))
!Squared Geometric Average
Me%SpecificHour%SquareGeomAverage (i, j, k) = &
(Me%SpecificHour%SquareGeomAverage (i, j, k) * &
Me%SpecificHour%RunPeriod + &
AuxValue**2 * DT) / (Me%SpecificHour%RunPeriod + DT)
!Geometric Standard Deviation
DX = Me%SpecificHour%SquareGeomAverage (i, j, k) - &
Me%SpecificHour%GeomAverage (i, j, k) ** 2
DX = abs(DX)
Me%SpecificHour%GeomStandardDeviation(i, j, k) = sqrt(DX)
endif
!guillaume juan
if (Me%Critical .and. (Value (i, j, k) > FillValueReal / 2.)) then
if (Value (i, j, k) < Me%CriticalValue) then
Me%SpecificHour%PercentBelowCriticalValue (i, j, k) = &
(Me%SpecificHour%PercentBelowCriticalValue (i, j, k) * &
Me%SpecificHour%RunPeriod + DT) / (Me%SpecificHour%RunPeriod + DT)
else
Me%SpecificHour%PercentBelowCriticalValue (i, j, k) = &
(Me%SpecificHour%PercentBelowCriticalValue (i, j, k) * &
Me%SpecificHour%RunPeriod + 0.) / (Me%SpecificHour%RunPeriod + DT)
endif
endif
endif
enddo
enddo
enddo
!Updates Time
Me%SpecificHour%RunPeriod = Me%SpecificHour%RunPeriod + DT
Me%SpecificHour%LastCalculation = Me%ExternalVar%Now
endif cd1
endif if1
end subroutine ModifySpecificHourStatistic
!--------------------------------------------------------------------------
! subroutine ModifySpecificHourStatistic_R4 (Value_R4, WaterPoints3D, KLB, KUB)
!
! !Arguments-------------------------------------------------------------
! real(4), dimension(:, :, :), pointer :: Value_R4
! integer, dimension(:, :, :), pointer :: WaterPoints3D
! integer :: KLB, KUB
!
! !Local-----------------------------------------------------------------
! integer :: ILB, IUB, i
! integer :: JLB, JUB, j
! integer :: k
! real :: DT
! real(4) :: DX, AuxValue
! real :: PresentHour
!
! !Shorten
! ILB = Me%ExternalVar%WorkSize%ILB
! IUB = Me%ExternalVar%WorkSize%IUB
! JLB = Me%ExternalVar%WorkSize%JLB
! JUB = Me%ExternalVar%WorkSize%JUB
!
! call ExtractDate (Time1=Me%ExternalVar%Now, Hour = PresentHour)
!
!if1: if (int(PresentHour) == int(Me%SpecificHourValue)) then
!
! !Time since last calculation
! DT = Me%ExternalVar%Now - Me%SpecificHour%LastCalculation
!
!cd1: if (DT>0) then
!
! !Loops
! do k = KLB, KUB
! do j = JLB, JUB
! do i = ILB, IUB
! if (WaterPoints3D(i, j, k) == WaterPoint) then
!
! !Minimum Value
! if (Value_R4 (i, j, k) < Me%SpecificHour%Minimum_R4 (i, j, k)) &
! Me%SpecificHour%Minimum_R4 (i, j, k) = Value_R4 (i, j, k)
!
! !Maximum Value
! if (Value_R4 (i, j, k) > Me%SpecificHour%Maximum_R4 (i, j, k)) &
! Me%SpecificHour%Maximum_R4 (i, j, k) = Value_R4 (i, j, k)
!
! !Average
! Me%SpecificHour%Average_R4 (i, j, k) = &
! (Me%SpecificHour%Average_R4 (i, j, k) * &
! Me%SpecificHour%RunPeriod + &
! Value_R4 (i, j, k) * DT) / (Me%SpecificHour%RunPeriod + DT)
!
! !Square Average
! Me%SpecificHour%SquareAverage_R4 (i, j, k) = &
! (Me%SpecificHour%SquareAverage_R4 (i, j, k) * &
! Me%SpecificHour%RunPeriod + &
! Value_R4 (i, j, k)**2 * DT) / (Me%SpecificHour%RunPeriod + DT)
!
! !Standard deviation
! DX = Me%SpecificHour%SquareAverage_R4 (i, j, k) - &
! Me%SpecificHour%Average_R4 (i, j, k) ** 2
!
! DX = abs(DX)
!
!
! Me%SpecificHour%StandardDeviation_R4(i, j, k) = sqrt(DX)
!
! !Accumulated Values
! if (Me%Accumulated) &
! Me%SpecificHour%Accumulated_R4 (i,j,k) = Me%SpecificHour%Accumulated_R4 (i,j,k) &
! + Value_R4 (i, j, k)
!
! if (Me%GeomMean) then !Geometric Average to calculate
! !Geometric Average
! AuxValue = Value_R4 (i, j, k)
! if (AuxValue == 0.0) then
! AuxValue = 1.0
! elseif (AuxValue .lt. 0.0) then
! write(*,*) 'Negative valued property.'
! write(*,*) 'Geometric Average cannot be calculated.'
! stop 'ModifySpecificHourStatistic_R4 - ModuleStatistic - ERR01'
! endif
!
! if (Me%SpecificHour%GeomAverage_R4 (i, j, k) == 0.0) then
! Me%SpecificHour%GeomAverage_R4 (i, j, k) = 1.0
! endif
!
! Me%SpecificHour%GeomAverage_R4 (i, j, k) = 10** &
! ((LOG10(Me%SpecificHour%GeomAverage_R4 (i, j, k)) * &
! Me%SpecificHour%RunPeriod + &
! LOG10(AuxValue) * DT) / (Me%SpecificHour%RunPeriod + DT))
!
! !Squared Geometric Average
! Me%SpecificHour%SquareGeomAverage_R4 (i, j, k) = &
! (Me%SpecificHour%SquareGeomAverage_R4 (i, j, k) * &
! Me%SpecificHour%RunPeriod + &
! AuxValue**2 * DT) / (Me%SpecificHour%RunPeriod + DT)
!
! !Geometric Standard Deviation
! DX = Me%SpecificHour%SquareGeomAverage_R4 (i, j, k) - &
! Me%SpecificHour%GeomAverage_R4 (i, j, k) ** 2
!
! DX = abs(DX)
!
! Me%SpecificHour%GeomStandardDeviation_R4(i, j, k) = sqrt(DX)
!
! endif
!
! !guillaume juan
! if (Me%Critical .and. (Value_R4 (i, j, k) > FillValueReal / 2.)) then
! if (Value_R4 (i, j, k) < Me%CriticalValue) then
! Me%SpecificHour%PercentBelowCriticalValue_R4 (i, j, k) = &
! (Me%SpecificHour%PercentBelowCriticalValue_R4 (i, j, k) * &
! Me%SpecificHour%RunPeriod + DT) / (Me%SpecificHour%RunPeriod + DT)
! else
! Me%SpecificHour%PercentBelowCriticalValue_R4 (i, j, k) = &
! (Me%SpecificHour%PercentBelowCriticalValue_R4 (i, j, k) * &
! Me%SpecificHour%RunPeriod + 0.) / (Me%SpecificHour%RunPeriod + DT)
! endif
! endif
!
! endif
! enddo
! enddo
! enddo
!
!
! !Updates Time
! Me%SpecificHour%RunPeriod = Me%SpecificHour%RunPeriod + DT
! Me%SpecificHour%LastCalculation = Me%ExternalVar%Now
!
! endif cd1
!
! endif if1
!
! end subroutine ModifySpecificHourStatistic_R4
!--------------------------------------------------------------------------
subroutine ModifyClassification (Value, WaterPoints3D, KLB, KUB)
!Arguments-------------------------------------------------------------
real, dimension(:, :, :), pointer :: Value
integer, dimension(:, :, :), pointer :: WaterPoints3D
integer :: KLB, KUB
!Local-----------------------------------------------------------------
integer :: ILB, IUB, i
integer :: JLB, JUB, j
integer :: k
integer :: iClass
real :: DT, Aux
!Shorten
ILB = Me%ExternalVar%WorkSize%ILB
IUB = Me%ExternalVar%WorkSize%IUB
JLB = Me%ExternalVar%WorkSize%JLB
JUB = Me%ExternalVar%WorkSize%JUB
!Time since last calculation
DT = Me%ExternalVar%Now - Me%Classification%LastCalculation
cd1: if (DT>0) then
!Loops
do k = KLB, KUB
do j = JLB, JUB
do i = ILB, IUB
if (WaterPoints3D(i, j, k) == WaterPoint) then
doClass: do iClass = 1, Me%Classification%nClasses
if (Value(i, j, k) >= Me%Classification%Classes(iClass, 1) .and. &
Value(i, j, k) < Me%Classification%Classes(iClass, 2)) then
Aux = DT
else
Aux = 0
endif
Me%Classification%Frequency (i, j, k, iClass) = &
(Me%Classification%Frequency(i, j, k, iClass) * &
Me%Classification%RunPeriod + Aux ) / &
(Me%Classification%RunPeriod + DT)
enddo doClass
endif
enddo
enddo
enddo
Me%Classification%RunPeriod = Me%Classification%RunPeriod + DT
Me%Classification%LastCalculation = Me%ExternalVar%Now
endif cd1
end subroutine ModifyClassification
!--------------------------------------------------------------------------
! subroutine ModifyClassification_R4 (Value_R4, WaterPoints3D, KLB, KUB)
!
! !Arguments-------------------------------------------------------------
! real(4), dimension(:, :, :), pointer :: Value_R4
! integer, dimension(:, :, :), pointer :: WaterPoints3D
! integer :: KLB, KUB
!
!
! !Local-----------------------------------------------------------------
! integer :: ILB, IUB, i
! integer :: JLB, JUB, j
! integer :: k
! integer :: iClass
! real :: DT
! real(4) :: Aux
!
! !Shorten
! ILB = Me%ExternalVar%WorkSize%ILB
! IUB = Me%ExternalVar%WorkSize%IUB
! JLB = Me%ExternalVar%WorkSize%JLB
! JUB = Me%ExternalVar%WorkSize%JUB
!
!
! !Time since last calculation
! DT = Me%ExternalVar%Now - Me%Classification%LastCalculation
!
!cd1: if (DT>0) then
!
! !Loops
! do k = KLB, KUB
! do j = JLB, JUB
! do i = ILB, IUB
! if (WaterPoints3D(i, j, k) == WaterPoint) then
!doClass: do iClass = 1, Me%Classification%nClasses
! if (Value_R4(i, j, k) >= Me%Classification%Classes(iClass, 1) .and. &
! Value_R4(i, j, k) < Me%Classification%Classes(iClass, 2)) then
! Aux = DT
! else
! Aux = 0
! endif
!
! Me%Classification%Frequency_R4 (i, j, k, iClass) = &
! (Me%Classification%Frequency_R4(i, j, k, iClass) * &
! Me%Classification%RunPeriod + Aux ) / &
! (Me%Classification%RunPeriod + DT)
! enddo doClass
! endif
! enddo
! enddo
! enddo
!
!
! Me%Classification%RunPeriod = Me%Classification%RunPeriod + DT
!
! Me%Classification%LastCalculation = Me%ExternalVar%Now
!
! endif cd1
!
! end subroutine ModifyClassification_R4
!--------------------------------------------------------------------------
subroutine ModifyGlobalStatistic2D (Value2D, WaterPoints2D)
!Arguments-------------------------------------------------------------
real, dimension(:,:), pointer :: Value2D
integer, dimension(:,:), pointer :: WaterPoints2D
!Local-----------------------------------------------------------------
integer :: ILB, IUB, i
integer :: JLB, JUB, j
real :: DT, DX, AuxValue
!Shorten
ILB = Me%ExternalVar%WorkSize%ILB
IUB = Me%ExternalVar%WorkSize%IUB
JLB = Me%ExternalVar%WorkSize%JLB
JUB = Me%ExternalVar%WorkSize%JUB
!Time since last calculation
DT = Me%ExternalVar%Now - Me%Global%LastCalculation
cd1: if (DT>0) then
!Loops
do j = JLB, JUB
do i = ILB, IUB
if (WaterPoints2D(i, j) == WaterPoint) then
!Minimum Value2D
if (Value2D (i, j) < Me%Global%Minimum2D (i, j)) &
Me%Global%Minimum2D (i, j) = Value2D (i, j)
!Maximum Value2D
if (Value2D (i, j) > Me%Global%Maximum2D (i, j)) &
Me%Global%Maximum2D (i, j) = Value2D (i, j)
!Average
Me%Global%Average2D (i, j) = &
(Me%Global%Average2D (i, j) * &
Me%Global%RunPeriod + &
Value2D (i, j) * DT) / (Me%Global%RunPeriod + DT)
!Square Average
Me%Global%SquareAverage2D (i, j) = &
(Me%Global%SquareAverage2D (i, j) * &
Me%Global%RunPeriod + &
Value2D (i, j)**2 * DT) / (Me%Global%RunPeriod + DT)
!Standard deviation
DX = Me%Global%SquareAverage2D (i, j) - &
Me%Global%Average2D (i, j) ** 2
DX = abs(DX)
Me%Global%StandardDeviation2D(i, j) = sqrt(DX)
!Accumulated Values
if (Me%Accumulated) &
Me%Global%Accumulated2D (i,j) = Me%Global%Accumulated2D (i,j) &
+ Value2D (i, j)
if (Me%GeomMean) then !Geometric Average to calculate
!Geometric Average
AuxValue = Value2D (i, j)
if (AuxValue == 0.0) then
AuxValue = 1.0
elseif (AuxValue .lt. 0.0) then
write(*,*) 'Negative valued property.'
write(*,*) 'Geometric Average cannot be calculated.'
stop 'ModifyGlobalStatistic2D - ModuleStatistic - ERR01'
endif
Me%Global%GeomAverage2D (i, j) = 10** &
((LOG10(Me%Global%GeomAverage2D(i, j)) * &
Me%Global%RunPeriod + &
LOG10(AuxValue) * DT) / (Me%Global%RunPeriod + DT))
!Squared Geometric Average
Me%Global%SquareGeomAverage2D (i, j) = &
(Me%Global%SquareGeomAverage2D (i, j) * &
Me%Global%RunPeriod + &
AuxValue**2 * DT) / (Me%Global%RunPeriod + DT)
!Geometric Standard Deviation
DX = Me%Global%SquareGeomAverage2D (i, j) - &
Me%Global%GeomAverage2D (i, j) ** 2
DX = abs(DX)
Me%Global%GeomStandardDeviation2D(i, j) = sqrt(DX)
endif
!guillaume juan
if (Me%Critical .and. (Value2D (i, j) > FillValueReal / 2.)) then
if (Value2D (i, j) < Me%CriticalValue) then
Me%Global%PercentBelowCriticalValue2D (i, j) = &
(Me%Global%PercentBelowCriticalValue2D (i, j) * &
Me%Global%RunPeriod + DT) / (Me%Global%RunPeriod + DT)
else
Me%Global%PercentBelowCriticalValue2D (i, j) = &
(Me%Global%PercentBelowCriticalValue2D (i, j) * &
Me%Global%RunPeriod + 0.) / (Me%Global%RunPeriod + DT)
endif
endif
endif
enddo
enddo
!Updates Time
Me%Global%RunPeriod = Me%Global%RunPeriod + DT
Me%Global%LastCalculation = Me%ExternalVar%Now
endif cd1
end subroutine ModifyGlobalStatistic2D
!--------------------------------------------------------------------------
subroutine ModifyDailyStatistic2D (Value2D, WaterPoints2D)
!Arguments-------------------------------------------------------------
real, dimension(:,:), pointer :: Value2D
integer, dimension(:,:), pointer :: WaterPoints2D
!Local-----------------------------------------------------------------
integer :: ILB, IUB, i
integer :: JLB, JUB, j
real :: DT, DX, AuxValue
! real :: OldDay, PresentDay
!Shorten
ILB = Me%ExternalVar%WorkSize%ILB
IUB = Me%ExternalVar%WorkSize%IUB
JLB = Me%ExternalVar%WorkSize%JLB
JUB = Me%ExternalVar%WorkSize%JUB
!Time since last calculation
DT = Me%ExternalVar%Now - Me%Daily%LastCalculation
cd1: if (DT>0) then
!Loops
do j = JLB, JUB
do i = ILB, IUB
if (WaterPoints2D(i, j) == WaterPoint) then
!Minimum Value2D
if (Value2D (i, j) < Me%Daily%Minimum2D (i, j)) &
Me%Daily%Minimum2D (i, j) = Value2D (i, j)
!Maximum Value2D
if (Value2D (i, j) > Me%Daily%Maximum2D (i, j)) &
Me%Daily%Maximum2D (i, j) = Value2D (i, j)
!Average
Me%Daily%Average2D (i, j) = &
(Me%Daily%Average2D (i, j) * &
Me%Daily%RunPeriod + &
Value2D (i, j) * DT) / (Me%Daily%RunPeriod + DT)
!Square Average
Me%Daily%SquareAverage2D (i, j) = &
(Me%Daily%SquareAverage2D (i, j) * &
Me%Daily%RunPeriod + &
Value2D (i, j)**2 * DT) / (Me%Daily%RunPeriod + DT)
!Standard deviation
DX = Me%Daily%SquareAverage2D (i, j) - &
Me%Daily%Average2D (i, j) ** 2
DX = abs(DX)
Me%Daily%StandardDeviation2D(i, j) = sqrt(DX)
!Accumulated Values
if (Me%Accumulated) &
Me%Daily%Accumulated2D (i,j) = Me%Daily%Accumulated2D (i,j) &
+ Value2D (i, j)
if (Me%GeomMean) then !Geometric Average to calculate
!Geometric Average
AuxValue = Value2D (i, j)
if (AuxValue == 0.0) then
AuxValue = 1.0
elseif (AuxValue .lt. 0.0) then
write(*,*) 'Negative valued property.'
write(*,*) 'Geometric Average cannot be calculated.'
stop 'ModifyDailyStatistic2D - ModuleStatistic - ERR01'
endif
if (Me%Daily%GeomAverage2D (i, j) == 0.0) then
Me%Daily%GeomAverage2D (i, j) = 1.0
endif
Me%Daily%GeomAverage2D (i, j) = 10** &
((LOG10(Me%Daily%GeomAverage2D (i, j)) * &
Me%Daily%RunPeriod + &
LOG10(AuxValue) * DT) / (Me%Daily%RunPeriod + DT))
!Squared Geometric Average
Me%Daily%SquareGeomAverage2D (i, j) = &
(Me%Daily%SquareGeomAverage2D (i, j) * &
Me%Daily%RunPeriod + &
AuxValue**2 * DT) / (Me%Daily%RunPeriod + DT)
!Geometric Standard Deviation
DX = Me%Daily%SquareGeomAverage2D (i, j) - &
Me%Daily%GeomAverage2D (i, j) ** 2
DX = abs(DX)
Me%Daily%GeomStandardDeviation2D(i, j) = sqrt(DX)
endif
!guillaume juan
if (Me%Critical .and. (Value2D (i, j) > FillValueReal / 2.)) then
if (Value2D (i, j) < Me%CriticalValue) then
Me%Daily%PercentBelowCriticalValue2D (i, j) = &
(Me%Daily%PercentBelowCriticalValue2D (i, j) * &
Me%Daily%RunPeriod + DT) / (Me%Daily%RunPeriod + DT)
else
Me%Daily%PercentBelowCriticalValue2D (i, j) = &
(Me%Daily%PercentBelowCriticalValue2D (i, j) * &
Me%Daily%RunPeriod + 0.) / (Me%Daily%RunPeriod + DT)
endif
endif
endif
enddo
enddo
!Verifies if the present time is a new output
! call ExtractDate (Me%ExternalVar%Now, Day = PresentDay)
! call ExtractDate (Me%Daily%LastCalculation, Day = OldDay)
! if (int(PresentDay) /= int(OldDay)) then
if (Me%ExternalVar%Now .GT. Me%Daily%NextOutputTime) then
call WriteValuesToFileHDF5 (.false., .true., .false., .false., .false.)
Me%Daily%Minimum2D = Value2D
Me%Daily%Maximum2D = Value2D
Me%Daily%Average2D = Value2D
Me%Daily%SquareAverage2D = Value2D
Me%Daily%StandardDeviation2D = Value2D
if (Me%Accumulated) Me%Daily%Accumulated2D = 0.0
if (Me%GeomMean) then !Geometric Average to calculate
Me%Daily%GeomAverage2D = Value2D
Me%Daily%SquareGeomAverage2D = Value2D
Me%Daily%GeomStandardDeviation2D = Value2D
endif
!guillaume juan
if (Me%Critical) then
do j = JLB, JUB
do i = ILB, IUB
if ( (WaterPoints2D(i, j) == WaterPoint) .and. (Value2D (i, j) > FillValueReal / 2.)) then
if(Value2D(i, j) < Me%CriticalValue) then
Me%Daily%PercentBelowCriticalValue2D (i, j) = 1.
else
Me%Daily%PercentBelowCriticalValue2D (i, j) = 0.
endif
endif
enddo
enddo
endif
Me%Daily%RunPeriod = 0.
Me%Daily%NextOutputTime = Me%Daily%NextOutputTime + 24*3600.
else
Me%Daily%RunPeriod = Me%Daily%RunPeriod + DT
endif
!Updates Time
Me%Daily%LastCalculation = Me%ExternalVar%Now
endif cd1
end subroutine ModifyDailyStatistic2D
!--------------------------------------------------------------------------
subroutine ModifyMonthlyStatistic2D (Value2D, WaterPoints2D)
!Arguments-------------------------------------------------------------
real, dimension(:,:), pointer :: Value2D
integer, dimension(:,:), pointer :: WaterPoints2D
!Local-----------------------------------------------------------------
integer :: ILB, IUB, i
integer :: JLB, JUB, j
real :: DT, DX, AuxValue
real :: OldMonth, PresentMonth
!Shorten
ILB = Me%ExternalVar%WorkSize%ILB
IUB = Me%ExternalVar%WorkSize%IUB
JLB = Me%ExternalVar%WorkSize%JLB
JUB = Me%ExternalVar%WorkSize%JUB
!Time since last calculation
DT = Me%ExternalVar%Now - Me%Monthly%LastCalculation
cd1: if (DT>0) then
!Loops
do j = JLB, JUB
do i = ILB, IUB
if (WaterPoints2D(i, j) == WaterPoint) then
!Minimum Value2D
if (Value2D (i, j) < Me%Monthly%Minimum2D (i, j)) &
Me%Monthly%Minimum2D (i, j) = Value2D (i, j)
!Maximum Value2D
if (Value2D (i, j) > Me%Monthly%Maximum2D (i, j)) &
Me%Monthly%Maximum2D (i, j) = Value2D (i, j)
!Average
Me%Monthly%Average2D (i, j) = &
(Me%Monthly%Average2D (i, j) * &
Me%Monthly%RunPeriod + &
Value2D (i, j) * DT) / (Me%Monthly%RunPeriod + DT)
!Square Average
Me%Monthly%SquareAverage2D (i, j) = &
(Me%Monthly%SquareAverage2D (i, j) * &
Me%Monthly%RunPeriod + &
Value2D (i, j)**2 * DT) / (Me%Monthly%RunPeriod + DT)
!Standard deviation
DX = Me%Monthly%SquareAverage2D (i, j) - &
Me%Monthly%Average2D (i, j) ** 2
DX = abs(DX)
Me%Monthly%StandardDeviation2D(i, j) = sqrt(DX)
!Accumulated Values
if (Me%Accumulated) &
Me%Monthly%Accumulated2D (i,j) = Me%Monthly%Accumulated2D (i,j) &
+ Value2D (i, j)
if (Me%GeomMean) then !Geometric Average to calculate
!Geometric Average
AuxValue = Value2D (i, j)
if (AuxValue == 0.0) then
AuxValue = 1.0
elseif (AuxValue .lt. 0.0) then
write(*,*) 'Negative valued property.'
write(*,*) 'Geometric Average cannot be calculated.'
stop 'ModifyMonthlyStatistic2D - ModuleStatistic - ERR01'
endif
if (Me%Monthly%GeomAverage2D (i, j) == 0.0) then
Me%Monthly%GeomAverage2D (i, j) = 1.0
endif
Me%Monthly%GeomAverage2D (i, j) = 10** &
((LOG10(Me%Monthly%GeomAverage2D (i, j)) * &
Me%Monthly%RunPeriod + &
LOG10(AuxValue) * DT) / (Me%Monthly%RunPeriod + DT))
!Squared Geometric Average
Me%Monthly%SquareGeomAverage2D (i, j) = &
(Me%Monthly%SquareGeomAverage2D (i, j) * &
Me%Monthly%RunPeriod + &
AuxValue**2 * DT) / (Me%Monthly%RunPeriod + DT)
!Geometric Standard Deviation
DX = Me%Monthly%SquareGeomAverage2D (i, j) - &
Me%Monthly%GeomAverage2D (i, j) ** 2
DX = abs(DX)
Me%Monthly%GeomStandardDeviation2D(i, j) = sqrt(DX)
endif
!guillaume juan
if (Me%Critical .and. (Value2D (i, j) > FillValueReal / 2.)) then
if (Value2D (i, j) < Me%CriticalValue) then
Me%Monthly%PercentBelowCriticalValue2D (i, j) = &
(Me%Monthly%PercentBelowCriticalValue2D (i, j) * &
Me%Monthly%RunPeriod + DT) / (Me%Monthly%RunPeriod + DT)
else
Me%Monthly%PercentBelowCriticalValue2D (i, j) = &
(Me%Monthly%PercentBelowCriticalValue2D (i, j) * &
Me%Monthly%RunPeriod + 0.) / (Me%Monthly%RunPeriod + DT)
endif
endif
endif
enddo
enddo
!Verifies if the present time is a new output
call ExtractDate (Time1=Me%ExternalVar%Now, Month = PresentMonth)
call ExtractDate (Time1=Me%Monthly%LastCalculation, Month = OldMonth)
if (int(PresentMonth) /= int(OldMonth)) then
call WriteValuesToFileHDF5 (.false., .false., .true., .false., .false.)
Me%Monthly%Minimum2D = Value2D
Me%Monthly%Maximum2D = Value2D
Me%Monthly%Average2D = Value2D
Me%Monthly%SquareAverage2D = Value2D
Me%Monthly%StandardDeviation2D = Value2D
if (Me%Accumulated) Me%Monthly%Accumulated2D = 0.0
if (Me%GeomMean) then !Geometric Average to calculate
Me%Monthly%GeomAverage2D = Value2D
Me%Monthly%SquareGeomAverage2D = Value2D
Me%Monthly%GeomStandardDeviation2D = Value2D
endif
!guillaume juan
if (Me%Critical) then
do j = JLB, JUB
do i = ILB, IUB
if ( (WaterPoints2D(i, j) == WaterPoint) .and. (Value2D (i, j) > FillValueReal / 2.)) then
if(Value2D(i, j) < Me%CriticalValue) then
Me%Monthly%PercentBelowCriticalValue2D (i, j) = 1.
else
Me%Monthly%PercentBelowCriticalValue2D (i, j) = 0.
endif
endif
enddo
enddo
endif
Me%Monthly%RunPeriod = 0.
else
Me%Monthly%RunPeriod = Me%Monthly%RunPeriod + DT
endif
!Updates Time
Me%Monthly%LastCalculation = Me%ExternalVar%Now
endif cd1
end subroutine ModifyMonthlyStatistic2D
!--------------------------------------------------------------------------
subroutine ModifySpecificHourStatistic2D (Value2D, WaterPoints2D)
!Arguments-------------------------------------------------------------
real, dimension(:,:), pointer :: Value2D
integer, dimension(:,:), pointer :: WaterPoints2D
!Local-----------------------------------------------------------------
integer :: ILB, IUB, i
integer :: JLB, JUB, j
real :: DT, DX, AuxValue
real :: PresentHour
!Shorten
ILB = Me%ExternalVar%WorkSize%ILB
IUB = Me%ExternalVar%WorkSize%IUB
JLB = Me%ExternalVar%WorkSize%JLB
JUB = Me%ExternalVar%WorkSize%JUB
call ExtractDate (Time1=Me%ExternalVar%Now, Hour = PresentHour)
if1: if (int(PresentHour) == int(Me%SpecificHourValue)) then
!Time since last calculation
DT = Me%ExternalVar%Now - Me%SpecificHour%LastCalculation
cd1: if (DT>0) then
!Loops
do j = JLB, JUB
do i = ILB, IUB
if (WaterPoints2D(i, j) == WaterPoint) then
!Minimum Value2D
if (Value2D (i, j) < Me%SpecificHour%Minimum2D (i, j)) &
Me%SpecificHour%Minimum2D (i, j) = Value2D (i, j)
!Maximum Value2D
if (Value2D (i, j) > Me%SpecificHour%Maximum2D (i, j)) &
Me%SpecificHour%Maximum2D (i, j) = Value2D (i, j)
!Average
Me%SpecificHour%Average2D (i, j) = &
(Me%SpecificHour%Average2D (i, j) * &
Me%SpecificHour%RunPeriod + &
Value2D (i, j) * DT) / (Me%SpecificHour%RunPeriod + DT)
!Square Average
Me%SpecificHour%SquareAverage2D (i, j) = &
(Me%SpecificHour%SquareAverage2D (i, j) * &
Me%SpecificHour%RunPeriod + &
Value2D (i, j)**2 * DT) / (Me%SpecificHour%RunPeriod + DT)
!Standard deviation
DX = Me%SpecificHour%SquareAverage2D (i, j) - &
Me%SpecificHour%Average2D (i, j) ** 2
DX = abs(DX)
Me%SpecificHour%StandardDeviation2D(i, j) = sqrt(DX)
!Accumulated Values
if (Me%Accumulated) &
Me%SpecificHour%Accumulated2D (i,j) = Me%SpecificHour%Accumulated2D (i,j) &
+ Value2D (i, j)
if (Me%GeomMean) then !Geometric Average to calculate
!Geometric Average
AuxValue = Value2D (i, j)
if (AuxValue == 0.0) then
AuxValue = 1.0
elseif (AuxValue .lt. 0.0) then
write(*,*) 'Negative valued property.'
write(*,*) 'Geometric Average cannot be calculated.'
stop 'ModifySpecificHourStatistic2D - ModuleStatistic - ERR01'
endif
Me%SpecificHour%GeomAverage2D (i, j) = 10** &
((LOG10(Me%SpecificHour%GeomAverage2D(i, j)) * &
Me%SpecificHour%RunPeriod + &
LOG10(AuxValue) * DT) / (Me%SpecificHour%RunPeriod + DT))
!Squared Geometric Average
Me%SpecificHour%SquareGeomAverage2D (i, j) = &
(Me%SpecificHour%SquareGeomAverage2D (i, j) * &
Me%SpecificHour%RunPeriod + &
AuxValue**2 * DT) / (Me%SpecificHour%RunPeriod + DT)
!Geometric Standard Deviation
DX = Me%SpecificHour%SquareGeomAverage2D (i, j) - &
Me%SpecificHour%GeomAverage2D (i, j) ** 2
DX = abs(DX)
Me%SpecificHour%GeomStandardDeviation2D(i, j) = sqrt(DX)
endif
!guillaume juan
if (Me%Critical .and. (Value2D (i, j) > FillValueReal / 2.)) then
if (Value2D (i, j) < Me%CriticalValue) then
Me%SpecificHour%PercentBelowCriticalValue2D (i, j) = &
(Me%SpecificHour%PercentBelowCriticalValue2D (i, j) * &
Me%SpecificHour%RunPeriod + DT) / (Me%SpecificHour%RunPeriod + DT)
else
Me%SpecificHour%PercentBelowCriticalValue2D (i, j) = &
(Me%SpecificHour%PercentBelowCriticalValue2D (i, j) * &
Me%SpecificHour%RunPeriod + 0.) / (Me%SpecificHour%RunPeriod + DT)
endif
endif
endif
enddo
enddo
!Updates Time
Me%SpecificHour%RunPeriod = Me%SpecificHour%RunPeriod + DT
Me%SpecificHour%LastCalculation = Me%ExternalVar%Now
endif cd1
endif if1
end subroutine ModifySpecificHourStatistic2D
!--------------------------------------------------------------------------
subroutine ModifyClassification2D (Value2D, WaterPoints2D)
!Arguments-------------------------------------------------------------
real, dimension(:,:), pointer :: Value2D
integer, dimension(:,:), pointer :: WaterPoints2D
!Local--------------