LCOV - code coverage report
Current view: top level - src - mo_errormeasures.f90 (source / functions) Hit Total Coverage
Test: forces coverage Lines: 1129 1195 94.5 %
Date: 2024-03-13 19:03:28 Functions: 66 66 100.0 %

          Line data    Source code
       1             : !> \file mo_errormeasures.f90
       2             : !> \brief \copybrief mo_errormeasures
       3             : !> \details \copydetails mo_errormeasures
       4             : 
       5             : !> \brief Calculation of error measures.
       6             : !> \details This module contains routines for the masked calculation of
       7             : !! error measures like MSE, RMSE, BIAS, SSE, NSE, ...
       8             : !> \note all except variance and standard deviation are population and not sample moments,
       9             : !!       i.e. they are normally divided by n and not (n-1)
      10             : !> \authors Mathias Zink
      11             : !> \date Aug 2012
      12             : !> \copyright Copyright 2005-\today, the CHS Developers, Sabine Attinger: All rights reserved.
      13             : !! FORCES is released under the LGPLv3+ license \license_note
      14             : MODULE mo_errormeasures
      15             : 
      16             :   USE mo_kind, ONLY : i4, sp, dp
      17             : 
      18             :   IMPLICIT NONE
      19             : 
      20             :   PUBLIC :: BIAS                         ! bias
      21             :   PUBLIC :: KGE                          ! Kling-Gupta efficiency measure
      22             :   PUBLIC :: KGEnocorr                    ! KGE without correlation
      23             :   PUBLIC :: LNNSE                        ! Logarithmic Nash Sutcliffe efficiency
      24             :   PUBLIC :: MAE                          ! Mean of absolute errors
      25             :   PUBLIC :: MSE                          ! Mean of squared errors
      26             :   PUBLIC :: NSE                          ! Nash Sutcliffe efficiency
      27             :   PUBLIC :: SSE                          ! Sum of squared errors
      28             :   PUBLIC :: SAE                          ! Sum of absolute errors
      29             :   PUBLIC :: RMSE                         ! Root mean squared error
      30             :   PUBLIC :: WNSE                         ! weighted NSE
      31             : 
      32             :   ! ------------------------------------------------------------------
      33             : 
      34             :   !>      \brief Calculates bias.
      35             : 
      36             :   !>      \details  Calculates the bias
      37             :   !!
      38             :   !!      \f[BIAS = \bar y - \bar x\f]
      39             :   !!
      40             :   !!      Where \f$ \bar y \f$ and \f$ \bar x \f$ are means of the data points.
      41             :   !!
      42             :   !!      If an optinal mask is given, the calculations are over those locations that correspond to true values in the mask.
      43             :   !!      \f$ x\f$ and \f$ y\f$ can be single or double precision. The result will have the same numerical precision.
      44             :   !!
      45             :   !!      \b Example
      46             :   !!
      47             :   !!      \code{.f90}
      48             :   !!         vec1 = (/ 1., 2, 3., -999., 5., 6. /)
      49             :   !!         vec2 = (/ 1., 2, 3., -999., 5., 6. /)
      50             :   !!         m   = BIAS(vec1, vec2, mask=(vec >= 0.))
      51             :   !!         --> m = 0.0
      52             :   !!      \endcode
      53             :   !!
      54             :   !!     See also example in test directory.
      55             : 
      56             : 
      57             :   !>      \param[in]  "real(sp/dp), dimension()     :: x, y"    1D/2D/3D-array with input numbers.
      58             :   !>      \param[in]  "logical, optional     :: mask"         1D/2D/Array-array of logical values with `size(x/y)`.
      59             :   !!         If present, only those locations in vec corresponding to the true values in mask are used.
      60             : 
      61             :   !>      \returns    "real(sp/dp) :: BIAS"                     Bias.
      62             : 
      63             :   !>     \note
      64             :   !!     Input values must be floating points.
      65             : 
      66             :   !>      \authors Matthias Zink
      67             :   !>      \date Sept 2012
      68             :   INTERFACE BIAS
      69             :     MODULE PROCEDURE BIAS_sp_1d, BIAS_dp_1d, BIAS_sp_2d, BIAS_dp_2d, BIAS_sp_3d, BIAS_dp_3d
      70             :   END INTERFACE BIAS
      71             : 
      72             :   ! ------------------------------------------------------------------
      73             : 
      74             :   !>        \brief Kling-Gupta-Efficiency measure.
      75             : 
      76             :   !>        \details
      77             :   !!        The Kling-Gupta model efficiency coefficient \f$ KGE \f$ is
      78             :   !!            \f[ KGE = 1 - \sqrt{( (1-r)^2 + (1-\alpha)^2 + (1-\beta)^2 )} \f]
      79             :   !!        where \n
      80             :   !!            \f$ r \f$      = Pearson product-moment correlation coefficient \n
      81             :   !!            \f$ \alpha \f$ = ratio of simulated mean to observed mean  \n
      82             :   !!            \f$ \beta  \f$ = ratio of simulated standard deviation to
      83             :   !!                             observed standard deviation \n
      84             :   !!        This three measures are calculated between two arrays (1d, 2d, or 3d).
      85             :   !!        Usually, one is an observation and the second is a modelled variable.\n
      86             :   !!
      87             :   !!        The higher the KGE the better the observation and simulation are matching.
      88             :   !!        The upper limit of KGE is 1.\n
      89             :   !!
      90             :   !!        Therefore, if you apply a minimization algorithm to calibrate regarding
      91             :   !!        KGE you have to use the objective function
      92             :   !!            \f[ obj\_value = 1.0 - KGE \f]
      93             :   !!        which has then the optimum at 0.0.
      94             :   !!        (Like for the NSE where you always optimize 1-NSE.)\n
      95             :   !!
      96             :   !!        \b Example
      97             :   !!
      98             :   !!        \code{.f90}
      99             :   !!        para = (/ 1., 2, 3., -999., 5., 6. /)
     100             :   !!        kge = kge(x,y,mask=mask)
     101             :   !!        \endcode
     102             :   !!
     103             :   !!        \b Literature
     104             :   !!
     105             :   !>        1. Gupta, Hoshin V., et al.
     106             :   !!           _"Decomposition of the mean squared error and NSE performance criteria:
     107             :   !!           Implications for improving hydrological modelling"_.
     108             :   !!           Journal of Hydrology 377.1 (2009): 80-91.
     109             :   !!
     110             :   !>        \param[in]  "real(sp/dp)   :: x, y"           1D/2D/3D-array with input numbers
     111             :   !>        \param[in]  "logical, optional     :: mask"   1D/2D/3D-array of logical values with size(x/y).
     112             :   !>        \retval     "real(sp/dp) ::kge"               Kling-Gupta-Efficiency (value less equal 1.0)
     113             : 
     114             :   !>       \note Input values must be floating points. \n
     115             : 
     116             :   !>        \author Rohini Kumar
     117             :   !>        \date August 2014
     118             : 
     119             :   !>        \author R. Kumar, J. Mai, & O. Rakovec
     120             :   !>        \date Sep. 2014
     121             :   !!          - remove double packing of input data (bug)
     122             :   !!          - KGE instead of 1.0-KGE
     123             :   !!          - 1d, 2d, 3d, version in sp and dp
     124             : 
     125             :   INTERFACE KGE
     126             :     MODULE PROCEDURE KGE_dp_1d, KGE_dp_2d, KGE_dp_3d, KGE_sp_1d, KGE_sp_2d, KGE_sp_3d
     127             :   END INTERFACE KGE
     128             : 
     129             :   ! ------------------------------------------------------------------
     130             : 
     131             :   !>        \brief Kling-Gupta-Efficiency measure without correlation
     132             : 
     133             :   !>        \details The modified Kling-Gupta model efficiency coefficient \f$ KGEnocorr \f$ is
     134             :   !!            \f[ KGEnocorr = 1 - \sqrt{( (1-\alpha)^2 + (1-\beta)^2 )} \f]
     135             :   !!        where \n
     136             :   !!            \f$ \alpha \f$ = ratio of simulated mean to observed mean  \n
     137             :   !!            \f$ \beta  \f$ = ratio of simulated standard deviation to
     138             :   !!                             observed standard deviation \n
     139             :   !!        This two measures are calculated between two arrays (1d, 2d, or 3d).
     140             :   !!        Usually, one is an observation and the second is a modelled variable.\n
     141             :   !!
     142             :   !!        The higher the KGEnocorr the better the observation and simulation are matching.
     143             :   !!        The upper limit of KGEnocorr is 1.\n
     144             :   !!
     145             :   !!        Therefore, if you apply a minimization algorithm to calibrate regarding
     146             :   !!        KGEnocorr you have to use the objective function
     147             :   !!            \f[ obj\_value = 1.0 - KGEnocorr \f]
     148             :   !!        which has then the optimum at 0.0.
     149             :   !!        (Like for the NSE where you always optimize 1-NSE.)\n
     150             :   !!
     151             :   !!        \b Example
     152             :   !!
     153             :   !!        \code{.f90}
     154             :   !!        para = (/ 1., 2, 3., -999., 5., 6. /)
     155             :   !!        kgenocorr = kgenocorr(x,y,mask=mask)
     156             :   !!        \endcode
     157             :   !!
     158             :   !!        \b Literature
     159             :   !!
     160             :   !!        1. Gupta, Hoshin V., et al.
     161             :   !!           _"Decomposition of the mean squared error and NSE performance criteria:
     162             :   !!           Implications for improving hydrological modelling"_.
     163             :   !!           Journal of Hydrology 377.1 (2009): 80-91.
     164             :   !!
     165             :   !>        \param[in]  "real(sp/dp)   :: x, y"           1D/2D/3D-array with input numbers
     166             :   !>        \param[in]  "logical, optional     :: mask"   1D/2D/3D-array of logical values with size(x/y).
     167             :   !>        \retval     "real(sp/dp) :: kgenocorr"        Kling-Gupta-Efficiency without correlation (value less equal 1.0)
     168             : 
     169             :   !>       \note Input values must be floating points. \n
     170             : 
     171             :   !>        \author Rohini Kumar
     172             :   !>        \date Aug 2014
     173             : 
     174             :   !>        \author M. Schroen
     175             :   !>        \date Jul 2017
     176             :   !!          - add KGEnocorr
     177             : 
     178             :   !>        \author R. Kumar, J. Mai, & O. Rakovec
     179             :   !>        \date Sep 2014
     180             :   !!          - remove double packing of input data (bug)
     181             :   !!          - KGE instead of 1.0-KGE
     182             :   !!          - 1d, 2d, 3d, version in sp and dp
     183             : 
     184             :   INTERFACE KGEnocorr
     185             :     MODULE PROCEDURE KGEnocorr_dp_1d, KGEnocorr_dp_2d, KGEnocorr_dp_3d, KGEnocorr_sp_1d, KGEnocorr_sp_2d, KGEnocorr_sp_3d
     186             :   END INTERFACE KGEnocorr
     187             : 
     188             : 
     189             :   ! ------------------------------------------------------------------
     190             : 
     191             :   !>        \brief Logarithmic Nash Sutcliffe Efficiency.
     192             : 
     193             :   !>        \details  Calculates the Logarithmic Nash Sutcliffe Efficiency
     194             :   !!
     195             :   !!        \f[LNNSE = \frac{\sum_i(\ln(y_i) - \ln(x_i))^2}  {\sum_i (\ln(x_i) - \ln(\bar x))^2 }\f]
     196             :   !!
     197             :   !!        where \f$ x\f$ is the observation and \f$ y\f$ is the modelled data.\n
     198             :   !!
     199             :   !!        If an optinal mask is given, the calculations are over those locations that correspond to true values in the mask.
     200             :   !!        Note that the mask is intent inout, since values which are less or equal zero will be masked additionally.
     201             :   !!        \f$ x \f$ and \f$ y\f$ can be single or double precision. The result will have the same numerical precision.
     202             :   !!
     203             :   !!        \b Example
     204             :   !!
     205             :   !!        \code{.f90}
     206             :   !!        vec1 = (/ 1., 2, 3., -999., 5., 6. /)
     207             :   !!        vec2 = (/ 1., 2, 3., -999., 5., 6. /)
     208             :   !!        m   = LNNSE(vec1, vec2, mask=(vec >= 0.))
     209             :   !!        --> m = 1.0
     210             :   !!        \endcode
     211             :   !!
     212             :   !!        See also example in test directory.
     213             : 
     214             :   !>      \param[in]  "real(sp/dp), dimension()     :: x, y"    1D/2D/3D-array with input numbers.
     215             :   !>      \param[in]  "logical, optional     :: mask"         1D/2D/Array-array of logical values with `size(x/y)`.
     216             :   !!         If present, only those locations in vec corresponding to the true values in mask are used.
     217             : 
     218             :   !>      \retval    "real(sp/dp) :: LNNSE"                     LNNSE.
     219             : 
     220             :   !>     \note
     221             :   !!     Input values must be floating points.
     222             : 
     223             :   !>      \author Juliane Mai
     224             :   !>      \date May 2013
     225             : 
     226             :   !>      \author Rohini Kumar
     227             :   !>      \date May 2013
     228             :   !!        - mean of logQ
     229             :   INTERFACE LNNSE
     230             :     MODULE PROCEDURE LNNSE_sp_1d, LNNSE_dp_1d, LNNSE_dp_2d, LNNSE_sp_2d, LNNSE_sp_3d, LNNSE_dp_3d
     231             :   END INTERFACE LNNSE
     232             : 
     233             :   ! ------------------------------------------------------------------
     234             : 
     235             :   !>        \brief Mean absolute error.
     236             : 
     237             :   !>        \details Calculates the mean absolute error,
     238             :   !!
     239             :   !!        \f[ MAE = \sum_i\frac{|y_i - x_i|}{N_\text{mask}} \f]
     240             :   !!
     241             :   !!        If an optinal mask is given, the calculations are over those locations that correspond to true values in the mask.
     242             :   !!        \f$ x\f$ and \f$ y\f$ can be single or double precision. The result will have the same numerical precision.
     243             :   !!
     244             :   !!        \b Example
     245             :   !!
     246             :   !!        \code{.f90}
     247             :   !!        vec1 = (/ 1., 2, 3., -999., 5., 6. /)
     248             :   !!        vec2 = (/ 1., 2, 3., -999., 5., 6. /)
     249             :   !!        m   = MAE(vec1, vec2, mask=(vec >= 0.))
     250             :   !!        --> m = 0.0
     251             :   !!        \endcode
     252             :   !!
     253             :   !!        See also example in test directory.
     254             : 
     255             :   !>      \param[in]  "real(sp/dp), dimension()     :: x, y"    1D/2D/3D-array with input numbers.
     256             :   !>      \param[in]  "logical, optional     :: mask"         1D/2D/Array-array of logical values with `size(x/y)`.
     257             :   !!         If present, only those locations in vec corresponding to the true values in mask are used.
     258             : 
     259             :   !>      \returns    "real(sp/dp) :: MAE"                     MAE.
     260             : 
     261             :   !>     \note
     262             :   !!     Input values must be floating points.
     263             : 
     264             :   !>      \authors Matthias Zink
     265             :   !>      \date Sept 2012
     266             : 
     267             :   ! ------------------------------------------------------------------
     268             : 
     269             :   INTERFACE MAE
     270             :     MODULE PROCEDURE MAE_sp_1d, MAE_dp_1d, MAE_sp_2d, MAE_dp_2d, MAE_sp_3d, MAE_dp_3d
     271             :   END INTERFACE MAE
     272             : 
     273             :   ! ------------------------------------------------------------------
     274             : 
     275             :   !>        \brief Mean squared error.
     276             : 
     277             :   !>        \details Calculates the mean squared error
     278             :   !!
     279             :   !!        \f[ MSE = \sum_i\frac{(y_i - x_i)^2}{N_\text{mask}} \f]
     280             :   !!
     281             :   !!        If an optional mask is given, the calculations are over those locations that correspond to true values in the mask.
     282             :   !!        x and y can be single or double precision. The result will have the same numerical precision.
     283             :   !!
     284             :   !!        \b Example
     285             :   !!
     286             :   !!        \code{.f90}
     287             :   !!        vec1 = (/ 1., 2, 3., -999., 5., 6. /)
     288             :   !!        vec2 = (/ 1., 2, 3., -999., 5., 6. /)
     289             :   !!        m   = MSE(vec1, vec2, mask=(vec >= 0.))
     290             :   !!        --> m = 0.0
     291             :   !!        \endcode
     292             :   !!
     293             :   !!        See also example in test directory.
     294             : 
     295             :   !>      \param[in]  "real(sp/dp), dimension()     :: x, y"    1D/2D/3D-array with input numbers.
     296             :   !>      \param[in]  "logical, optional     :: mask"         1D/2D/Array-array of logical values with `size(x/y)`.
     297             :   !!         If present, only those locations in vec corresponding to the true values in mask are used.
     298             :   !>      \retval    "real(sp/dp) :: MSE"                     MSE.
     299             : 
     300             :   !>     \note
     301             :   !!     Input values must be floating points.
     302             : 
     303             :   !>        \authors Matthias Zink
     304             :   !>        \date Sept 2012
     305             : 
     306             :   ! ------------------------------------------------------------------
     307             : 
     308             :   INTERFACE MSE
     309             :     MODULE PROCEDURE MSE_sp_1d, MSE_dp_1d, MSE_sp_2d, MSE_dp_2d, MSE_sp_3d, MSE_dp_3d
     310             :   END INTERFACE MSE
     311             : 
     312             :   ! ------------------------------------------------------------------
     313             : 
     314             :   !>        \brief Nash Sutcliffe Efficiency.
     315             : 
     316             :   !>        \details Calculates the Nash Sutcliffe Efficiency
     317             :   !!
     318             :   !!        \f[NSE = \frac{\sum_i(y_i - x_i)^2}  {\sum_i (x_i - \bar x)^2 }\f]
     319             :   !!
     320             :   !!        where \f$ x\f$ is the observation and \f$ y\f$ is the modelled data.
     321             :   !!
     322             :   !!        If an optinal mask is given, the calculations are over those locations that correspond to true values in the mask.
     323             :   !!        \f$ x\f$ and \f$ y\f$ can be single or double precision. The result will have the same numerical precision.
     324             :   !!
     325             :   !!        \b Example
     326             :   !!
     327             :   !!        \code{.f90}
     328             :   !!        vec1 = (/ 1., 2, 3., -999., 5., 6. /)
     329             :   !!        vec2 = (/ 1., 2, 3., -999., 5., 6. /)
     330             :   !!        m   = NSE(vec1, vec2, mask=(vec >= 0.))
     331             :   !!        --> m = 1.0
     332             :   !!        \endcode
     333             :   !!
     334             :   !!        See also example in test directory.
     335             :   !!
     336             :   !!        \b Literature
     337             :   !!
     338             :   !!        1. Nash, J., & Sutcliffe, J. (1970). _River flow forecasting through conceptual models part I: A discussion of
     339             :   !!           principles_. Journal of Hydrology, 10(3), 282-290. doi:10.1016/0022-1694(70)90255-6
     340             :   !!
     341             :   !>      \param[in]  "real(sp/dp), dimension()     :: x, y"    1D/2D/3D-array with input numbers.
     342             :   !>      \param[in]  "logical, optional     :: mask"         1D/2D/Array-array of logical values with `size(x/y)`.
     343             :   !!         If present, only those locations in vec corresponding to the true values in mask are used.
     344             :   !>      \retval    "real(sp/dp) :: NSE"                     NSE.
     345             : 
     346             :   !>     \note
     347             :   !!     Input values must be floating points.
     348             : 
     349             :   !>        \authors Matthias Zink
     350             :   !>        \date Sept 2012
     351             : 
     352             :   ! ------------------------------------------------------------------
     353             : 
     354             :   INTERFACE NSE
     355             :     MODULE PROCEDURE NSE_sp_1d, NSE_dp_1d, NSE_dp_2d, NSE_sp_2d, NSE_sp_3d, NSE_dp_3d
     356             :   END INTERFACE NSE
     357             : 
     358             :   ! ------------------------------------------------------------------
     359             : 
     360             :   !>     \brief  Sum of absolute errors.
     361             : 
     362             :   !>     \details Calculates the sum of absolute errors
     363             :   !!
     364             :   !!     \f[ SAE = \sum_i|y_i - x_i| \f]
     365             :   !!
     366             :   !!     If an optional mask is given, the calculations are over those locations that correspond to true values in the mask.
     367             :   !!     \f$ x\f$ and \f$ y\f$ can be single or double precision. The result will have the same numerical precision.
     368             :   !!
     369             :   !!     \b Example
     370             :   !!
     371             :   !!     \code{.f90}
     372             :   !!     vec1 = (/ 1., 2, 3., -999., 5., 6. /)
     373             :   !!     vec2 = (/ 1., 2, 3., -999., 5., 6. /)
     374             :   !!     m = NSE(vec1, vec2, mask=(vec >= 0.))
     375             :   !!     --> m = 0.0
     376             :   !!     \endcode
     377             :   !!
     378             :   !!     See also example in test directory.
     379             : 
     380             :   !>     \param[in]  "real(sp/dp), dimension()     :: x, y"    1D/2D/3D-array with input numbers.
     381             :   !>     \param[in]  "logical, optional     :: mask"         1D/2D/Array-array of logical values with `size(x/y)`.
     382             :   !!         If present, only those locations in vec corresponding to the true values in mask are used.
     383             : 
     384             :   !>     \returns    "real(sp/dp) :: NSE"                     NSE.
     385             : 
     386             :   !>     \note
     387             :   !!     Input values must be floating points.
     388             : 
     389             :   !>     \authors Matthias Zink
     390             :   !>     \date Sept 2012
     391             : 
     392             :   ! ------------------------------------------------------------------
     393             : 
     394             :   INTERFACE SAE
     395             :     MODULE PROCEDURE SAE_sp_1d, SAE_dp_1d, SAE_sp_2d, SAE_dp_2d, SAE_sp_3d, SAE_dp_3d
     396             :   END INTERFACE SAE
     397             : 
     398             :   ! ------------------------------------------------------------------
     399             : 
     400             :   !>    \brief Sum of squared errors
     401             : 
     402             :   !>    \details Calculates the sum of squared errors
     403             :   !!
     404             :   !!    \f[ SSE = \sum_i(y_i - x_i)^2 \f]
     405             :   !!
     406             :   !!    If an optional mask is given, the calculations are over those locations that correspond to true values in the mask.
     407             :   !!    \f$ x\f$ and \f$ y\f$ can be single or double precision. The result will have the same numerical precision.
     408             :   !!
     409             :   !!    \b Example
     410             :   !!
     411             :   !!    \code{.f90}
     412             :   !!    vec1 = (/ 1., 2, 3., -999., 5., 6. /)
     413             :   !!    vec2 = (/ 1., 2, 3., -999., 5., 6. /)
     414             :   !!    m   = SSE(vec1, vec2, mask=(vec >= 0.))
     415             :   !!    --> m = 0.0
     416             :   !!    \endcode
     417             :   !!
     418             :   !!    See also example in test directory.
     419             : 
     420             :   !>    \param[in]  "real(sp/dp), dimension()     :: x, y"    1D/2D/3D-array with input numbers.
     421             :   !>    \param[in]  "logical, optional     :: mask"         1D/2D/Array-array of logical values with `size(x/y)`.
     422             :   !!    If present, only those locations in vec corresponding to the true values in mask are used.
     423             :   !>    \retval    "real(sp/dp) :: SSE"                     SSE.
     424             : 
     425             :   !>    \note
     426             :   !!    Input values must be floating points.
     427             : 
     428             :   !>    \authors Matthias Zink
     429             :   !>    \date Sept 2012
     430             : 
     431             :   ! ------------------------------------------------------------------
     432             : 
     433             :   INTERFACE SSE
     434             :     MODULE PROCEDURE SSE_sp_1d, SSE_dp_1d, SSE_sp_2d, SSE_dp_2d, SSE_sp_3d, SSE_dp_3d
     435             :   END INTERFACE SSE
     436             : 
     437             :   ! ------------------------------------------------------------------
     438             : 
     439             :   !>    \brief RMS Error.
     440             : 
     441             :   !>    \details Calculates the root-mean-square error
     442             :   !!
     443             :   !!    \f[ RMSE = \sqrt{\frac{\sum_i{(y_i - x_i)^2}} {{N_\text{count}}}} \f]
     444             :   !!
     445             :   !!    If an optional mask is given, the calculations are over those locations that correspond to true values in the mask.
     446             :   !!    \f$ x\f$ and \f$ y\f$ can be single or double precision. The result will have the same numerical precision.
     447             :   !!
     448             :   !!    \b Example
     449             :   !!
     450             :   !!    \code{.f90}
     451             :   !!    vec1 = (/ 1., 2, 3., -999., 5., 6. /)
     452             :   !!    vec2 = (/ 1., 2, 3., -999., 5., 6. /)
     453             :   !!    m   = RMSE(vec1, vec2, mask=(vec >= 0.))
     454             :   !!    --> m = 0.0
     455             :   !!    \endcode
     456             :   !!
     457             :   !!    See also example in test directory.
     458             : 
     459             :   !>    \param[in]  "real(sp/dp), dimension()     :: x, y"    1D/2D/3D-array with input numbers.
     460             :   !>    \param[in]  "logical, optional     :: mask"         1D/2D/Array-array of logical values with `size(x/y)`.
     461             :   !!    If present, only those locations in vec corresponding to the true values in mask are used.
     462             :   !>    \retval    "real(sp/dp) :: RMSE"                     RMSE.
     463             : 
     464             :   !>    \note
     465             :   !!    Input values must be floating points.
     466             : 
     467             :   !>    \authors Matthias Zink
     468             :   !>    \date Sept 2012
     469             : 
     470             :   ! ------------------------------------------------------------------
     471             : 
     472             :   INTERFACE RMSE
     473             :     MODULE PROCEDURE RMSE_sp_1d, RMSE_dp_1d, RMSE_sp_2d, RMSE_dp_2d, RMSE_sp_3d, RMSE_dp_3d
     474             :   END INTERFACE RMSE
     475             : 
     476             :   ! ------------------------------------------------------------------
     477             : 
     478             :   !>    \brief weighted Nash Sutcliffe Efficiency.
     479             : 
     480             :   !>    Calculates the weighted Nash Sutcliffe Efficiency
     481             :   !!
     482             :   !!    \f[ wNSE = \frac{\sum_i {x_i (y_i - x_i)^2}} {\sum_i{ x_i (x_i - \bar x)^2}} \f]
     483             :   !!
     484             :   !!    where \f$ x\f$ is the observation and \f$ y\f$ is the modelled data.
     485             :   !!    This objective function is introduced in Hundecha and Bardossy, 2004.
     486             :   !!
     487             :   !!    If an optinal mask is given, the calculations are over those locations that correspond to true values in the mask.
     488             :   !!    \f$ x\f$ and \f$ y\f$ can be single or double precision. The result will have the same numerical precision.
     489             :   !!
     490             :   !!    \b Example
     491             :   !!
     492             :   !!    \code{.f90}
     493             :   !!    vec1 = (/ 1., 2, 3., -999., 5., 6. /)
     494             :   !!    vec2 = (/ 1., 2, 3., -999., 5., 6. /)
     495             :   !!    m   = wNSE(vec1, vec2, mask=(vec >= 0.))
     496             :   !!    --> m = 1.0
     497             :   !!    \endcode
     498             :   !!
     499             :   !!    See also example in test directory.
     500             :   !!
     501             :   !!    \b Literature
     502             :   !!
     503             :   !!    1.   Nash, J., & Sutcliffe, J. (1970). _River flow forecasting through conceptual models part I: A discussion of
     504             :   !!                principles_. Journal of Hydrology, 10(3), 282-290. doi:10.1016/0022-1694(70)90255-6\n
     505             :   !!    2.   Hundecha and Bardossy (2004). _Modeling of the effect of land use changes on the runoff generation of a river
     506             :   !!                domain through parameter regionalization of a watershed model_. Journal of Hydrology, 292, 281-295
     507             :   !!
     508             :   !>    \param[in]  "real(sp/dp), dimension()     :: x, y"    1D/2D/3D-array with input numbers.
     509             :   !>    \param[in]  "logical, optional     :: mask"         1D/2D/Array-array of logical values with `size(x/y)`.
     510             :   !!    If present, only those locations in vec corresponding to the true values in mask are used.
     511             :   !>    \retval    "real(sp/dp) :: wNSE"                     wNSE.
     512             : 
     513             :   !>    \note
     514             :   !!    Input values must be floating points.
     515             : 
     516             :   !>    \authors Matthias Zink & Bjoern Guse
     517             :   !>    \date May 2018
     518             : 
     519             :   ! ------------------------------------------------------------------
     520             : 
     521             :   INTERFACE wNSE
     522             :      MODULE PROCEDURE wNSE_sp_1d, wNSE_dp_1d, wNSE_dp_2d, wNSE_sp_2d, wNSE_sp_3d, wNSE_dp_3d
     523             :   END INTERFACE wNSE
     524             : 
     525             :   ! ------------------------------------------------------------------
     526             : 
     527             :   PRIVATE
     528             : 
     529             :   ! ------------------------------------------------------------------
     530             : 
     531             : CONTAINS
     532             : 
     533             :   ! ------------------------------------------------------------------
     534             : 
     535           4 :   FUNCTION BIAS_sp_1d(x, y, mask)
     536             : 
     537             :     USE mo_moment, ONLY : average
     538             : 
     539             :     IMPLICIT NONE
     540             : 
     541             :     REAL(sp), DIMENSION(:), INTENT(IN) :: x, y
     542             :     LOGICAL, DIMENSION(:), OPTIONAL, INTENT(IN) :: mask
     543             :     REAL(sp) :: BIAS_sp_1d
     544             : 
     545             :     INTEGER(i4) :: n
     546             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
     547           2 :     LOGICAL, DIMENSION(size(x)) :: maske
     548             : 
     549           2 :     if (present(mask)) then
     550           2 :       shapemask = shape(mask)
     551             :     else
     552           2 :       shapemask = shape(x)
     553             :     end if
     554             :     !
     555           8 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
     556           0 :             stop 'BIAS_sp_1d: shapes of inputs(x,y) or mask are not matching'
     557             :     !
     558           2 :     if (present(mask)) then
     559          56 :       maske = mask
     560          55 :       n = count(maske)
     561             :     else
     562          55 :       maske = .true.
     563           1 :       n = size(x)
     564             :     end if
     565             :     !
     566           2 :     if (n .LE. 1_i4) stop 'BIAS_sp_1d: number of arguments must be at least 2'
     567             :     !
     568           2 :     BIAS_sp_1d = average(y, mask = maske) - average(x, mask = maske)
     569             : 
     570           2 :   END FUNCTION BIAS_sp_1d
     571             : 
     572           4 :   FUNCTION BIAS_dp_1d(x, y, mask)
     573             : 
     574           2 :     USE mo_moment, ONLY : average
     575             : 
     576             :     IMPLICIT NONE
     577             : 
     578             :     REAL(dp), DIMENSION(:), INTENT(IN) :: x, y
     579             :     LOGICAL, DIMENSION(:), OPTIONAL, INTENT(IN) :: mask
     580             :     REAL(dp) :: BIAS_dp_1d
     581             : 
     582             :     INTEGER(i4) :: n
     583             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
     584           2 :     LOGICAL, DIMENSION(size(x)) :: maske
     585             : 
     586           2 :     if (present(mask)) then
     587           2 :       shapemask = shape(mask)
     588             :     else
     589           2 :       shapemask = shape(x)
     590             :     end if
     591             :     !
     592           8 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
     593           0 :             stop 'BIAS_dp_1d: shapes of inputs(x,y) or mask are not matching'
     594             :     !
     595           2 :     if (present(mask)) then
     596          56 :       maske = mask
     597          55 :       n = count(maske)
     598             :     else
     599          55 :       maske = .true.
     600           1 :       n = size(x)
     601             :     end if
     602           2 :     if (n .LE. 1_i4) stop 'BIAS_dp_1d: number of arguments must be at least 2'
     603             :     !
     604           2 :     BIAS_dp_1d = average(y, mask = maske) - average(x, mask = maske)
     605             : 
     606           2 :   END FUNCTION BIAS_dp_1d
     607             : 
     608           4 :   FUNCTION BIAS_sp_2d(x, y, mask)
     609             : 
     610           2 :     USE mo_moment, ONLY : average
     611             : 
     612             :     IMPLICIT NONE
     613             : 
     614             :     REAL(sp), DIMENSION(:, :), INTENT(IN) :: x, y
     615             :     LOGICAL, DIMENSION(:, :), OPTIONAL, INTENT(IN) :: mask
     616             :     REAL(sp) :: BIAS_sp_2d
     617             : 
     618             :     INTEGER(i4) :: n
     619             : 
     620             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
     621           2 :     LOGICAL, DIMENSION(size(x, dim = 1), size(x, dim = 2)) :: maske
     622             : 
     623           2 :     if (present(mask)) then
     624           3 :       shapemask = shape(mask)
     625             :     else
     626           3 :       shapemask = shape(x)
     627             :     end if
     628             :     !
     629          12 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
     630           0 :             stop 'BIAS_sp_2d: shapes of inputs(x,y) or mask are not matching'
     631             :     !
     632           2 :     if (present(mask)) then
     633          62 :       maske = mask
     634          61 :       n = count(maske)
     635             :     else
     636          61 :       maske = .true.
     637           1 :       n = size(x, dim = 1) * size(x, dim = 2)
     638             :     end if
     639             :     !
     640           2 :     if (n .LE. 1_i4) stop 'BIAS_sp_2d: number of arguments must be at least 2'
     641             :     !
     642             :     BIAS_sp_2d = average(reshape(y, (/size(y, dim = 1) * size(y, dim = 2)/)), &
     643             :             mask = reshape(maske, (/size(y, dim = 1) * size(y, dim = 2)/))) - &
     644             :             average(reshape(x, (/size(x, dim = 1) * size(x, dim = 2)/)), &
     645          10 :                     mask = reshape(maske, (/size(x, dim = 1) * size(x, dim = 2)/)))
     646             :     !
     647           2 :   END FUNCTION BIAS_sp_2d
     648             : 
     649           4 :   FUNCTION BIAS_dp_2d(x, y, mask)
     650             : 
     651           2 :     USE mo_moment, ONLY : average
     652             : 
     653             :     IMPLICIT NONE
     654             : 
     655             :     REAL(dp), DIMENSION(:, :), INTENT(IN) :: x, y
     656             :     LOGICAL, DIMENSION(:, :), OPTIONAL, INTENT(IN) :: mask
     657             :     REAL(dp) :: BIAS_dp_2d
     658             : 
     659             :     INTEGER(i4) :: n
     660             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
     661           2 :     LOGICAL, DIMENSION(size(x, dim = 1), size(x, dim = 2)) :: maske
     662             : 
     663           2 :     if (present(mask)) then
     664           3 :       shapemask = shape(mask)
     665             :     else
     666           3 :       shapemask = shape(x)
     667             :     end if
     668             :     !
     669          12 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
     670           0 :             stop 'BIAS_dp_2d: shapes of inputs(x,y) or mask are not matching'
     671             :     !
     672           2 :     if (present(mask)) then
     673          62 :       maske = mask
     674          61 :       n = count(maske)
     675             :     else
     676          61 :       maske = .true.
     677           1 :       n = size(x, dim = 1) * size(x, dim = 2)
     678             :     end if
     679             :     !
     680           2 :     if (n .LE. 1_i4) stop 'BIAS_dp_2d: number of arguments must be at least 2'
     681             :     !
     682             :     BIAS_dp_2d = average(reshape(y, (/size(y, dim = 1) * size(y, dim = 2)/)), &
     683             :             mask = reshape(maske, (/size(y, dim = 1) * size(y, dim = 2)/))) - &
     684             :             average(reshape(x, (/size(x, dim = 1) * size(x, dim = 2)/)), &
     685          10 :                     mask = reshape(maske, (/size(x, dim = 1) * size(x, dim = 2)/)))
     686             :     !
     687           2 :   END FUNCTION BIAS_dp_2d
     688             : 
     689           4 :   FUNCTION BIAS_sp_3d(x, y, mask)
     690             : 
     691           2 :     USE mo_moment, ONLY : average
     692             : 
     693             :     IMPLICIT NONE
     694             : 
     695             :     REAL(sp), DIMENSION(:, :, :), INTENT(IN) :: x, y
     696             :     LOGICAL, DIMENSION(:, :, :), OPTIONAL, INTENT(IN) :: mask
     697             :     REAL(sp) :: BIAS_sp_3d
     698             : 
     699             :     INTEGER(i4) :: n
     700             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
     701             :     LOGICAL, DIMENSION(size(x, dim = 1), &
     702           2 :             size(x, dim = 2), size(x, dim = 3)) :: maske
     703             : 
     704           2 :     if (present(mask)) then
     705           4 :       shapemask = shape(mask)
     706             :     else
     707           4 :       shapemask = shape(x)
     708             :     end if
     709             :     !
     710          16 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
     711           0 :             stop 'BIAS_sp_3d: shapes of inputs(x,y) or mask are not matching'
     712             :     !
     713           2 :     if (present(mask)) then
     714         490 :       maske = mask
     715         489 :       n = count(maske)
     716             :     else
     717         489 :       maske = .true.
     718           1 :       n = size(x, dim = 1) * size(x, dim = 2) * size(x, dim = 3)
     719             :     end if
     720             :     !
     721             :     ! not really sopisticated, it has to be checked if the 3 numbers of x and y are matching in arry position
     722           2 :     if (n .LE. 1_i4) stop 'BIAS_sp_3d: number of arguments must be at least 2'
     723             :     !
     724             :     BIAS_sp_3d = average(reshape(y, (/size(y, dim = 1) * size(y, dim = 2) * size(y, dim = 3)/)), &
     725             :             mask = reshape(maske, (/size(y, dim = 1) * size(y, dim = 2) * size(y, dim = 3)/))) - &
     726             :             average(reshape(x, (/size(x, dim = 1) * size(x, dim = 2) * size(x, dim = 3)/)), &
     727          10 :                     mask = reshape(maske, (/size(x, dim = 1) * size(x, dim = 2) * size(x, dim = 3)/)))
     728             :     !
     729           2 :   END FUNCTION BIAS_sp_3d
     730             : 
     731           4 :   FUNCTION BIAS_dp_3d(x, y, mask)
     732             : 
     733           2 :     USE mo_moment, ONLY : average
     734             : 
     735             :     IMPLICIT NONE
     736             : 
     737             :     REAL(dp), DIMENSION(:, :, :), INTENT(IN) :: x, y
     738             :     LOGICAL, DIMENSION(:, :, :), OPTIONAL, INTENT(IN) :: mask
     739             :     REAL(dp) :: BIAS_dp_3d
     740             : 
     741             :     INTEGER(i4) :: n
     742             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
     743             :     LOGICAL, DIMENSION(size(x, dim = 1), &
     744           2 :             size(x, dim = 2), size(x, dim = 3)) :: maske
     745             : 
     746           2 :     if (present(mask)) then
     747           4 :       shapemask = shape(mask)
     748             :     else
     749           4 :       shapemask = shape(x)
     750             :     end if
     751             :     !
     752          16 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
     753           0 :             stop 'BIAS_dp_3d: shapes of inputs(x,y) or mask are not matching'
     754             :     !
     755           2 :     if (present(mask)) then
     756         490 :       maske = mask
     757         489 :       n = count(maske)
     758             :     else
     759         489 :       maske = .true.
     760           1 :       n = size(x, dim = 1) * size(x, dim = 2) * size(x, dim = 3)
     761             :     end if
     762             :     !
     763             :     ! not really sopisticated, it has to be checked if the 3 numbers of x and y are matching in arry position
     764           2 :     if (n .LE. 1_i4) stop 'BIAS_dp_3d: number of arguments must be at least 2'
     765             :     !
     766             :     BIAS_dp_3d = average(reshape(y, (/size(y, dim = 1) * size(y, dim = 2) * size(y, dim = 3)/)), &
     767             :             mask = reshape(maske, (/size(y, dim = 1) * size(y, dim = 2) * size(y, dim = 3)/))) - &
     768             :             average(reshape(x, (/size(x, dim = 1) * size(x, dim = 2) * size(x, dim = 3)/)), &
     769          10 :                     mask = reshape(maske, (/size(x, dim = 1) * size(x, dim = 2) * size(x, dim = 3)/)))
     770             :     !
     771           2 :   END FUNCTION BIAS_dp_3d
     772             : 
     773             :   ! ------------------------------------------------------------------
     774             : 
     775           4 :   FUNCTION KGE_sp_1d(x, y, mask)
     776             : 
     777           2 :     USE mo_moment, ONLY : average, stddev, correlation
     778             : 
     779             :     IMPLICIT NONE
     780             : 
     781             :     REAL(sp), DIMENSION(:), INTENT(IN) :: x, y
     782             :     LOGICAL, DIMENSION(:), OPTIONAL, INTENT(IN) :: mask
     783             :     REAL(sp) :: KGE_sp_1d
     784             : 
     785             :     ! local variables
     786             :     INTEGER(i4) :: n
     787             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
     788           2 :     LOGICAL, DIMENSION(size(x)) :: maske
     789             : 
     790           2 :     REAL(sp) :: mu_Obs, mu_Sim       ! Mean          of x and y
     791           2 :     REAL(sp) :: sigma_Obs, sigma_Sim ! Standard dev. of x and y
     792           2 :     REAL(sp) :: pearson_coor         ! Pearson Corr. of x and y
     793             : 
     794           2 :     if (present(mask)) then
     795           2 :       shapemask = shape(mask)
     796             :     else
     797           2 :       shapemask = shape(x)
     798             :     end if
     799           8 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
     800           0 :             stop 'KGE_sp_1d: shapes of inputs(x,y) or mask are not matching'
     801             :     !
     802           2 :     if (present(mask)) then
     803          55 :       maske = mask
     804          55 :       n = count(maske)
     805             :     else
     806          55 :       maske = .true.
     807           1 :       n = size(x)
     808             :     end if
     809           2 :     if (n .LE. 1_i4) stop 'KGE_sp_1d: sample size must be at least 2'
     810             : 
     811             :     ! Mean
     812           2 :     mu_Obs = average(x, mask = maske)
     813           2 :     mu_Sim = average(y, mask = maske)
     814             :     ! Standard Deviation
     815           2 :     sigma_Obs = stddev(x, mask = maske)
     816           2 :     sigma_Sim = stddev(y, mask = maske)
     817             :     ! Pearson product-moment correlation coefficient is with (N-1) not N
     818           2 :     pearson_coor = correlation(x, y, mask = maske) * real(n, sp) / real(n - 1, sp)
     819             :     !
     820             :     KGE_sp_1d = 1.0 - SQRT(&
     821             :             (1.0_sp - (mu_Sim / mu_Obs))**2 + &
     822             :                     (1.0_sp - (sigma_Sim / sigma_Obs))**2 + &
     823             :                     (1.0_sp - pearson_coor)**2             &
     824           2 :             )
     825             : 
     826           2 :   END FUNCTION KGE_sp_1d
     827             : 
     828           4 :   FUNCTION KGE_sp_2d(x, y, mask)
     829             : 
     830           2 :     USE mo_moment, ONLY : average, stddev, correlation
     831             : 
     832             :     IMPLICIT NONE
     833             : 
     834             :     REAL(sp), DIMENSION(:, :), INTENT(IN) :: x, y
     835             :     LOGICAL, DIMENSION(:, :), OPTIONAL, INTENT(IN) :: mask
     836             :     REAL(sp) :: KGE_sp_2d
     837             : 
     838             :     ! local variables
     839             :     INTEGER(i4) :: n
     840             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
     841           2 :     LOGICAL, DIMENSION(size(x, dim = 1), size(x, dim = 2)) :: maske
     842           2 :     REAL(sp) :: mu_Obs, mu_Sim       ! Mean          of x and y
     843           2 :     REAL(sp) :: sigma_Obs, sigma_Sim ! Standard dev. of x and y
     844           2 :     REAL(sp) :: pearson_coor         ! Pearson Corr. of x and y
     845             : 
     846           2 :     if (present(mask)) then
     847           3 :       shapemask = shape(mask)
     848             :     else
     849           3 :       shapemask = shape(x)
     850             :     end if
     851          12 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
     852           0 :             stop 'KGE_sp_2d: shapes of inputs(x,y) or mask are not matching'
     853             :     !
     854           2 :     if (present(mask)) then
     855          61 :       maske = mask
     856          61 :       n = count(maske)
     857             :     else
     858          61 :       maske = .true.
     859           3 :       n = size(x)
     860             :     end if
     861           2 :     if (n .LE. 1_i4) stop 'KGE_sp_2d: sample size must be at least 2'
     862             : 
     863             :     ! Mean
     864             :     mu_Obs = average(&
     865             :             reshape(x(:, :), (/size(x, dim = 1) * size(x, dim = 2)/)), &
     866           6 :             mask = reshape(maske(:, :), (/size(x, dim = 1) * size(x, dim = 2)/)))
     867             :     mu_Sim = average(&
     868             :             reshape(y(:, :), (/size(y, dim = 1) * size(y, dim = 2)/)), &
     869           6 :             mask = reshape(maske(:, :), (/size(y, dim = 1) * size(y, dim = 2)/)))
     870             :     ! Standard Deviation
     871             :     sigma_Obs = stddev(&
     872             :             reshape(x(:, :), (/size(x, dim = 1) * size(x, dim = 2)/)), &
     873           6 :             mask = reshape(maske(:, :), (/size(x, dim = 1) * size(x, dim = 2)/)))
     874             :     sigma_Sim = stddev(&
     875             :             reshape(y(:, :), (/size(y, dim = 1) * size(y, dim = 2)/)), &
     876           6 :             mask = reshape(maske(:, :), (/size(y, dim = 1) * size(y, dim = 2)/)))
     877             :     ! Pearson product-moment correlation coefficient is with (N-1) not N
     878             :     pearson_coor = correlation(&
     879             :             reshape(x(:, :), (/size(x, dim = 1) * size(x, dim = 2)/)), &
     880             :             reshape(y(:, :), (/size(y, dim = 1) * size(y, dim = 2)/)), &
     881             :             mask = reshape(maske(:, :), (/size(y, dim = 1) * size(y, dim = 2)/))) * &
     882           8 :             real(n, sp) / real(n - 1, sp)
     883             :     !
     884             :     KGE_sp_2d = 1.0 - SQRT(&
     885             :             (1.0_sp - (mu_Sim / mu_Obs))**2 + &
     886             :                     (1.0_sp - (sigma_Sim / sigma_Obs))**2 + &
     887             :                     (1.0_sp - pearson_coor)**2             &
     888           2 :             )
     889             : 
     890           2 :   END FUNCTION KGE_sp_2d
     891             : 
     892           4 :   FUNCTION KGE_sp_3d(x, y, mask)
     893             : 
     894           2 :     USE mo_moment, ONLY : average, stddev, correlation
     895             : 
     896             :     IMPLICIT NONE
     897             : 
     898             :     REAL(sp), DIMENSION(:, :, :), INTENT(IN) :: x, y
     899             :     LOGICAL, DIMENSION(:, :, :), OPTIONAL, INTENT(IN) :: mask
     900             :     REAL(sp) :: KGE_sp_3d
     901             : 
     902             :     ! local variables
     903             :     INTEGER(i4) :: n
     904             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
     905           2 :     LOGICAL, DIMENSION(size(x, dim = 1), size(x, dim = 2), size(x, dim = 3)) :: maske
     906           2 :     REAL(sp) :: mu_Obs, mu_Sim       ! Mean          of x and y
     907           2 :     REAL(sp) :: sigma_Obs, sigma_Sim ! Standard dev. of x and y
     908           2 :     REAL(sp) :: pearson_coor         ! Pearson Corr. of x and y
     909             : 
     910           2 :     if (present(mask)) then
     911           4 :       shapemask = shape(mask)
     912             :     else
     913           4 :       shapemask = shape(x)
     914             :     end if
     915          16 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
     916           0 :             stop 'KGE_sp_3d: shapes of inputs(x,y) or mask are not matching'
     917             :     !
     918           2 :     if (present(mask)) then
     919         489 :       maske = mask
     920         489 :       n = count(maske)
     921             :     else
     922         489 :       maske = .true.
     923           4 :       n = size(x)
     924             :     end if
     925           2 :     if (n .LE. 1_i4) stop 'KGE_sp_3d: sample size must be at least 2'
     926             : 
     927             :     ! Mean
     928             :     mu_Obs = average(&
     929             :             reshape(x(:, :, :), (/size(x, dim = 1) * size(x, dim = 2) * size(x, dim = 3)/)), &
     930           6 :             mask = reshape(maske(:, :, :), (/size(x, dim = 1) * size(x, dim = 2) * size(x, dim = 3)/)))
     931             :     mu_Sim = average(&
     932             :             reshape(y(:, :, :), (/size(y, dim = 1) * size(y, dim = 2) * size(y, dim = 3)/)), &
     933           6 :             mask = reshape(maske(:, :, :), (/size(y, dim = 1) * size(y, dim = 2) * size(y, dim = 3)/)))
     934             :     ! Standard Deviation
     935             :     sigma_Obs = stddev(&
     936             :             reshape(x(:, :, :), (/size(x, dim = 1) * size(x, dim = 2) * size(x, dim = 3)/)), &
     937           6 :             mask = reshape(maske(:, :, :), (/size(x, dim = 1) * size(x, dim = 2) * size(x, dim = 3)/)))
     938             :     sigma_Sim = stddev(&
     939             :             reshape(y(:, :, :), (/size(y, dim = 1) * size(y, dim = 2) * size(y, dim = 3)/)), &
     940           6 :             mask = reshape(maske(:, :, :), (/size(y, dim = 1) * size(y, dim = 2) * size(y, dim = 3)/)))
     941             :     ! Pearson product-moment correlation coefficient is with (N-1) not N
     942             :     pearson_coor = correlation(&
     943             :             reshape(x(:, :, :), (/size(x, dim = 1) * size(x, dim = 2) * size(x, dim = 3)/)), &
     944             :             reshape(y(:, :, :), (/size(y, dim = 1) * size(y, dim = 2) * size(y, dim = 3)/)), &
     945             :             mask = reshape(maske(:, :, :), (/size(y, dim = 1) * size(y, dim = 2) * size(y, dim = 3)/))) * &
     946           8 :             real(n, sp) / real(n - 1, sp)
     947             :     !
     948             :     KGE_sp_3d = 1.0 - SQRT(&
     949             :             (1.0_sp - (mu_Sim / mu_Obs))**2 + &
     950             :                     (1.0_sp - (sigma_Sim / sigma_Obs))**2 + &
     951             :                     (1.0_sp - pearson_coor)**2             &
     952           2 :             )
     953             : 
     954           2 :   END FUNCTION KGE_sp_3d
     955             : 
     956           4 :   FUNCTION KGE_dp_1d(x, y, mask)
     957             : 
     958           2 :     USE mo_moment, ONLY : average, stddev, correlation
     959             : 
     960             :     IMPLICIT NONE
     961             : 
     962             :     REAL(dp), DIMENSION(:), INTENT(IN) :: x, y
     963             :     LOGICAL, DIMENSION(:), OPTIONAL, INTENT(IN) :: mask
     964             :     REAL(dp) :: KGE_dp_1d
     965             : 
     966             :     ! local variables
     967             :     INTEGER(i4) :: n
     968             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
     969           2 :     LOGICAL, DIMENSION(size(x)) :: maske
     970             : 
     971           2 :     REAL(dp) :: mu_Obs, mu_Sim       ! Mean          of x and y
     972           2 :     REAL(dp) :: sigma_Obs, sigma_Sim ! Standard dev. of x and y
     973           2 :     REAL(dp) :: pearson_coor         ! Pearson Corr. of x and y
     974             : 
     975           2 :     if (present(mask)) then
     976           2 :       shapemask = shape(mask)
     977             :     else
     978           2 :       shapemask = shape(x)
     979             :     end if
     980           8 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
     981           0 :             stop 'KGE_dp_1d: shapes of inputs(x,y) or mask are not matching'
     982             :     !
     983           2 :     if (present(mask)) then
     984          55 :       maske = mask
     985          55 :       n = count(maske)
     986             :     else
     987          55 :       maske = .true.
     988           1 :       n = size(x)
     989             :     end if
     990           2 :     if (n .LE. 1_i4) stop 'KGE_dp_1d: sample size must be at least 2'
     991             : 
     992             :     ! Mean
     993           2 :     mu_Obs = average(x, mask = maske)
     994           2 :     mu_Sim = average(y, mask = maske)
     995             :     ! Standard Deviation
     996           2 :     sigma_Obs = stddev(x, mask = maske)
     997           2 :     sigma_Sim = stddev(y, mask = maske)
     998             :     ! Pearson product-moment correlation coefficient is with (N-1) not N
     999           2 :     pearson_coor = correlation(x, y, mask = maske) * real(n, dp) / real(n - 1, dp)
    1000             :     !
    1001             :     KGE_dp_1d = 1.0 - SQRT(&
    1002             :             (1.0_dp - (mu_Sim / mu_Obs))**2 + &
    1003             :                     (1.0_dp - (sigma_Sim / sigma_Obs))**2 + &
    1004             :                     (1.0_dp - pearson_coor)**2             &
    1005           2 :             )
    1006             : 
    1007           2 :   END FUNCTION KGE_dp_1d
    1008             : 
    1009           4 :   FUNCTION KGE_dp_2d(x, y, mask)
    1010             : 
    1011           2 :     USE mo_moment, ONLY : average, stddev, correlation
    1012             : 
    1013             :     IMPLICIT NONE
    1014             : 
    1015             :     REAL(dp), DIMENSION(:, :), INTENT(IN) :: x, y
    1016             :     LOGICAL, DIMENSION(:, :), OPTIONAL, INTENT(IN) :: mask
    1017             :     REAL(dp) :: KGE_dp_2d
    1018             : 
    1019             :     ! local variables
    1020             :     INTEGER(i4) :: n
    1021             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    1022           2 :     LOGICAL, DIMENSION(size(x, dim = 1), size(x, dim = 2)) :: maske
    1023           2 :     REAL(dp) :: mu_Obs, mu_Sim       ! Mean          of x and y
    1024           2 :     REAL(dp) :: sigma_Obs, sigma_Sim ! Standard dev. of x and y
    1025           2 :     REAL(dp) :: pearson_coor         ! Pearson Corr. of x and y
    1026             : 
    1027           2 :     if (present(mask)) then
    1028           3 :       shapemask = shape(mask)
    1029             :     else
    1030           3 :       shapemask = shape(x)
    1031             :     end if
    1032          12 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    1033           0 :             stop 'KGE_dp_2d: shapes of inputs(x,y) or mask are not matching'
    1034             :     !
    1035           2 :     if (present(mask)) then
    1036          61 :       maske = mask
    1037          61 :       n = count(maske)
    1038             :     else
    1039          61 :       maske = .true.
    1040           3 :       n = size(x)
    1041             :     end if
    1042           2 :     if (n .LE. 1_i4) stop 'KGE_dp_2d: sample size must be at least 2'
    1043             : 
    1044             :     ! Mean
    1045             :     mu_Obs = average(&
    1046             :             reshape(x(:, :), (/size(x, dim = 1) * size(x, dim = 2)/)), &
    1047           6 :             mask = reshape(maske(:, :), (/size(x, dim = 1) * size(x, dim = 2)/)))
    1048             :     mu_Sim = average(&
    1049             :             reshape(y(:, :), (/size(y, dim = 1) * size(y, dim = 2)/)), &
    1050           6 :             mask = reshape(maske(:, :), (/size(y, dim = 1) * size(y, dim = 2)/)))
    1051             :     ! Standard Deviation
    1052             :     sigma_Obs = stddev(&
    1053             :             reshape(x(:, :), (/size(x, dim = 1) * size(x, dim = 2)/)), &
    1054           6 :             mask = reshape(maske(:, :), (/size(x, dim = 1) * size(x, dim = 2)/)))
    1055             :     sigma_Sim = stddev(&
    1056             :             reshape(y(:, :), (/size(y, dim = 1) * size(y, dim = 2)/)), &
    1057           6 :             mask = reshape(maske(:, :), (/size(y, dim = 1) * size(y, dim = 2)/)))
    1058             :     ! Pearson product-moment correlation coefficient is with (N-1) not N
    1059             :     pearson_coor = correlation(&
    1060             :             reshape(x(:, :), (/size(x, dim = 1) * size(x, dim = 2)/)), &
    1061             :             reshape(y(:, :), (/size(y, dim = 1) * size(y, dim = 2)/)), &
    1062             :             mask = reshape(maske(:, :), (/size(y, dim = 1) * size(y, dim = 2)/))) * &
    1063           8 :             real(n, dp) / real(n - 1, dp)
    1064             :     !
    1065             :     KGE_dp_2d = 1.0 - SQRT(&
    1066             :             (1.0_dp - (mu_Sim / mu_Obs))**2 + &
    1067             :                     (1.0_dp - (sigma_Sim / sigma_Obs))**2 + &
    1068             :                     (1.0_dp - pearson_coor)**2             &
    1069           2 :             )
    1070             : 
    1071           2 :   END FUNCTION KGE_dp_2d
    1072             : 
    1073           4 :   FUNCTION KGE_dp_3d(x, y, mask)
    1074             : 
    1075           2 :     USE mo_moment, ONLY : average, stddev, correlation
    1076             : 
    1077             :     IMPLICIT NONE
    1078             : 
    1079             :     REAL(dp), DIMENSION(:, :, :), INTENT(IN) :: x, y
    1080             :     LOGICAL, DIMENSION(:, :, :), OPTIONAL, INTENT(IN) :: mask
    1081             :     REAL(dp) :: KGE_dp_3d
    1082             : 
    1083             :     ! local variables
    1084             :     INTEGER(i4) :: n
    1085             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    1086           2 :     LOGICAL, DIMENSION(size(x, dim = 1), size(x, dim = 2), size(x, dim = 3)) :: maske
    1087           2 :     REAL(dp) :: mu_Obs, mu_Sim       ! Mean          of x and y
    1088           2 :     REAL(dp) :: sigma_Obs, sigma_Sim ! Standard dev. of x and y
    1089           2 :     REAL(dp) :: pearson_coor         ! Pearson Corr. of x and y
    1090             : 
    1091           2 :     if (present(mask)) then
    1092           4 :       shapemask = shape(mask)
    1093             :     else
    1094           4 :       shapemask = shape(x)
    1095             :     end if
    1096          16 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    1097           0 :             stop 'KGE_dp_3d: shapes of inputs(x,y) or mask are not matching'
    1098             :     !
    1099           2 :     if (present(mask)) then
    1100         489 :       maske = mask
    1101         489 :       n = count(maske)
    1102             :     else
    1103         489 :       maske = .true.
    1104           4 :       n = size(x)
    1105             :     end if
    1106           2 :     if (n .LE. 1_i4) stop 'KGE_dp_3d: sample size must be at least 2'
    1107             : 
    1108             :     ! Mean
    1109             :     mu_Obs = average(&
    1110             :             reshape(x(:, :, :), (/size(x, dim = 1) * size(x, dim = 2) * size(x, dim = 3)/)), &
    1111           6 :             mask = reshape(maske(:, :, :), (/size(x, dim = 1) * size(x, dim = 2) * size(x, dim = 3)/)))
    1112             :     mu_Sim = average(&
    1113             :             reshape(y(:, :, :), (/size(y, dim = 1) * size(y, dim = 2) * size(y, dim = 3)/)), &
    1114           6 :             mask = reshape(maske(:, :, :), (/size(y, dim = 1) * size(y, dim = 2) * size(y, dim = 3)/)))
    1115             :     ! Standard Deviation
    1116             :     sigma_Obs = stddev(&
    1117             :             reshape(x(:, :, :), (/size(x, dim = 1) * size(x, dim = 2) * size(x, dim = 3)/)), &
    1118           6 :             mask = reshape(maske(:, :, :), (/size(x, dim = 1) * size(x, dim = 2) * size(x, dim = 3)/)))
    1119             :     sigma_Sim = stddev(&
    1120             :             reshape(y(:, :, :), (/size(y, dim = 1) * size(y, dim = 2) * size(y, dim = 3)/)), &
    1121           6 :             mask = reshape(maske(:, :, :), (/size(y, dim = 1) * size(y, dim = 2) * size(y, dim = 3)/)))
    1122             :     ! Pearson product-moment correlation coefficient is with (N-1) not N
    1123             :     pearson_coor = correlation(&
    1124             :             reshape(x(:, :, :), (/size(x, dim = 1) * size(x, dim = 2) * size(x, dim = 3)/)), &
    1125             :             reshape(y(:, :, :), (/size(y, dim = 1) * size(y, dim = 2) * size(y, dim = 3)/)), &
    1126             :             mask = reshape(maske(:, :, :), (/size(y, dim = 1) * size(y, dim = 2) * size(y, dim = 3)/))) * &
    1127           8 :             real(n, dp) / real(n - 1, dp)
    1128             :     !
    1129             :     KGE_dp_3d = 1.0 - SQRT(&
    1130             :             (1.0_dp - (mu_Sim / mu_Obs))**2 + &
    1131             :                     (1.0_dp - (sigma_Sim / sigma_Obs))**2 + &
    1132             :                     (1.0_dp - pearson_coor)**2             &
    1133           2 :             )
    1134             : 
    1135           2 :   END FUNCTION KGE_dp_3d
    1136             : 
    1137             :   ! ------------------------------------------------------------------
    1138             : 
    1139           4 :   FUNCTION KGEnocorr_sp_1d(x, y, mask)
    1140             : 
    1141           2 :     USE mo_moment, ONLY : average, stddev
    1142             : 
    1143             :     IMPLICIT NONE
    1144             : 
    1145             :     REAL(sp), DIMENSION(:), INTENT(IN) :: x, y
    1146             :     LOGICAL, DIMENSION(:), OPTIONAL, INTENT(IN) :: mask
    1147             :     REAL(sp) :: KGEnocorr_sp_1d
    1148             : 
    1149             :     ! local variables
    1150             :     INTEGER(i4) :: n
    1151             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    1152           2 :     LOGICAL, DIMENSION(size(x)) :: maske
    1153             : 
    1154           2 :     REAL(sp) :: mu_Obs, mu_Sim       ! Mean          of x and y
    1155           2 :     REAL(sp) :: sigma_Obs, sigma_Sim ! Standard dev. of x and y
    1156             : 
    1157           2 :     if (present(mask)) then
    1158           2 :       shapemask = shape(mask)
    1159             :     else
    1160           2 :       shapemask = shape(x)
    1161             :     end if
    1162           8 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    1163           0 :             stop 'KGEnocorr_sp_1d: shapes of inputs(x,y) or mask are not matching'
    1164             :     !
    1165           2 :     if (present(mask)) then
    1166          55 :       maske = mask
    1167          55 :       n = count(maske)
    1168             :     else
    1169          55 :       maske = .true.
    1170           1 :       n = size(x)
    1171             :     end if
    1172           2 :     if (n .LE. 1_i4) stop 'KGEnocorr_sp_1d: sample size must be at least 2'
    1173             : 
    1174             :     ! Mean
    1175           2 :     mu_Obs = average(x, mask = maske)
    1176           2 :     mu_Sim = average(y, mask = maske)
    1177             :     ! Standard Deviation
    1178           2 :     sigma_Obs = stddev(x, mask = maske)
    1179           2 :     sigma_Sim = stddev(y, mask = maske)
    1180             : 
    1181             :     !
    1182             :     KGEnocorr_sp_1d = 1.0 - SQRT(&
    1183             :             (1.0_sp - (mu_Sim / mu_Obs))**2 + &
    1184             :                     (1.0_sp - (sigma_Sim / sigma_Obs))**2 &
    1185           2 :             )
    1186             : 
    1187           2 :   END FUNCTION KGEnocorr_sp_1d
    1188             : 
    1189           4 :   FUNCTION KGEnocorr_sp_2d(x, y, mask)
    1190             : 
    1191           2 :     USE mo_moment, ONLY : average, stddev
    1192             : 
    1193             :     IMPLICIT NONE
    1194             : 
    1195             :     REAL(sp), DIMENSION(:, :), INTENT(IN) :: x, y
    1196             :     LOGICAL, DIMENSION(:, :), OPTIONAL, INTENT(IN) :: mask
    1197             :     REAL(sp) :: KGEnocorr_sp_2d
    1198             : 
    1199             :     ! local variables
    1200             :     INTEGER(i4) :: n
    1201             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    1202           2 :     LOGICAL, DIMENSION(size(x, dim = 1), size(x, dim = 2)) :: maske
    1203           2 :     REAL(sp) :: mu_Obs, mu_Sim       ! Mean          of x and y
    1204           2 :     REAL(sp) :: sigma_Obs, sigma_Sim ! Standard dev. of x and y
    1205             : 
    1206           2 :     if (present(mask)) then
    1207           3 :       shapemask = shape(mask)
    1208             :     else
    1209           3 :       shapemask = shape(x)
    1210             :     end if
    1211          12 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    1212           0 :             stop 'KGEnocorr_sp_2d: shapes of inputs(x,y) or mask are not matching'
    1213             :     !
    1214           2 :     if (present(mask)) then
    1215          61 :       maske = mask
    1216          61 :       n = count(maske)
    1217             :     else
    1218          61 :       maske = .true.
    1219           3 :       n = size(x)
    1220             :     end if
    1221           2 :     if (n .LE. 1_i4) stop 'KGEnocorr_sp_2d: sample size must be at least 2'
    1222             : 
    1223             :     ! Mean
    1224             :     mu_Obs = average(&
    1225             :             reshape(x(:, :), (/size(x, dim = 1) * size(x, dim = 2)/)), &
    1226           6 :             mask = reshape(maske(:, :), (/size(x, dim = 1) * size(x, dim = 2)/)))
    1227             :     mu_Sim = average(&
    1228             :             reshape(y(:, :), (/size(y, dim = 1) * size(y, dim = 2)/)), &
    1229           6 :             mask = reshape(maske(:, :), (/size(y, dim = 1) * size(y, dim = 2)/)))
    1230             :     ! Standard Deviation
    1231             :     sigma_Obs = stddev(&
    1232             :             reshape(x(:, :), (/size(x, dim = 1) * size(x, dim = 2)/)), &
    1233           6 :             mask = reshape(maske(:, :), (/size(x, dim = 1) * size(x, dim = 2)/)))
    1234             :     sigma_Sim = stddev(&
    1235             :             reshape(y(:, :), (/size(y, dim = 1) * size(y, dim = 2)/)), &
    1236           6 :             mask = reshape(maske(:, :), (/size(y, dim = 1) * size(y, dim = 2)/)))
    1237             :     !
    1238             :     KGEnocorr_sp_2d = 1.0 - SQRT(&
    1239             :             (1.0_sp - (mu_Sim / mu_Obs))**2 + &
    1240             :                     (1.0_sp - (sigma_Sim / sigma_Obs))**2 &
    1241           2 :             )
    1242             : 
    1243           2 :   END FUNCTION KGEnocorr_sp_2d
    1244             : 
    1245           4 :   FUNCTION KGEnocorr_sp_3d(x, y, mask)
    1246             : 
    1247           2 :     USE mo_moment, ONLY : average, stddev
    1248             : 
    1249             :     IMPLICIT NONE
    1250             : 
    1251             :     REAL(sp), DIMENSION(:, :, :), INTENT(IN) :: x, y
    1252             :     LOGICAL, DIMENSION(:, :, :), OPTIONAL, INTENT(IN) :: mask
    1253             :     REAL(sp) :: KGEnocorr_sp_3d
    1254             : 
    1255             :     ! local variables
    1256             :     INTEGER(i4) :: n
    1257             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    1258           2 :     LOGICAL, DIMENSION(size(x, dim = 1), size(x, dim = 2), size(x, dim = 3)) :: maske
    1259           2 :     REAL(sp) :: mu_Obs, mu_Sim       ! Mean          of x and y
    1260           2 :     REAL(sp) :: sigma_Obs, sigma_Sim ! Standard dev. of x and y
    1261             : 
    1262           2 :     if (present(mask)) then
    1263           4 :       shapemask = shape(mask)
    1264             :     else
    1265           4 :       shapemask = shape(x)
    1266             :     end if
    1267          16 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    1268           0 :             stop 'KGEnocorr_sp_3d: shapes of inputs(x,y) or mask are not matching'
    1269             :     !
    1270           2 :     if (present(mask)) then
    1271         489 :       maske = mask
    1272         489 :       n = count(maske)
    1273             :     else
    1274         489 :       maske = .true.
    1275           4 :       n = size(x)
    1276             :     end if
    1277           2 :     if (n .LE. 1_i4) stop 'KGEnocorr_sp_3d: sample size must be at least 2'
    1278             : 
    1279             :     ! Mean
    1280             :     mu_Obs = average(&
    1281             :             reshape(x(:, :, :), (/size(x, dim = 1) * size(x, dim = 2) * size(x, dim = 3)/)), &
    1282           6 :             mask = reshape(maske(:, :, :), (/size(x, dim = 1) * size(x, dim = 2) * size(x, dim = 3)/)))
    1283             :     mu_Sim = average(&
    1284             :             reshape(y(:, :, :), (/size(y, dim = 1) * size(y, dim = 2) * size(y, dim = 3)/)), &
    1285           6 :             mask = reshape(maske(:, :, :), (/size(y, dim = 1) * size(y, dim = 2) * size(y, dim = 3)/)))
    1286             :     ! Standard Deviation
    1287             :     sigma_Obs = stddev(&
    1288             :             reshape(x(:, :, :), (/size(x, dim = 1) * size(x, dim = 2) * size(x, dim = 3)/)), &
    1289           6 :             mask = reshape(maske(:, :, :), (/size(x, dim = 1) * size(x, dim = 2) * size(x, dim = 3)/)))
    1290             :     sigma_Sim = stddev(&
    1291             :             reshape(y(:, :, :), (/size(y, dim = 1) * size(y, dim = 2) * size(y, dim = 3)/)), &
    1292           6 :             mask = reshape(maske(:, :, :), (/size(y, dim = 1) * size(y, dim = 2) * size(y, dim = 3)/)))
    1293             : 
    1294             :     !
    1295             :     KGEnocorr_sp_3d = 1.0 - SQRT(&
    1296             :             (1.0_sp - (mu_Sim / mu_Obs))**2 + &
    1297             :                     (1.0_sp - (sigma_Sim / sigma_Obs))**2 &
    1298           2 :             )
    1299             : 
    1300           2 :   END FUNCTION KGEnocorr_sp_3d
    1301             : 
    1302           4 :   FUNCTION KGEnocorr_dp_1d(x, y, mask)
    1303             : 
    1304           2 :     USE mo_moment, ONLY : average, stddev
    1305             : 
    1306             :     IMPLICIT NONE
    1307             : 
    1308             :     REAL(dp), DIMENSION(:), INTENT(IN) :: x, y
    1309             :     LOGICAL, DIMENSION(:), OPTIONAL, INTENT(IN) :: mask
    1310             :     REAL(dp) :: KGEnocorr_dp_1d
    1311             : 
    1312             :     ! local variables
    1313             :     INTEGER(i4) :: n
    1314             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    1315           2 :     LOGICAL, DIMENSION(size(x)) :: maske
    1316             : 
    1317           2 :     REAL(dp) :: mu_Obs, mu_Sim       ! Mean          of x and y
    1318           2 :     REAL(dp) :: sigma_Obs, sigma_Sim ! Standard dev. of x and y
    1319             : 
    1320           2 :     if (present(mask)) then
    1321           2 :       shapemask = shape(mask)
    1322             :     else
    1323           2 :       shapemask = shape(x)
    1324             :     end if
    1325           8 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    1326           0 :             stop 'KGEnocorr_dp_1d: shapes of inputs(x,y) or mask are not matching'
    1327             :     !
    1328           2 :     if (present(mask)) then
    1329          55 :       maske = mask
    1330          55 :       n = count(maske)
    1331             :     else
    1332          55 :       maske = .true.
    1333           1 :       n = size(x)
    1334             :     end if
    1335           2 :     if (n .LE. 1_i4) stop 'KGEnocorr_dp_1d: sample size must be at least 2'
    1336             : 
    1337             :     ! Mean
    1338           2 :     mu_Obs = average(x, mask = maske)
    1339           2 :     mu_Sim = average(y, mask = maske)
    1340             :     ! Standard Deviation
    1341           2 :     sigma_Obs = stddev(x, mask = maske)
    1342           2 :     sigma_Sim = stddev(y, mask = maske)
    1343             : 
    1344             :     !
    1345             :     KGEnocorr_dp_1d = 1.0 - SQRT(&
    1346             :             (1.0_dp - (mu_Sim / mu_Obs))**2 + &
    1347             :                     (1.0_dp - (sigma_Sim / sigma_Obs))**2 &
    1348           2 :             )
    1349             : 
    1350           2 :   END FUNCTION KGEnocorr_dp_1d
    1351             : 
    1352           4 :   FUNCTION KGEnocorr_dp_2d(x, y, mask)
    1353             : 
    1354           2 :     USE mo_moment, ONLY : average, stddev
    1355             : 
    1356             :     IMPLICIT NONE
    1357             : 
    1358             :     REAL(dp), DIMENSION(:, :), INTENT(IN) :: x, y
    1359             :     LOGICAL, DIMENSION(:, :), OPTIONAL, INTENT(IN) :: mask
    1360             :     REAL(dp) :: KGEnocorr_dp_2d
    1361             : 
    1362             :     ! local variables
    1363             :     INTEGER(i4) :: n
    1364             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    1365           2 :     LOGICAL, DIMENSION(size(x, dim = 1), size(x, dim = 2)) :: maske
    1366           2 :     REAL(dp) :: mu_Obs, mu_Sim       ! Mean          of x and y
    1367           2 :     REAL(dp) :: sigma_Obs, sigma_Sim ! Standard dev. of x and y
    1368             : 
    1369           2 :     if (present(mask)) then
    1370           3 :       shapemask = shape(mask)
    1371             :     else
    1372           3 :       shapemask = shape(x)
    1373             :     end if
    1374          12 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    1375           0 :             stop 'KGEnocorr_dp_2d: shapes of inputs(x,y) or mask are not matching'
    1376             :     !
    1377           2 :     if (present(mask)) then
    1378          61 :       maske = mask
    1379          61 :       n = count(maske)
    1380             :     else
    1381          61 :       maske = .true.
    1382           3 :       n = size(x)
    1383             :     end if
    1384           2 :     if (n .LE. 1_i4) stop 'KGEnocorr_dp_2d: sample size must be at least 2'
    1385             : 
    1386             :     ! Mean
    1387             :     mu_Obs = average(&
    1388             :             reshape(x(:, :), (/size(x, dim = 1) * size(x, dim = 2)/)), &
    1389           6 :             mask = reshape(maske(:, :), (/size(x, dim = 1) * size(x, dim = 2)/)))
    1390             :     mu_Sim = average(&
    1391             :             reshape(y(:, :), (/size(y, dim = 1) * size(y, dim = 2)/)), &
    1392           6 :             mask = reshape(maske(:, :), (/size(y, dim = 1) * size(y, dim = 2)/)))
    1393             :     ! Standard Deviation
    1394             :     sigma_Obs = stddev(&
    1395             :             reshape(x(:, :), (/size(x, dim = 1) * size(x, dim = 2)/)), &
    1396           6 :             mask = reshape(maske(:, :), (/size(x, dim = 1) * size(x, dim = 2)/)))
    1397             :     sigma_Sim = stddev(&
    1398             :             reshape(y(:, :), (/size(y, dim = 1) * size(y, dim = 2)/)), &
    1399           6 :             mask = reshape(maske(:, :), (/size(y, dim = 1) * size(y, dim = 2)/)))
    1400             :     !
    1401             :     KGEnocorr_dp_2d = 1.0 - SQRT(&
    1402             :             (1.0_dp - (mu_Sim / mu_Obs))**2 + &
    1403             :                     (1.0_dp - (sigma_Sim / sigma_Obs))**2 &
    1404           2 :             )
    1405             : 
    1406           2 :   END FUNCTION KGEnocorr_dp_2d
    1407             : 
    1408           4 :   FUNCTION KGEnocorr_dp_3d(x, y, mask)
    1409             : 
    1410           2 :     USE mo_moment, ONLY : average, stddev
    1411             : 
    1412             :     IMPLICIT NONE
    1413             : 
    1414             :     REAL(dp), DIMENSION(:, :, :), INTENT(IN) :: x, y
    1415             :     LOGICAL, DIMENSION(:, :, :), OPTIONAL, INTENT(IN) :: mask
    1416             :     REAL(dp) :: KGEnocorr_dp_3d
    1417             : 
    1418             :     ! local variables
    1419             :     INTEGER(i4) :: n
    1420             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    1421           2 :     LOGICAL, DIMENSION(size(x, dim = 1), size(x, dim = 2), size(x, dim = 3)) :: maske
    1422           2 :     REAL(dp) :: mu_Obs, mu_Sim       ! Mean          of x and y
    1423           2 :     REAL(dp) :: sigma_Obs, sigma_Sim ! Standard dev. of x and y
    1424             : 
    1425           2 :     if (present(mask)) then
    1426           4 :       shapemask = shape(mask)
    1427             :     else
    1428           4 :       shapemask = shape(x)
    1429             :     end if
    1430          16 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    1431           0 :             stop 'KGEnocorr_dp_3d: shapes of inputs(x,y) or mask are not matching'
    1432             :     !
    1433           2 :     if (present(mask)) then
    1434         489 :       maske = mask
    1435         489 :       n = count(maske)
    1436             :     else
    1437         489 :       maske = .true.
    1438           4 :       n = size(x)
    1439             :     end if
    1440           2 :     if (n .LE. 1_i4) stop 'KGEnocorr_dp_3d: sample size must be at least 2'
    1441             : 
    1442             :     ! Mean
    1443             :     mu_Obs = average(&
    1444             :             reshape(x(:, :, :), (/size(x, dim = 1) * size(x, dim = 2) * size(x, dim = 3)/)), &
    1445           6 :             mask = reshape(maske(:, :, :), (/size(x, dim = 1) * size(x, dim = 2) * size(x, dim = 3)/)))
    1446             :     mu_Sim = average(&
    1447             :             reshape(y(:, :, :), (/size(y, dim = 1) * size(y, dim = 2) * size(y, dim = 3)/)), &
    1448           6 :             mask = reshape(maske(:, :, :), (/size(y, dim = 1) * size(y, dim = 2) * size(y, dim = 3)/)))
    1449             :     ! Standard Deviation
    1450             :     sigma_Obs = stddev(&
    1451             :             reshape(x(:, :, :), (/size(x, dim = 1) * size(x, dim = 2) * size(x, dim = 3)/)), &
    1452           6 :             mask = reshape(maske(:, :, :), (/size(x, dim = 1) * size(x, dim = 2) * size(x, dim = 3)/)))
    1453             :     sigma_Sim = stddev(&
    1454             :             reshape(y(:, :, :), (/size(y, dim = 1) * size(y, dim = 2) * size(y, dim = 3)/)), &
    1455           6 :             mask = reshape(maske(:, :, :), (/size(y, dim = 1) * size(y, dim = 2) * size(y, dim = 3)/)))
    1456             : 
    1457             :     !
    1458             :     KGEnocorr_dp_3d = 1.0 - SQRT(&
    1459             :             (1.0_dp - (mu_Sim / mu_Obs))**2 + &
    1460             :                     (1.0_dp - (sigma_Sim / sigma_Obs))**2 &
    1461           2 :             )
    1462             : 
    1463           2 :   END FUNCTION KGEnocorr_dp_3d
    1464             : 
    1465             : 
    1466             :   ! ------------------------------------------------------------------
    1467             : 
    1468           4 :   FUNCTION LNNSE_sp_1d(x, y, mask)
    1469             : 
    1470           2 :     USE mo_moment, ONLY : average
    1471             : 
    1472             :     IMPLICIT NONE
    1473             : 
    1474             :     REAL(sp), DIMENSION(:), INTENT(IN) :: x, y
    1475             :     LOGICAL, DIMENSION(:), OPTIONAL, INTENT(INOUT) :: mask
    1476             :     REAL(sp) :: LNNSE_sp_1d
    1477             : 
    1478             :     INTEGER(i4) :: n
    1479             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    1480         110 :     REAL(sp) :: xmean
    1481         436 :     REAL(sp), DIMENSION(size(x)) :: logx, logy, v1, v2
    1482           2 :     LOGICAL, DIMENSION(size(x)) :: maske
    1483             : 
    1484           2 :     if (present(mask)) then
    1485           2 :       shapemask = shape(mask)
    1486             :     else
    1487           2 :       shapemask = shape(x)
    1488             :     end if
    1489           8 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    1490           0 :             stop 'LNNSE_sp_1d: shapes of inputs(x,y) or mask are not matching'
    1491             :     !
    1492           2 :     if (present(mask)) then
    1493          55 :       maske = mask
    1494             :     else
    1495          55 :       maske = .true.
    1496             :     end if
    1497             : 
    1498             :     ! mask all negative and zero entries
    1499         110 :     where (x .lt. tiny(1.0_sp) .or. y .lt. tiny(1.0_sp))
    1500             :       maske = .false.
    1501             :     end where
    1502         110 :     n = count(maske)
    1503           2 :     if (n .LE. 1_i4) stop 'LNNSE_sp_1d: number of arguments must be at least 2'
    1504             : 
    1505             :     ! logarithms
    1506         110 :     logx = 0.0_sp
    1507         110 :     logy = 0.0_sp
    1508         326 :     where (maske)
    1509             :       logx = log(x)
    1510             :       logy = log(y)
    1511             :     end where
    1512             : 
    1513             :     ! mean of x
    1514           2 :     xmean = average(logx, mask = maske)
    1515             : 
    1516             :     ! NSE
    1517         110 :     v1 = merge(logy - logx, 0.0_sp, maske)
    1518         110 :     v2 = merge(logx - xmean, 0.0_sp, maske)
    1519         218 :     LNNSE_sp_1d = 1.0_sp - dot_product(v1, v1) / dot_product(v2, v2)
    1520             : 
    1521           2 :   END FUNCTION LNNSE_sp_1d
    1522             : 
    1523             :   ! ------------------------------------------------------------------
    1524             : 
    1525           4 :   FUNCTION LNNSE_dp_1d(x, y, mask)
    1526             : 
    1527           2 :     USE mo_moment, ONLY : average
    1528             : 
    1529             :     IMPLICIT NONE
    1530             : 
    1531             :     REAL(dp), DIMENSION(:), INTENT(IN) :: x, y
    1532             :     LOGICAL, DIMENSION(:), OPTIONAL, INTENT(INOUT) :: mask
    1533             :     REAL(dp) :: LNNSE_dp_1d
    1534             : 
    1535             :     INTEGER(i4) :: n
    1536             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    1537         110 :     REAL(dp) :: xmean
    1538         436 :     REAL(dp), DIMENSION(size(x)) :: logx, logy, v1, v2
    1539           2 :     LOGICAL, DIMENSION(size(x)) :: maske
    1540             : 
    1541           2 :     if (present(mask)) then
    1542           2 :       shapemask = shape(mask)
    1543             :     else
    1544           2 :       shapemask = shape(x)
    1545             :     end if
    1546           8 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    1547           0 :             stop 'LNNSE_dp_1d: shapes of inputs(x,y) or mask are not matching'
    1548             :     !
    1549           2 :     if (present(mask)) then
    1550          55 :       maske = mask
    1551             :     else
    1552          55 :       maske = .true.
    1553             :     end if
    1554             : 
    1555             :     ! mask all negative and zero entries
    1556         110 :     where (x .lt. tiny(1.0_dp) .or. y .lt. tiny(1.0_dp))
    1557             :       maske = .false.
    1558             :     end where
    1559         110 :     n = count(maske)
    1560           2 :     if (n .LE. 1_i4) stop 'LNNSE_dp_1d: number of arguments must be at least 2'
    1561             : 
    1562             :     ! logarithms
    1563         110 :     logx = 0.0_dp
    1564         110 :     logy = 0.0_dp
    1565         326 :     where (maske)
    1566             :       logx = log(x)
    1567             :       logy = log(y)
    1568             :     end where
    1569             : 
    1570             :     ! mean of x
    1571           2 :     xmean = average(logx, mask = maske)
    1572             : 
    1573             :     ! NSE
    1574         110 :     v1 = merge(logy - logx, 0.0_dp, maske)
    1575         110 :     v2 = merge(logx - xmean, 0.0_dp, maske)
    1576         218 :     LNNSE_dp_1d = 1.0_dp - sum(v1 * v1, mask = maske) / sum(v2 * v2, mask = maske)
    1577             : 
    1578           2 :   END FUNCTION LNNSE_dp_1d
    1579             : 
    1580             :   ! ------------------------------------------------------------------
    1581             : 
    1582           4 :   FUNCTION LNNSE_sp_2d(x, y, mask)
    1583             : 
    1584           2 :     USE mo_moment, ONLY : average
    1585             : 
    1586             :     IMPLICIT NONE
    1587             : 
    1588             :     REAL(sp), DIMENSION(:, :), INTENT(IN) :: x, y
    1589             :     LOGICAL, DIMENSION(:, :), OPTIONAL, INTENT(INOUT) :: mask
    1590             :     REAL(sp) :: LNNSE_sp_2d
    1591             : 
    1592             :     INTEGER(i4) :: n
    1593             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    1594           2 :     REAL(sp) :: xmean
    1595         484 :     REAL(sp), DIMENSION(size(x, dim = 1), size(x, dim = 2)) :: logx, logy, v1, v2
    1596           2 :     LOGICAL, DIMENSION(size(x, dim = 1), size(x, dim = 2)) :: maske
    1597             : 
    1598           2 :     if (present(mask)) then
    1599           3 :       shapemask = shape(mask)
    1600             :     else
    1601           3 :       shapemask = shape(x)
    1602             :     end if
    1603          12 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    1604           0 :             stop 'LNNSE_sp_2d: shapes of inputs(x,y) or mask are not matching'
    1605             :     !
    1606           2 :     if (present(mask)) then
    1607          61 :       maske = mask
    1608             :     else
    1609          61 :       maske = .true.
    1610             :     end if
    1611             : 
    1612             :     ! mask all negative and zero entries
    1613         122 :     where (x .lt. tiny(1.0_sp) .or. y .lt. tiny(1.0_sp))
    1614             :       maske = .false.
    1615             :     end where
    1616         122 :     n = count(maske)
    1617           2 :     if (n .LE. 1_i4) stop 'LNNSE_sp_2d: number of arguments must be at least 2'
    1618             : 
    1619             :     ! logarithms
    1620         122 :     logx = 0.0_sp
    1621         122 :     logy = 0.0_sp
    1622         362 :     where (maske)
    1623             :       logx = log(x)
    1624             :       logy = log(y)
    1625             :     end where
    1626             : 
    1627             :     ! mean of x
    1628           2 :     xmean = average(pack(logx, maske))
    1629             : 
    1630             :     ! NSE
    1631         122 :     v1 = merge(logy - logx, 0.0_sp, maske)
    1632         122 :     v2 = merge(logx - xmean, 0.0_sp, maske)
    1633         242 :     LNNSE_sp_2d = 1.0_sp - sum(v1 * v1, mask = maske) / sum(v2 * v2, mask = maske)
    1634             : 
    1635           2 :   END FUNCTION LNNSE_sp_2d
    1636             : 
    1637             :   ! ------------------------------------------------------------------
    1638             : 
    1639           4 :   FUNCTION LNNSE_dp_2d(x, y, mask)
    1640             : 
    1641           2 :     USE mo_moment, ONLY : average
    1642             : 
    1643             :     IMPLICIT NONE
    1644             : 
    1645             :     REAL(dp), DIMENSION(:, :), INTENT(IN) :: x, y
    1646             :     LOGICAL, DIMENSION(:, :), OPTIONAL, INTENT(INOUT) :: mask
    1647             :     REAL(dp) :: LNNSE_dp_2d
    1648             : 
    1649             :     INTEGER(i4) :: n
    1650             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    1651           2 :     REAL(dp) :: xmean
    1652         484 :     REAL(dp), DIMENSION(size(x, dim = 1), size(x, dim = 2)) :: logx, logy, v1, v2
    1653           2 :     LOGICAL, DIMENSION(size(x, dim = 1), size(x, dim = 2)) :: maske
    1654             : 
    1655           2 :     if (present(mask)) then
    1656           3 :       shapemask = shape(mask)
    1657             :     else
    1658           3 :       shapemask = shape(x)
    1659             :     end if
    1660          12 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    1661           0 :             stop 'LNNSE_dp_2d: shapes of inputs(x,y) or mask are not matching'
    1662             :     !
    1663           2 :     if (present(mask)) then
    1664          61 :       maske = mask
    1665             :     else
    1666          61 :       maske = .true.
    1667             :     end if
    1668             : 
    1669             :     ! mask all negative and zero entries
    1670         122 :     where (x .lt. tiny(1.0_dp) .or. y .lt. tiny(1.0_dp))
    1671             :       maske = .false.
    1672             :     end where
    1673         122 :     n = count(maske)
    1674           2 :     if (n .LE. 1_i4) stop 'LNNSE_dp_2d: number of arguments must be at least 2'
    1675             : 
    1676             :     ! logarithms
    1677         122 :     logx = 0.0_dp
    1678         122 :     logy = 0.0_dp
    1679         362 :     where (maske)
    1680             :       logx = log(x)
    1681             :       logy = log(y)
    1682             :     end where
    1683             : 
    1684             :     ! mean of x
    1685           2 :     xmean = average(pack(logx, maske))
    1686             : 
    1687             :     ! NSE
    1688         122 :     v1 = merge(logy - logx, 0.0_dp, maske)
    1689         122 :     v2 = merge(logx - xmean, 0.0_dp, maske)
    1690         242 :     LNNSE_dp_2d = 1.0_dp - sum(v1 * v1, mask = maske) / sum(v2 * v2, mask = maske)
    1691             : 
    1692           2 :   END FUNCTION LNNSE_dp_2d
    1693             : 
    1694             :   ! ------------------------------------------------------------------
    1695             : 
    1696           4 :   FUNCTION LNNSE_sp_3d(x, y, mask)
    1697             : 
    1698           2 :     USE mo_moment, ONLY : average
    1699             : 
    1700             :     IMPLICIT NONE
    1701             : 
    1702             :     REAL(sp), DIMENSION(:, :, :), INTENT(IN) :: x, y
    1703             :     LOGICAL, DIMENSION(:, :, :), OPTIONAL, INTENT(INOUT) :: mask
    1704             :     REAL(sp) :: LNNSE_sp_3d
    1705             : 
    1706             :     INTEGER(i4) :: n
    1707             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    1708           2 :     REAL(sp) :: xmean
    1709        3908 :     REAL(sp), DIMENSION(size(x, dim = 1), size(x, dim = 2), size(x, dim = 3)) :: logx, logy, v1, v2
    1710           2 :     LOGICAL, DIMENSION(size(x, dim = 1), size(x, dim = 2), size(x, dim = 3)) :: maske
    1711             : 
    1712           2 :     if (present(mask)) then
    1713           4 :       shapemask = shape(mask)
    1714             :     else
    1715           4 :       shapemask = shape(x)
    1716             :     end if
    1717          16 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    1718           0 :             stop 'LNNSE_sp_3d: shapes of inputs(x,y) or mask are not matching'
    1719             :     !
    1720           2 :     if (present(mask)) then
    1721         489 :       maske = mask
    1722             :     else
    1723         489 :       maske = .true.
    1724             :     end if
    1725             : 
    1726             :     ! mask all negative and zero entries
    1727         978 :     where (x .lt. tiny(1.0_sp) .or. y .lt. tiny(1.0_sp))
    1728             :       maske = .false.
    1729             :     end where
    1730         978 :     n = count(maske)
    1731           2 :     if (n .LE. 1_i4) stop 'LNNSE_sp_3d: number of arguments must be at least 2'
    1732             : 
    1733             :     ! logarithms
    1734         978 :     logx = 0.0_sp
    1735         978 :     logy = 0.0_sp
    1736        2930 :     where (maske)
    1737             :       logx = log(x)
    1738             :       logy = log(y)
    1739             :     end where
    1740             : 
    1741             :     ! mean of x
    1742           2 :     xmean = average(pack(logx, maske))
    1743             : 
    1744             :     ! NSE
    1745         978 :     v1 = merge(logy - logx, 0.0_sp, maske)
    1746         978 :     v2 = merge(logx - xmean, 0.0_sp, maske)
    1747        1954 :     LNNSE_sp_3d = 1.0_sp - sum(v1 * v1, mask = maske) / sum(v2 * v2, mask = maske)
    1748             : 
    1749           2 :   END FUNCTION LNNSE_sp_3d
    1750             : 
    1751             :   ! ------------------------------------------------------------------
    1752             : 
    1753           4 :   FUNCTION LNNSE_dp_3d(x, y, mask)
    1754             : 
    1755           2 :     USE mo_moment, ONLY : average
    1756             : 
    1757             :     IMPLICIT NONE
    1758             : 
    1759             :     REAL(dp), DIMENSION(:, :, :), INTENT(IN) :: x, y
    1760             :     LOGICAL, DIMENSION(:, :, :), OPTIONAL, INTENT(INOUT) :: mask
    1761             :     REAL(dp) :: LNNSE_dp_3d
    1762             : 
    1763             :     INTEGER(i4) :: n
    1764             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    1765           2 :     REAL(dp) :: xmean
    1766        3908 :     REAL(dp), DIMENSION(size(x, dim = 1), size(x, dim = 2), size(x, dim = 3)) :: logx, logy, v1, v2
    1767           2 :     LOGICAL, DIMENSION(size(x, dim = 1), size(x, dim = 2), size(x, dim = 3)) :: maske
    1768             : 
    1769           2 :     if (present(mask)) then
    1770           4 :       shapemask = shape(mask)
    1771             :     else
    1772           4 :       shapemask = shape(x)
    1773             :     end if
    1774          16 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    1775           0 :             stop 'LNNSE_dp_3d: shapes of inputs(x,y) or mask are not matching'
    1776             :     !
    1777           2 :     if (present(mask)) then
    1778         489 :       maske = mask
    1779             :     else
    1780         489 :       maske = .true.
    1781             :     end if
    1782             : 
    1783             :     ! mask all negative and zero entries
    1784         978 :     where (x .lt. tiny(1.0_dp) .or. y .lt. tiny(1.0_dp))
    1785             :       maske = .false.
    1786             :     end where
    1787         978 :     n = count(maske)
    1788           2 :     if (n .LE. 1_i4) stop 'LNNSE_dp_3d: number of arguments must be at least 2'
    1789             : 
    1790             :     ! logarithms
    1791         978 :     logx = 0.0_dp
    1792         978 :     logy = 0.0_dp
    1793        2930 :     where (maske)
    1794             :       logx = log(x)
    1795             :       logy = log(y)
    1796             :     end where
    1797             : 
    1798             :     ! mean of x
    1799           2 :     xmean = average(pack(logx, maske))
    1800             : 
    1801             :     ! NSE
    1802         978 :     v1 = merge(logy - logx, 0.0_dp, maske)
    1803         978 :     v2 = merge(logx - xmean, 0.0_dp, maske)
    1804        1954 :     LNNSE_dp_3d = 1.0_dp - sum(v1 * v1, mask = maske) / sum(v2 * v2, mask = maske)
    1805             : 
    1806           2 :   END FUNCTION LNNSE_dp_3d
    1807             : 
    1808             :   ! ------------------------------------------------------------------
    1809             : 
    1810           4 :   FUNCTION MAE_sp_1d(x, y, mask)
    1811             : 
    1812             :     IMPLICIT NONE
    1813             : 
    1814             :     REAL(sp), DIMENSION(:), INTENT(IN) :: x, y
    1815             :     LOGICAL, DIMENSION(:), OPTIONAL, INTENT(IN) :: mask
    1816             :     REAL(sp) :: MAE_sp_1d
    1817             : 
    1818             :     INTEGER(i4) :: n
    1819             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    1820           2 :     LOGICAL, DIMENSION(size(x, dim = 1)) :: maske
    1821             : 
    1822           2 :     if (present(mask)) then
    1823           2 :       shapemask = shape(mask)
    1824             :     else
    1825           2 :       shapemask = shape(x)
    1826             :     end if
    1827           8 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    1828           0 :             stop 'MAE_sp_1d: shapes of inputs(x,y) or mask are not matching'
    1829             :     !
    1830           2 :     if (present(mask)) then
    1831          56 :       maske = mask
    1832          55 :       n = count(maske)
    1833             :     else
    1834          55 :       maske = .true.
    1835           1 :       n = size(x, dim = 1)
    1836             :     end if
    1837           2 :     if (n .LE. 1_i4) stop 'MAE_sp_1d: number of arguments must be at least 2'
    1838             :     !
    1839           2 :     MAE_sp_1d = SAE_sp_1d(x, y, mask = maske) / real(n, sp)
    1840             : 
    1841           2 :   END FUNCTION MAE_sp_1d
    1842             : 
    1843           4 :   FUNCTION MAE_dp_1d(x, y, mask)
    1844             : 
    1845             :     IMPLICIT NONE
    1846             : 
    1847             :     REAL(dp), DIMENSION(:), INTENT(IN) :: x, y
    1848             :     LOGICAL, DIMENSION(:), OPTIONAL, INTENT(IN) :: mask
    1849             :     REAL(dp) :: MAE_dp_1d
    1850             : 
    1851             :     INTEGER(i4) :: n
    1852             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    1853           2 :     LOGICAL, DIMENSION(size(x, dim = 1)) :: maske
    1854             : 
    1855           2 :     if (present(mask)) then
    1856           2 :       shapemask = shape(mask)
    1857             :     else
    1858           2 :       shapemask = shape(x)
    1859             :     end if
    1860           8 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    1861           0 :             stop 'MAE_dp_1d: shapes of inputs(x,y) or mask are not matching'
    1862             :     !
    1863           2 :     if (present(mask)) then
    1864          56 :       maske = mask
    1865          55 :       n = count(maske)
    1866             :     else
    1867          55 :       maske = .true.
    1868           1 :       n = size(x, dim = 1)
    1869             :     end if
    1870           2 :     if (n .LE. 1_i4) stop 'MAE_dp_1d: number of arguments must be at least 2'
    1871             :     !
    1872           2 :     MAE_dp_1d = SAE_dp_1d(x, y, mask = maske) / real(n, dp)
    1873             : 
    1874           2 :   END FUNCTION MAE_dp_1d
    1875             : 
    1876           4 :   FUNCTION MAE_sp_2d(x, y, mask)
    1877             : 
    1878             :     IMPLICIT NONE
    1879             : 
    1880             :     REAL(sp), DIMENSION(:, :), INTENT(IN) :: x, y
    1881             :     LOGICAL, DIMENSION(:, :), OPTIONAL, INTENT(IN) :: mask
    1882             :     REAL(sp) :: MAE_sp_2d
    1883             : 
    1884             :     INTEGER(i4) :: n
    1885             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    1886           2 :     LOGICAL, DIMENSION(size(x, dim = 1), size(x, dim = 2)) :: maske
    1887             : 
    1888           2 :     if (present(mask)) then
    1889           3 :       shapemask = shape(mask)
    1890             :     else
    1891           3 :       shapemask = shape(x)
    1892             :     end if
    1893          12 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    1894           0 :             stop 'MAE_sp_2d: shapes of inputs(x,y) or mask are not matching'
    1895             :     !
    1896           2 :     if (present(mask)) then
    1897          62 :       maske = mask
    1898          61 :       n = count(maske)
    1899             :     else
    1900          61 :       maske = .true.
    1901           1 :       n = size(x, dim = 1) * size(x, dim = 2)
    1902             :     end if
    1903           2 :     if (n .LE. 1_i4) stop 'MAE_sp_2d: number of arguments must be at least 2'
    1904             :     !
    1905           2 :     MAE_sp_2d = SAE_sp_2d(x, y, mask = maske) / real(n, sp)
    1906             : 
    1907           2 :   END FUNCTION MAE_sp_2d
    1908             : 
    1909           4 :   FUNCTION MAE_dp_2d(x, y, mask)
    1910             : 
    1911             :     IMPLICIT NONE
    1912             : 
    1913             :     REAL(dp), DIMENSION(:, :), INTENT(IN) :: x, y
    1914             :     LOGICAL, DIMENSION(:, :), OPTIONAL, INTENT(IN) :: mask
    1915             :     REAL(dp) :: MAE_dp_2d
    1916             : 
    1917             :     INTEGER(i4) :: n
    1918             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    1919           2 :     LOGICAL, DIMENSION(size(x, dim = 1), size(x, dim = 2)) :: maske
    1920             : 
    1921           2 :     if (present(mask)) then
    1922           3 :       shapemask = shape(mask)
    1923             :     else
    1924           3 :       shapemask = shape(x)
    1925             :     end if
    1926          12 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    1927           0 :             stop 'MAE_dp_2d: shapes of inputs(x,y) or mask are not matching'
    1928             :     !
    1929           2 :     if (present(mask)) then
    1930          62 :       maske = mask
    1931          61 :       n = count(maske)
    1932             :     else
    1933          61 :       maske = .true.
    1934           1 :       n = size(x, dim = 1) * size(x, dim = 2)
    1935             :     end if
    1936           2 :     if (n .LE. 1_i4) stop 'MAE_dp_2d: number of arguments must be at least 2'
    1937             :     !
    1938           2 :     MAE_dp_2d = SAE_dp_2d(x, y, mask = maske) / real(n, dp)
    1939             : 
    1940           2 :   END FUNCTION MAE_dp_2d
    1941             : 
    1942           4 :   FUNCTION MAE_sp_3d(x, y, mask)
    1943             : 
    1944             :     IMPLICIT NONE
    1945             : 
    1946             :     REAL(sp), DIMENSION(:, :, :), INTENT(IN) :: x, y
    1947             :     LOGICAL, DIMENSION(:, :, :), OPTIONAL, INTENT(IN) :: mask
    1948             :     REAL(sp) :: MAE_sp_3d
    1949             : 
    1950             :     INTEGER(i4) :: n
    1951             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    1952             :     LOGICAL, DIMENSION(size(x, dim = 1), size(x, dim = 2), &
    1953           2 :             size(x, dim = 3)) :: maske
    1954             : 
    1955           2 :     if (present(mask)) then
    1956           4 :       shapemask = shape(mask)
    1957             :     else
    1958           4 :       shapemask = shape(x)
    1959             :     end if
    1960          16 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    1961           0 :             stop 'MAE_sp_3d: shapes of inputs(x,y) or mask are not matching'
    1962             :     !
    1963           2 :     if (present(mask)) then
    1964         490 :       maske = mask
    1965         489 :       n = count(maske)
    1966             :     else
    1967         489 :       maske = .true.
    1968           1 :       n = size(x, dim = 1) * size(x, dim = 2) * size(x, dim = 3)
    1969             :     end if
    1970           2 :     if (n .LE. 1_i4) stop 'MAE_sp_3d: number of arguments must be at least 2'
    1971             :     !
    1972           2 :     MAE_sp_3d = SAE_sp_3d(x, y, mask = maske) / real(n, sp)
    1973             : 
    1974           2 :   END FUNCTION MAE_sp_3d
    1975             : 
    1976           4 :   FUNCTION MAE_dp_3d(x, y, mask)
    1977             : 
    1978             :     IMPLICIT NONE
    1979             : 
    1980             :     REAL(dp), DIMENSION(:, :, :), INTENT(IN) :: x, y
    1981             :     LOGICAL, DIMENSION(:, :, :), OPTIONAL, INTENT(IN) :: mask
    1982             :     REAL(dp) :: MAE_dp_3d
    1983             : 
    1984             :     INTEGER(i4) :: n
    1985             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    1986             :     LOGICAL, DIMENSION(size(x, dim = 1), size(x, dim = 2), &
    1987           2 :             size(x, dim = 3)) :: maske
    1988             : 
    1989           2 :     if (present(mask)) then
    1990           4 :       shapemask = shape(mask)
    1991             :     else
    1992           4 :       shapemask = shape(x)
    1993             :     end if
    1994          16 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    1995           0 :             stop 'MAE_dp_3d: shapes of inputs(x,y) or mask are not matching'
    1996             :     !
    1997           2 :     if (present(mask)) then
    1998         490 :       maske = mask
    1999         489 :       n = count(maske)
    2000             :     else
    2001         489 :       maske = .true.
    2002           1 :       n = size(x, dim = 1) * size(x, dim = 2) * size(x, dim = 3)
    2003             :     end if
    2004           2 :     if (n .LE. 1_i4) stop 'MAE_dp_3d: number of arguments must be at least 2'
    2005             :     !
    2006           2 :     MAE_dp_3d = SAE_dp_3d(x, y, mask = maske) / real(n, dp)
    2007             : 
    2008           2 :   END FUNCTION MAE_dp_3d
    2009             : 
    2010             :   ! ------------------------------------------------------------------
    2011             : 
    2012           8 :   FUNCTION MSE_sp_1d(x, y, mask)
    2013             : 
    2014             :     IMPLICIT NONE
    2015             : 
    2016             :     REAL(sp), DIMENSION(:), INTENT(IN) :: x, y
    2017             :     LOGICAL, DIMENSION(:), OPTIONAL, INTENT(IN) :: mask
    2018             :     REAL(sp) :: MSE_sp_1d
    2019             : 
    2020             :     INTEGER(i4) :: n
    2021             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    2022           4 :     LOGICAL, DIMENSION(size(x, dim = 1)) :: maske
    2023             : 
    2024           4 :     if (present(mask)) then
    2025           6 :       shapemask = shape(mask)
    2026             :     else
    2027           2 :       shapemask = shape(x)
    2028             :     end if
    2029          16 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    2030           0 :             stop 'MSE_sp_1d: shapes of inputs(x,y) or mask are not matching'
    2031             :     !
    2032           4 :     if (present(mask)) then
    2033         168 :       maske = mask
    2034         165 :       n = count(maske)
    2035             :     else
    2036          55 :       maske = .true.
    2037           1 :       n = size(x, dim = 1)
    2038             :     end if
    2039           4 :     if (n .LE. 1_i4) stop 'MSE_sp_1d: number of arguments must be at least 2'
    2040             :     !
    2041           4 :     MSE_sp_1d = SSE_sp_1d(x, y, mask = maske) / real(n, sp)
    2042             : 
    2043           2 :   END FUNCTION MSE_sp_1d
    2044             : 
    2045           8 :   FUNCTION MSE_dp_1d(x, y, mask)
    2046             : 
    2047             :     IMPLICIT NONE
    2048             : 
    2049             :     REAL(dp), DIMENSION(:), INTENT(IN) :: x, y
    2050             :     LOGICAL, DIMENSION(:), OPTIONAL, INTENT(IN) :: mask
    2051             :     REAL(dp) :: MSE_dp_1d
    2052             : 
    2053             :     INTEGER(i4) :: n
    2054             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    2055           4 :     LOGICAL, DIMENSION(size(x, dim = 1)) :: maske
    2056             : 
    2057           4 :     if (present(mask)) then
    2058           6 :       shapemask = shape(mask)
    2059             :     else
    2060           2 :       shapemask = shape(x)
    2061             :     end if
    2062          16 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    2063           0 :             stop 'MSE_dp_1d: shapes of inputs(x,y) or mask are not matching'
    2064             :     !
    2065           4 :     if (present(mask)) then
    2066         168 :       maske = mask
    2067         165 :       n = count(maske)
    2068             :     else
    2069          55 :       maske = .true.
    2070           1 :       n = size(x, dim = 1)
    2071             :     end if
    2072           4 :     if (n .LE. 1_i4) stop 'MSE_dp_1d: number of arguments must be at least 2'
    2073             :     !
    2074           4 :     MSE_dp_1d = SSE_dp_1d(x, y, mask = maske) / real(n, dp)
    2075             : 
    2076           4 :   END FUNCTION MSE_dp_1d
    2077             : 
    2078           8 :   FUNCTION MSE_sp_2d(x, y, mask)
    2079             : 
    2080             :     IMPLICIT NONE
    2081             : 
    2082             :     REAL(sp), DIMENSION(:, :), INTENT(IN) :: x, y
    2083             :     LOGICAL, DIMENSION(:, :), OPTIONAL, INTENT(IN) :: mask
    2084             :     REAL(sp) :: MSE_sp_2d
    2085             : 
    2086             :     INTEGER(i4) :: n
    2087             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    2088           4 :     LOGICAL, DIMENSION(size(x, dim = 1), size(x, dim = 2)) :: maske
    2089             : 
    2090           4 :     if (present(mask)) then
    2091           9 :       shapemask = shape(mask)
    2092             :     else
    2093           3 :       shapemask = shape(x)
    2094             :     end if
    2095          24 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    2096           0 :             stop 'MSE_sp_2d: shapes of inputs(x,y) or mask are not matching'
    2097             :     !
    2098           4 :     if (present(mask)) then
    2099         186 :       maske = mask
    2100         183 :       n = count(maske)
    2101             :     else
    2102          61 :       maske = .true.
    2103           1 :       n = size(x, dim = 1) * size(x, dim = 2)
    2104             :     end if
    2105           4 :     if (n .LE. 1_i4) stop 'MSE_sp_2d: number of arguments must be at least 2'
    2106             :     !
    2107           4 :     MSE_sp_2d = SSE_sp_2d(x, y, mask = maske) / real(n, sp)
    2108             : 
    2109           4 :   END FUNCTION MSE_sp_2d
    2110             : 
    2111           8 :   FUNCTION MSE_dp_2d(x, y, mask)
    2112             : 
    2113             :     IMPLICIT NONE
    2114             : 
    2115             :     REAL(dp), DIMENSION(:, :), INTENT(IN) :: x, y
    2116             :     LOGICAL, DIMENSION(:, :), OPTIONAL, INTENT(IN) :: mask
    2117             :     REAL(dp) :: MSE_dp_2d
    2118             : 
    2119             :     INTEGER(i4) :: n
    2120             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    2121           4 :     LOGICAL, DIMENSION(size(x, dim = 1), size(x, dim = 2)) :: maske
    2122             : 
    2123           4 :     if (present(mask)) then
    2124           9 :       shapemask = shape(mask)
    2125             :     else
    2126           3 :       shapemask = shape(x)
    2127             :     end if
    2128          24 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    2129           0 :             stop 'MSE_dp_2d: shapes of inputs(x,y) or mask are not matching'
    2130             :     !
    2131           4 :     if (present(mask)) then
    2132         186 :       maske = mask
    2133         183 :       n = count(maske)
    2134             :     else
    2135          61 :       maske = .true.
    2136           1 :       n = size(x, dim = 1) * size(x, dim = 2)
    2137             :     end if
    2138           4 :     if (n .LE. 1_i4) stop 'MSE_dp_2d: number of arguments must be at least 2'
    2139             :     !
    2140           4 :     MSE_dp_2d = SSE_dp_2d(x, y, mask = maske) / real(n, dp)
    2141             : 
    2142           4 :   END FUNCTION MSE_dp_2d
    2143             : 
    2144           8 :   FUNCTION MSE_sp_3d(x, y, mask)
    2145             : 
    2146             :     IMPLICIT NONE
    2147             : 
    2148             :     REAL(sp), DIMENSION(:, :, :), INTENT(IN) :: x, y
    2149             :     LOGICAL, DIMENSION(:, :, :), OPTIONAL, INTENT(IN) :: mask
    2150             :     REAL(sp) :: MSE_sp_3d
    2151             : 
    2152             :     INTEGER(i4) :: n
    2153             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    2154             :     LOGICAL, DIMENSION(size(x, dim = 1), size(x, dim = 2), &
    2155           4 :             size(x, dim = 3)) :: maske
    2156             : 
    2157           4 :     if (present(mask)) then
    2158          12 :       shapemask = shape(mask)
    2159             :     else
    2160           4 :       shapemask = shape(x)
    2161             :     end if
    2162          32 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    2163           0 :             stop 'MSE_sp_3d: shapes of inputs(x,y) or mask are not matching'
    2164             :     !
    2165           4 :     if (present(mask)) then
    2166        1470 :       maske = mask
    2167        1467 :       n = count(maske)
    2168             :     else
    2169         489 :       maske = .true.
    2170           1 :       n = size(x, dim = 1) * size(x, dim = 2) * size(x, dim = 3)
    2171             :     end if
    2172           4 :     if (n .LE. 1_i4) stop 'MSE_sp_3d: number of arguments must be at least 2'
    2173             :     !
    2174           4 :     MSE_sp_3d = SSE_sp_3d(x, y, mask = maske) / real(n, sp)
    2175             : 
    2176           4 :   END FUNCTION MSE_sp_3d
    2177             : 
    2178           8 :   FUNCTION MSE_dp_3d(x, y, mask)
    2179             : 
    2180             :     IMPLICIT NONE
    2181             : 
    2182             :     REAL(dp), DIMENSION(:, :, :), INTENT(IN) :: x, y
    2183             :     LOGICAL, DIMENSION(:, :, :), OPTIONAL, INTENT(IN) :: mask
    2184             :     REAL(dp) :: MSE_dp_3d
    2185             : 
    2186             :     INTEGER(i4) :: n
    2187             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    2188             :     LOGICAL, DIMENSION(size(x, dim = 1), size(x, dim = 2), &
    2189           4 :             size(x, dim = 3)) :: maske
    2190             : 
    2191           4 :     if (present(mask)) then
    2192          12 :       shapemask = shape(mask)
    2193             :     else
    2194           4 :       shapemask = shape(x)
    2195             :     end if
    2196          32 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    2197           0 :             stop 'MSE_dp_3d: shapes of inputs(x,y) or mask are not matching'
    2198             :     !
    2199           4 :     if (present(mask)) then
    2200        1470 :       maske = mask
    2201        1467 :       n = count(maske)
    2202             :     else
    2203         489 :       maske = .true.
    2204           1 :       n = size(x, dim = 1) * size(x, dim = 2) * size(x, dim = 3)
    2205             :     end if
    2206           4 :     if (n .LE. 1_i4) stop 'MSE_dp_3d: number of arguments must be at least 2'
    2207             :     !
    2208           4 :     MSE_dp_3d = SSE_dp_3d(x, y, mask = maske) / real(n, dp)
    2209             : 
    2210           4 :   END FUNCTION MSE_dp_3d
    2211             : 
    2212             :   ! ------------------------------------------------------------------
    2213             : 
    2214           4 :   FUNCTION NSE_sp_1d(x, y, mask)
    2215             : 
    2216           4 :     USE mo_moment, ONLY : average
    2217             : 
    2218             :     IMPLICIT NONE
    2219             : 
    2220             :     REAL(sp), DIMENSION(:), INTENT(IN) :: x, y
    2221             :     LOGICAL, DIMENSION(:), OPTIONAL, INTENT(IN) :: mask
    2222             :     REAL(sp) :: NSE_sp_1d
    2223             : 
    2224             :     INTEGER(i4) :: n
    2225             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    2226         110 :     REAL(sp) :: xmean
    2227         222 :     REAL(sp), DIMENSION(size(x)) :: v1, v2
    2228           2 :     LOGICAL, DIMENSION(size(x)) :: maske
    2229             : 
    2230           2 :     if (present(mask)) then
    2231           2 :       shapemask = shape(mask)
    2232             :     else
    2233           2 :       shapemask = shape(x)
    2234             :     end if
    2235           8 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    2236           0 :             stop 'NSE_sp_1d: shapes of inputs(x,y) or mask are not matching'
    2237             :     !
    2238           2 :     if (present(mask)) then
    2239          56 :       maske = mask
    2240          55 :       n = count(maske)
    2241             :     else
    2242          55 :       maske = .true.
    2243           1 :       n = size(x)
    2244             :     end if
    2245           2 :     if (n .LE. 1_i4) stop 'NSE_sp_1d: number of arguments must be at least 2'
    2246             :     ! mean of x
    2247           2 :     xmean = average(x, mask = maske)
    2248             :     !
    2249         112 :     v1 = merge(y - x, 0.0_sp, maske)
    2250         112 :     v2 = merge(x - xmean, 0.0_sp, maske)
    2251             :     !
    2252         220 :     NSE_sp_1d = 1.0_sp - dot_product(v1, v1) / dot_product(v2, v2)
    2253             : 
    2254           2 :   END FUNCTION NSE_sp_1d
    2255             : 
    2256           4 :   FUNCTION NSE_dp_1d(x, y, mask)
    2257             : 
    2258           2 :     USE mo_moment, ONLY : average
    2259             : 
    2260             :     IMPLICIT NONE
    2261             : 
    2262             :     REAL(dp), DIMENSION(:), INTENT(IN) :: x, y
    2263             :     LOGICAL, DIMENSION(:), OPTIONAL, INTENT(IN) :: mask
    2264             :     REAL(dp) :: NSE_dp_1d
    2265             : 
    2266             :     INTEGER(i4) :: n
    2267             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    2268         110 :     REAL(dp) :: xmean
    2269         222 :     REAL(dp), DIMENSION(size(x)) :: v1, v2
    2270           2 :     LOGICAL, DIMENSION(size(x)) :: maske
    2271             : 
    2272           2 :     if (present(mask)) then
    2273           2 :       shapemask = shape(mask)
    2274             :     else
    2275           2 :       shapemask = shape(x)
    2276             :     end if
    2277           8 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    2278           0 :             stop 'NSE_dp_1d: shapes of inputs(x,y) or mask are not matching'
    2279             :     !
    2280           2 :     if (present(mask)) then
    2281          56 :       maske = mask
    2282          55 :       n = count(maske)
    2283             :     else
    2284          55 :       maske = .true.
    2285           1 :       n = size(x)
    2286             :     end if
    2287           2 :     if (n .LE. 1_i4) stop 'NSE_dp_1d: number of arguments must be at least 2'
    2288             :     ! mean of x
    2289           2 :     xmean = average(x, mask = maske)
    2290             :     !
    2291         112 :     v1 = merge(y - x, 0.0_dp, maske)
    2292         112 :     v2 = merge(x - xmean, 0.0_dp, maske)
    2293             :     !
    2294         220 :     NSE_dp_1d = 1.0_dp - dot_product(v1, v1) / dot_product(v2, v2)
    2295             : 
    2296           2 :   END FUNCTION NSE_dp_1d
    2297             : 
    2298           4 :   FUNCTION NSE_sp_2d(x, y, mask)
    2299             : 
    2300           2 :     USE mo_moment, ONLY : average
    2301             : 
    2302             :     IMPLICIT NONE
    2303             : 
    2304             :     REAL(sp), DIMENSION(:, :), INTENT(IN) :: x, y
    2305             :     LOGICAL, DIMENSION(:, :), OPTIONAL, INTENT(IN) :: mask
    2306             :     REAL(sp) :: NSE_sp_2d
    2307             : 
    2308             :     INTEGER(i4) :: n
    2309             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    2310           2 :     REAL(sp) :: xmean
    2311           2 :     LOGICAL, DIMENSION(size(x, dim = 1), size(x, dim = 2)) :: maske
    2312             : 
    2313           2 :     if (present(mask)) then
    2314           3 :       shapemask = shape(mask)
    2315             :     else
    2316           3 :       shapemask = shape(x)
    2317             :     end if
    2318          12 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    2319           0 :             stop 'NSE_sp_2d: shapes of inputs(x,y) or mask are not matching'
    2320             :     !
    2321           2 :     if (present(mask)) then
    2322          62 :       maske = mask
    2323          61 :       n = count(maske)
    2324             :     else
    2325          61 :       maske = .true.
    2326           1 :       n = size(x, dim = 1) * size(x, dim = 2)
    2327             :     end if
    2328             :     !
    2329           2 :     if (n .LE. 1_i4) stop 'NSE_sp_2d: number of arguments must be at least 2'
    2330             :     ! mean of x
    2331             :     xmean = average(reshape(x(:, :), (/size(x, dim = 1) * size(x, dim = 2)/)), &
    2332           6 :             mask = reshape(maske(:, :), (/size(x, dim = 1) * size(x, dim = 2)/)))
    2333             :     !
    2334         246 :     NSE_sp_2d = 1.0_sp - sum((y - x) * (y - x), mask = maske) / sum((x - xmean) * (x - xmean), mask = maske)
    2335             :     !
    2336           2 :   END FUNCTION NSE_sp_2d
    2337             : 
    2338           4 :   FUNCTION NSE_dp_2d(x, y, mask)
    2339             : 
    2340           2 :     USE mo_moment, ONLY : average
    2341             : 
    2342             :     IMPLICIT NONE
    2343             : 
    2344             :     REAL(dp), DIMENSION(:, :), INTENT(IN) :: x, y
    2345             :     LOGICAL, DIMENSION(:, :), OPTIONAL, INTENT(IN) :: mask
    2346             :     REAL(dp) :: NSE_dp_2d
    2347             : 
    2348             :     INTEGER(i4) :: n
    2349             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    2350           2 :     REAL(dp) :: xmean
    2351           2 :     LOGICAL, DIMENSION(size(x, dim = 1), size(x, dim = 2)) :: maske
    2352             : 
    2353           2 :     if (present(mask)) then
    2354           3 :       shapemask = shape(mask)
    2355             :     else
    2356           3 :       shapemask = shape(x)
    2357             :     end if
    2358          12 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    2359           0 :             stop 'NSE_dp_2d: shapes of inputs(x,y) or mask are not matching'
    2360             :     !
    2361           2 :     if (present(mask)) then
    2362          62 :       maske = mask
    2363          61 :       n = count(maske)
    2364             :     else
    2365          61 :       maske = .true.
    2366           1 :       n = size(x, dim = 1) * size(x, dim = 2)
    2367             :     end if
    2368             :     !
    2369           2 :     if (n .LE. 1_i4) stop 'NSE_dp_2d: number of arguments must be at least 2'
    2370             :     ! mean of x
    2371             :     xmean = average(reshape(x(:, :), (/size(x, dim = 1) * size(x, dim = 2)/)), &
    2372           6 :             mask = reshape(maske(:, :), (/size(x, dim = 1) * size(x, dim = 2)/)))
    2373             :     !
    2374         246 :     NSE_dp_2d = 1.0_dp - sum((y - x) * (y - x), mask = maske) / sum((x - xmean) * (x - xmean), mask = maske)
    2375             :     !
    2376           2 :   END FUNCTION NSE_dp_2d
    2377             : 
    2378           4 :   FUNCTION NSE_sp_3d(x, y, mask)
    2379             : 
    2380           2 :     USE mo_moment, ONLY : average
    2381             : 
    2382             :     IMPLICIT NONE
    2383             : 
    2384             :     REAL(sp), DIMENSION(:, :, :), INTENT(IN) :: x, y
    2385             :     LOGICAL, DIMENSION(:, :, :), OPTIONAL, INTENT(IN) :: mask
    2386             :     REAL(sp) :: NSE_sp_3d
    2387             : 
    2388             :     INTEGER(i4) :: n
    2389             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    2390           2 :     REAL(sp) :: xmean
    2391             :     LOGICAL, DIMENSION(size(x, dim = 1), &
    2392           2 :             size(x, dim = 2), size(x, dim = 3)) :: maske
    2393             : 
    2394           2 :     if (present(mask)) then
    2395           4 :       shapemask = shape(mask)
    2396             :     else
    2397           4 :       shapemask = shape(x)
    2398             :     end if
    2399          16 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    2400           0 :             stop 'NSE_sp_3d: shapes of inputs(x,y) or mask are not matching'
    2401             :     !
    2402           2 :     if (present(mask)) then
    2403         490 :       maske = mask
    2404         489 :       n = count(maske)
    2405             :     else
    2406         489 :       maske = .true.
    2407           1 :       n = size(x, dim = 1) * size(x, dim = 2) * size(x, dim = 3)
    2408             :     end if
    2409             :     !
    2410           2 :     if (n .LE. 1_i4) stop 'NSE_sp_3d: number of arguments must be at least 2'
    2411             :     ! mean of x
    2412             :     xmean = average(reshape(x(:, :, :), (/size(x, dim = 1) * size(x, dim = 2) * size(x, dim = 3)/)), &
    2413           6 :             mask = reshape(maske(:, :, :), (/size(x, dim = 1) * size(x, dim = 2) * size(x, dim = 3)/)))
    2414             :     !
    2415        1958 :     NSE_sp_3d = 1.0_sp - sum((y - x) * (y - x), mask = maske) / sum((x - xmean) * (x - xmean), mask = maske)
    2416             :     !
    2417           2 :   END FUNCTION NSE_sp_3d
    2418             : 
    2419           4 :   FUNCTION NSE_dp_3d(x, y, mask)
    2420             : 
    2421           2 :     USE mo_moment, ONLY : average
    2422             : 
    2423             :     IMPLICIT NONE
    2424             : 
    2425             :     REAL(dp), DIMENSION(:, :, :), INTENT(IN) :: x, y
    2426             :     LOGICAL, DIMENSION(:, :, :), OPTIONAL, INTENT(IN) :: mask
    2427             :     REAL(dp) :: NSE_dp_3d
    2428             : 
    2429             :     INTEGER(i4) :: n
    2430             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    2431           2 :     REAL(dp) :: xmean
    2432             :     LOGICAL, DIMENSION(size(x, dim = 1), &
    2433           2 :             size(x, dim = 2), size(x, dim = 3)) :: maske
    2434             : 
    2435           2 :     if (present(mask)) then
    2436           4 :       shapemask = shape(mask)
    2437             :     else
    2438           4 :       shapemask = shape(x)
    2439             :     end if
    2440          16 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    2441           0 :             stop 'NSE_dp_3d: shapes of inputs(x,y) or mask are not matching'
    2442             :     !
    2443           2 :     if (present(mask)) then
    2444         490 :       maske = mask
    2445         489 :       n = count(maske)
    2446             :     else
    2447         489 :       maske = .true.
    2448           1 :       n = size(x, dim = 1) * size(x, dim = 2) * size(x, dim = 3)
    2449             :     end if
    2450             :     !
    2451           2 :     if (n .LE. 1_i4) stop 'NSE_dp_3d: number of arguments must be at least 2'
    2452             :     ! Average of x
    2453             :     xmean = average(reshape(x(:, :, :), (/size(x, dim = 1) * size(x, dim = 2) * size(x, dim = 3)/)), &
    2454           6 :             mask = reshape(maske(:, :, :), (/size(x, dim = 1) * size(x, dim = 2) * size(x, dim = 3)/)))
    2455             :     !
    2456        1958 :     NSE_dp_3d = 1.0_dp - sum((y - x) * (y - x), mask = maske) / sum((x - xmean) * (x - xmean), mask = maske)
    2457             :     !
    2458           2 :   END FUNCTION NSE_dp_3d
    2459             : 
    2460             : 
    2461             :   ! ------------------------------------------------------------------
    2462             : 
    2463          24 :   FUNCTION SAE_sp_1d(x, y, mask)
    2464             : 
    2465             :     IMPLICIT NONE
    2466             : 
    2467             :     REAL(sp), DIMENSION(:), INTENT(IN) :: x, y
    2468             :     LOGICAL, DIMENSION(:), OPTIONAL, INTENT(IN) :: mask
    2469             :     REAL(sp) :: SAE_sp_1d
    2470             : 
    2471             :     INTEGER(i4) :: n
    2472             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    2473          12 :     LOGICAL, DIMENSION(size(x)) :: maske
    2474             : 
    2475          12 :     if (present(mask)) then
    2476          22 :       shapemask = shape(mask)
    2477             :     else
    2478           2 :       shapemask = shape(x)
    2479             :     end if
    2480             :     !
    2481          48 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    2482           0 :             stop 'SAE_sp_1d: shapes of inputs(x,y) or mask are not matching'
    2483             :     !
    2484          12 :     if (present(mask)) then
    2485        2128 :       maske = mask
    2486        2117 :       n = count(maske)
    2487             :     else
    2488          55 :       maske = .true.
    2489           1 :       n = size(x)
    2490             :     end if
    2491          12 :     if (n .LE. 1_i4) stop 'SAE_sp_1d: number of arguments must be at least 2'
    2492             :     !
    2493        2184 :     SAE_sp_1d = sum(abs(y - x), mask = maske)
    2494             : 
    2495           2 :   END FUNCTION SAE_sp_1d
    2496             : 
    2497          24 :   FUNCTION SAE_dp_1d(x, y, mask)
    2498             : 
    2499             :     IMPLICIT NONE
    2500             : 
    2501             :     REAL(dp), DIMENSION(:), INTENT(IN) :: x, y
    2502             :     LOGICAL, DIMENSION(:), OPTIONAL, INTENT(IN) :: mask
    2503             :     REAL(dp) :: SAE_dp_1d
    2504             : 
    2505             :     INTEGER(i4) :: n
    2506             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    2507          12 :     LOGICAL, DIMENSION(size(x)) :: maske
    2508             : 
    2509          12 :     if (present(mask)) then
    2510          22 :       shapemask = shape(mask)
    2511             :     else
    2512           2 :       shapemask = shape(x)
    2513             :     end if
    2514             :     !
    2515          48 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    2516           0 :             stop 'SAE_dp_1d: shapes of inputs(x,y) or mask are not matching'
    2517             :     !
    2518          12 :     if (present(mask)) then
    2519        2128 :       maske = mask
    2520        2117 :       n = count(maske)
    2521             :     else
    2522          55 :       maske = .true.
    2523           1 :       n = size(x)
    2524             :     end if
    2525          12 :     if (n .LE. 1_i4) stop 'SAE_dp_1d: number of arguments must be at least 2'
    2526             :     !
    2527        2184 :     SAE_dp_1d = sum(abs(y - x), mask = maske)
    2528             : 
    2529          12 :   END FUNCTION SAE_dp_1d
    2530             : 
    2531           8 :   FUNCTION SAE_sp_2d(x, y, mask)
    2532             : 
    2533             :     IMPLICIT NONE
    2534             : 
    2535             :     REAL(sp), DIMENSION(:, :), INTENT(IN) :: x, y
    2536             :     LOGICAL, DIMENSION(:, :), OPTIONAL, INTENT(IN) :: mask
    2537             :     REAL(sp) :: SAE_sp_2d
    2538             : 
    2539             :     INTEGER(i4) :: n
    2540             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    2541           4 :     LOGICAL, DIMENSION(size(x, dim = 1), size(x, dim = 2)) :: maske
    2542             : 
    2543           4 :     if (present(mask)) then
    2544           9 :       shapemask = shape(mask)
    2545             :     else
    2546           3 :       shapemask = shape(x)
    2547             :     end if
    2548             :     !
    2549          24 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    2550           0 :             stop 'SAE_sp_2d: shapes of inputs(x,y) or mask are not matching'
    2551             :     !
    2552           4 :     if (present(mask)) then
    2553         186 :       maske = mask
    2554         183 :       n = count(maske)
    2555             :     else
    2556          61 :       maske = .true.
    2557           1 :       n = size(x, dim = 1) * size(x, dim = 2)
    2558             :     end if
    2559           4 :     if (n .LE. 1_i4) stop 'SAE_sp_2d: number of arguments must be at least 2'
    2560             :     !
    2561             :     SAE_sp_2d = SAE_sp_1d(reshape(x, (/size(x, dim = 1) * size(x, dim = 2)/)), &
    2562             :             reshape(y, (/size(y, dim = 1) * size(y, dim = 2)/)), &
    2563          16 :             mask = reshape(maske, (/size(maske, dim = 1) * size(maske, dim = 2)/)))
    2564             : 
    2565          12 :   END FUNCTION SAE_sp_2d
    2566             : 
    2567           8 :   FUNCTION SAE_dp_2d(x, y, mask)
    2568             : 
    2569             :     IMPLICIT NONE
    2570             : 
    2571             :     REAL(dp), DIMENSION(:, :), INTENT(IN) :: x, y
    2572             :     LOGICAL, DIMENSION(:, :), OPTIONAL, INTENT(IN) :: mask
    2573             :     REAL(dp) :: SAE_dp_2d
    2574             : 
    2575             :     INTEGER(i4) :: n
    2576             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    2577           4 :     LOGICAL, DIMENSION(size(x, dim = 1), size(x, dim = 2)) :: maske
    2578             : 
    2579           4 :     if (present(mask)) then
    2580           9 :       shapemask = shape(mask)
    2581             :     else
    2582           3 :       shapemask = shape(x)
    2583             :     end if
    2584             :     !
    2585          24 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    2586           0 :             stop 'SAE_dp_2d: shapes of inputs(x,y) or mask are not matching'
    2587             :     !
    2588           4 :     if (present(mask)) then
    2589         186 :       maske = mask
    2590         183 :       n = count(maske)
    2591             :     else
    2592          61 :       maske = .true.
    2593           1 :       n = size(x, dim = 1) * size(x, dim = 2)
    2594             :     end if
    2595           4 :     if (n .LE. 1_i4) stop 'SAE_dp_2d: number of arguments must be at least 2'
    2596             :     !
    2597             :     SAE_dp_2d = SAE_dp_1d(reshape(x, (/size(x, dim = 1) * size(x, dim = 2)/)), &
    2598             :             reshape(y, (/size(y, dim = 1) * size(y, dim = 2)/)), &
    2599          16 :             mask = reshape(maske, (/size(maske, dim = 1) * size(maske, dim = 2)/)))
    2600             : 
    2601           4 :   END FUNCTION SAE_dp_2d
    2602             : 
    2603           8 :   FUNCTION SAE_sp_3d(x, y, mask)
    2604             : 
    2605             :     IMPLICIT NONE
    2606             : 
    2607             :     REAL(sp), DIMENSION(:, :, :), INTENT(IN) :: x, y
    2608             :     LOGICAL, DIMENSION(:, :, :), OPTIONAL, INTENT(IN) :: mask
    2609             :     REAL(sp) :: SAE_sp_3d
    2610             : 
    2611             :     INTEGER(i4) :: n
    2612             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    2613             :     LOGICAL, DIMENSION(size(x, dim = 1), size(x, dim = 2), &
    2614           4 :             size(x, dim = 3)) :: maske
    2615             : 
    2616           4 :     if (present(mask)) then
    2617          12 :       shapemask = shape(mask)
    2618             :     else
    2619           4 :       shapemask = shape(x)
    2620             :     end if
    2621             :     !
    2622          32 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    2623           0 :             stop 'SAE_sp_3d: shapes of inputs(x,y) or mask are not matching'
    2624             :     !
    2625           4 :     if (present(mask)) then
    2626        1470 :       maske = mask
    2627        1467 :       n = count(maske)
    2628             :     else
    2629         489 :       maske = .true.
    2630           1 :       n = size(x, dim = 1) * size(x, dim = 2) * size(x, dim = 3)
    2631             :     end if
    2632           4 :     if (n .LE. 1_i4) stop 'SAE_sp_3d: number of arguments must be at least 2'
    2633             :     !
    2634             :     SAE_sp_3d = SAE_sp_1d(reshape(x, (/size(x, dim = 1) * size(x, dim = 2) * size(x, dim = 3)/)), &
    2635             :             reshape(y, (/size(y, dim = 1) * size(y, dim = 2) * size(x, dim = 3)/)), &
    2636             :             mask = reshape(maske, (/size(maske, dim = 1) * size(maske, dim = 2) &
    2637          16 :                     * size(maske, dim = 3)/)))
    2638             : 
    2639           4 :   END FUNCTION SAE_sp_3d
    2640             : 
    2641           8 :   FUNCTION SAE_dp_3d(x, y, mask)
    2642             : 
    2643             :     IMPLICIT NONE
    2644             : 
    2645             :     REAL(dp), DIMENSION(:, :, :), INTENT(IN) :: x, y
    2646             :     LOGICAL, DIMENSION(:, :, :), OPTIONAL, INTENT(IN) :: mask
    2647             :     REAL(dp) :: SAE_dp_3d
    2648             : 
    2649             :     INTEGER(i4) :: n
    2650             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    2651             :     LOGICAL, DIMENSION(size(x, dim = 1), size(x, dim = 2), &
    2652           4 :             size(x, dim = 3)) :: maske
    2653             : 
    2654           4 :     if (present(mask)) then
    2655          12 :       shapemask = shape(mask)
    2656             :     else
    2657           4 :       shapemask = shape(x)
    2658             :     end if
    2659             :     !
    2660          32 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    2661           0 :             stop 'SAE_dp_3d: shapes of inputs(x,y) or mask are not matching'
    2662             :     !
    2663           4 :     if (present(mask)) then
    2664        1470 :       maske = mask
    2665        1467 :       n = count(maske)
    2666             :     else
    2667         489 :       maske = .true.
    2668           1 :       n = size(x, dim = 1) * size(x, dim = 2) * size(x, dim = 3)
    2669             :     end if
    2670           4 :     if (n .LE. 1_i4) stop 'SAE_dp_3d: number of arguments must be at least 2'
    2671             :     !
    2672             :     SAE_dp_3d = SAE_dp_1d(reshape(x, (/size(x, dim = 1) * size(x, dim = 2) * size(x, dim = 3)/)), &
    2673             :             reshape(y, (/size(y, dim = 1) * size(y, dim = 2) * size(x, dim = 3)/)), &
    2674             :             mask = reshape(maske, (/size(maske, dim = 1) * size(maske, dim = 2) &
    2675          16 :                     * size(maske, dim = 3)/)))
    2676             : 
    2677           4 :   END FUNCTION SAE_dp_3d
    2678             : 
    2679             :   ! ------------------------------------------------------------------
    2680             : 
    2681          36 :   FUNCTION SSE_sp_1d(x, y, mask)
    2682             : 
    2683             :     IMPLICIT NONE
    2684             : 
    2685             :     REAL(sp), DIMENSION(:), INTENT(IN) :: x, y
    2686             :     LOGICAL, DIMENSION(:), OPTIONAL, INTENT(IN) :: mask
    2687             :     REAL(sp) :: SSE_sp_1d
    2688             : 
    2689             :     INTEGER(i4) :: n
    2690             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    2691          18 :     LOGICAL, DIMENSION(size(x)) :: maske
    2692             : 
    2693          18 :     if (present(mask)) then
    2694          34 :       shapemask = shape(mask)
    2695             :     else
    2696           2 :       shapemask = shape(x)
    2697             :     end if
    2698             :     !
    2699          72 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    2700           0 :             stop 'SSE_sp_1d: shapes of inputs(x,y) or mask are not matching'
    2701             :     !
    2702          18 :     if (present(mask)) then
    2703        3220 :       maske = mask
    2704        3203 :       n = count(maske)
    2705             :     else
    2706          55 :       maske = .true.
    2707           1 :       n = size(x)
    2708             :     end if
    2709          18 :     if (n .LE. 1_i4) stop 'SSE_sp_1d: number of arguments must be at least 2'
    2710             :     !
    2711        3276 :     SSE_sp_1d = sum((y - x)**2_i4, mask = maske)
    2712             : 
    2713           4 :   END FUNCTION SSE_sp_1d
    2714             : 
    2715          36 :   FUNCTION SSE_dp_1d(x, y, mask)
    2716             : 
    2717             :     IMPLICIT NONE
    2718             : 
    2719             :     REAL(dp), DIMENSION(:), INTENT(IN) :: x, y
    2720             :     LOGICAL, DIMENSION(:), OPTIONAL, INTENT(IN) :: mask
    2721             :     REAL(dp) :: SSE_dp_1d
    2722             : 
    2723             :     INTEGER(i4) :: n
    2724             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    2725          18 :     LOGICAL, DIMENSION(size(x)) :: maske
    2726             : 
    2727          18 :     if (present(mask)) then
    2728          34 :       shapemask = shape(mask)
    2729             :     else
    2730           2 :       shapemask = shape(x)
    2731             :     end if
    2732          72 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    2733           0 :             stop 'SSE_dp_1d: shapes of inputs(x,y) or mask are not matching'
    2734             :     !
    2735          18 :     if (present(mask)) then
    2736        3220 :       maske = mask
    2737        3203 :       n = count(maske)
    2738             :     else
    2739          55 :       maske = .true.
    2740           1 :       n = size(x)
    2741             :     end if
    2742          18 :     if (n .LE. 1_i4) stop 'SSE_dp_1d: number of arguments must be at least 2'
    2743             :     !
    2744        3276 :     SSE_dp_1d = sum((y - x)**2_i4, mask = maske)
    2745             : 
    2746          18 :   END FUNCTION SSE_dp_1d
    2747             : 
    2748          12 :   FUNCTION SSE_sp_2d(x, y, mask)
    2749             : 
    2750             :     IMPLICIT NONE
    2751             : 
    2752             :     REAL(sp), DIMENSION(:, :), INTENT(IN) :: x, y
    2753             :     LOGICAL, DIMENSION(:, :), OPTIONAL, INTENT(IN) :: mask
    2754             :     REAL(sp) :: SSE_sp_2d
    2755             : 
    2756             :     INTEGER(i4) :: n
    2757             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    2758           6 :     LOGICAL, DIMENSION(size(x, dim = 1), size(x, dim = 2)) :: maske
    2759             : 
    2760           6 :     if (present(mask)) then
    2761          15 :       shapemask = shape(mask)
    2762             :     else
    2763           3 :       shapemask = shape(x)
    2764             :     end if
    2765          36 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    2766           0 :             stop 'SSE_sp_2d: shapes of inputs(x,y) or mask are not matching'
    2767             :     !
    2768           6 :     if (present(mask)) then
    2769         310 :       maske = mask
    2770         305 :       n = count(maske)
    2771             :     else
    2772          61 :       maske = .true.
    2773           3 :       n = size(x)
    2774             :     end if
    2775           6 :     if (n .LE. 1_i4) stop 'SSE_sp_2d: number of arguments must be at least 2'
    2776             :     !
    2777             :     SSE_sp_2d = SSE_sp_1d(reshape(x, (/size(x, dim = 1) * size(x, dim = 2)/)), &
    2778             :             reshape(y, (/size(y, dim = 1) * size(y, dim = 2)/)), &
    2779          24 :             mask = reshape(maske, (/size(maske, dim = 1) * size(maske, dim = 2)/)))
    2780             : 
    2781          18 :   END FUNCTION SSE_sp_2d
    2782             : 
    2783          12 :   FUNCTION SSE_dp_2d(x, y, mask)
    2784             : 
    2785             :     IMPLICIT NONE
    2786             : 
    2787             :     REAL(dp), DIMENSION(:, :), INTENT(IN) :: x, y
    2788             :     LOGICAL, DIMENSION(:, :), OPTIONAL, INTENT(IN) :: mask
    2789             :     REAL(dp) :: SSE_dp_2d
    2790             : 
    2791             :     INTEGER(i4) :: n
    2792             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    2793           6 :     LOGICAL, DIMENSION(size(x, dim = 1), size(x, dim = 2)) :: maske
    2794             : 
    2795           6 :     if (present(mask)) then
    2796          15 :       shapemask = shape(mask)
    2797             :     else
    2798           3 :       shapemask = shape(x)
    2799             :     end if
    2800          36 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    2801           0 :             stop 'SSE_dp_2d: shapes of inputs(x,y) or mask are not matching'
    2802             :     !
    2803           6 :     if (present(mask)) then
    2804         310 :       maske = mask
    2805         305 :       n = count(maske)
    2806             :     else
    2807          61 :       maske = .true.
    2808           3 :       n = size(x)
    2809             :     end if
    2810           6 :     if (n .LE. 1_i4) stop 'SSE_dp_2d: number of arguments must be at least 2'
    2811             :     !
    2812             :     SSE_dp_2d = SSE_dp_1d(reshape(x, (/size(x, dim = 1) * size(x, dim = 2)/)), &
    2813             :             reshape(y, (/size(y, dim = 1) * size(y, dim = 2)/)), &
    2814          24 :             mask = reshape(maske, (/size(maske, dim = 1) * size(maske, dim = 2)/)))
    2815             : 
    2816           6 :   END FUNCTION SSE_dp_2d
    2817             : 
    2818          12 :   FUNCTION SSE_sp_3d(x, y, mask)
    2819             : 
    2820             :     IMPLICIT NONE
    2821             : 
    2822             :     REAL(sp), DIMENSION(:, :, :), INTENT(IN) :: x, y
    2823             :     LOGICAL, DIMENSION(:, :, :), OPTIONAL, INTENT(IN) :: mask
    2824             :     REAL(sp) :: SSE_sp_3d
    2825             : 
    2826             :     INTEGER(i4) :: n
    2827             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    2828             :     LOGICAL, DIMENSION(size(x, dim = 1), size(x, dim = 2), &
    2829           6 :             size(x, dim = 3)) :: maske
    2830             : 
    2831           6 :     if (present(mask)) then
    2832          20 :       shapemask = shape(mask)
    2833             :     else
    2834           4 :       shapemask = shape(x)
    2835             :     end if
    2836          48 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    2837           0 :             stop 'SSE_sp_3d: shapes of inputs(x,y) or mask are not matching'
    2838             :     !
    2839           6 :     if (present(mask)) then
    2840        2450 :       maske = mask
    2841        2445 :       n = count(maske)
    2842             :     else
    2843         489 :       maske = .true.
    2844           4 :       n = size(x)
    2845             :     end if
    2846           6 :     if (n .LE. 1_i4) stop 'SSE_sp_3d: number of arguments must be at least 2'
    2847             :     !
    2848             :     SSE_sp_3d = SSE_sp_1d(reshape(x, (/size(x, dim = 1) * size(x, dim = 2) * size(x, dim = 3)/)), &
    2849             :             reshape(y, (/size(y, dim = 1) * size(y, dim = 2) * size(x, dim = 3)/)), &
    2850             :             mask = reshape(maske, (/size(maske, dim = 1) * size(maske, dim = 2)   &
    2851          24 :                     * size(maske, dim = 3)/)))
    2852             : 
    2853           6 :   END FUNCTION SSE_sp_3d
    2854             : 
    2855          12 :   FUNCTION SSE_dp_3d(x, y, mask)
    2856             : 
    2857             :     IMPLICIT NONE
    2858             : 
    2859             :     REAL(dp), DIMENSION(:, :, :), INTENT(IN) :: x, y
    2860             :     LOGICAL, DIMENSION(:, :, :), OPTIONAL, INTENT(IN) :: mask
    2861             :     REAL(dp) :: SSE_dp_3d
    2862             : 
    2863             :     INTEGER(i4) :: n
    2864             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    2865             :     LOGICAL, DIMENSION(size(x, dim = 1), size(x, dim = 2), &
    2866           6 :             size(x, dim = 3)) :: maske
    2867             : 
    2868           6 :     if (present(mask)) then
    2869          20 :       shapemask = shape(mask)
    2870             :     else
    2871           4 :       shapemask = shape(x)
    2872             :     end if
    2873          48 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    2874           0 :             stop 'SSE_dp_3d: shapes of inputs(x,y) or mask are not matching'
    2875             :     !
    2876           6 :     if (present(mask)) then
    2877        2450 :       maske = mask
    2878        2445 :       n = count(maske)
    2879             :     else
    2880         489 :       maske = .true.
    2881           4 :       n = size(x)
    2882             :     end if
    2883           6 :     if (n .LE. 1_i4) stop 'SSE_dp_3d: number of arguments must be at least 2'
    2884             :     !
    2885             :     SSE_dp_3d = SSE_dp_1d(reshape(x, (/size(x, dim = 1) * size(x, dim = 2) * size(x, dim = 3)/)), &
    2886             :             reshape(y, (/size(y, dim = 1) * size(y, dim = 2) * size(x, dim = 3)/)), &
    2887             :             mask = reshape(maske, (/size(maske, dim = 1) * size(maske, dim = 2)   &
    2888          24 :                     * size(maske, dim = 3)/)))
    2889             : 
    2890           6 :   END FUNCTION SSE_dp_3d
    2891             : 
    2892             :   ! ------------------------------------------------------------------
    2893             : 
    2894           4 :   FUNCTION RMSE_sp_1d(x, y, mask)
    2895             : 
    2896             :     IMPLICIT NONE
    2897             : 
    2898             :     REAL(sp), DIMENSION(:), INTENT(IN) :: x, y
    2899             :     LOGICAL, DIMENSION(:), OPTIONAL, INTENT(IN) :: mask
    2900             :     REAL(sp) :: RMSE_sp_1d
    2901             : 
    2902             :     INTEGER(i4) :: n
    2903             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    2904           2 :     LOGICAL, DIMENSION(size(x, dim = 1)) :: maske
    2905             : 
    2906           2 :     if (present(mask)) then
    2907           2 :       shapemask = shape(mask)
    2908             :     else
    2909           2 :       shapemask = shape(x)
    2910             :     end if
    2911           8 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    2912           0 :             stop 'RMSE_sp_1d: shapes of inputs(x,y) or mask are not matching'
    2913             :     !
    2914           2 :     if (present(mask)) then
    2915          55 :       maske = mask
    2916          55 :       n = count(maske)
    2917             :     else
    2918          55 :       maske = .true.
    2919           1 :       n = size(x, dim = 1)
    2920             :     end if
    2921           2 :     if (n .LE. 1_i4) stop 'RMSE_sp_1d: number of arguments must be at least 2'
    2922             :     !
    2923           2 :     RMSE_sp_1d = sqrt(MSE_sp_1d(x, y, mask = maske))
    2924             : 
    2925           6 :   END FUNCTION RMSE_sp_1d
    2926             : 
    2927           4 :   FUNCTION RMSE_dp_1d(x, y, mask)
    2928             : 
    2929             :     IMPLICIT NONE
    2930             : 
    2931             :     REAL(dp), DIMENSION(:), INTENT(IN) :: x, y
    2932             :     LOGICAL, DIMENSION(:), OPTIONAL, INTENT(IN) :: mask
    2933             :     REAL(dp) :: RMSE_dp_1d
    2934             : 
    2935             :     INTEGER(i4) :: n
    2936             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    2937           2 :     LOGICAL, DIMENSION(size(x, dim = 1)) :: maske
    2938             : 
    2939           2 :     if (present(mask)) then
    2940           2 :       shapemask = shape(mask)
    2941             :     else
    2942           2 :       shapemask = shape(x)
    2943             :     end if
    2944           8 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    2945           0 :             stop 'RMSE_dp_1d: shapes of inputs(x,y) or mask are not matching'
    2946             :     !
    2947           2 :     if (present(mask)) then
    2948          55 :       maske = mask
    2949          55 :       n = count(maske)
    2950             :     else
    2951          55 :       maske = .true.
    2952           1 :       n = size(x, dim = 1)
    2953             :     end if
    2954           2 :     if (n .LE. 1_i4) stop 'RMSE_dp_1d: number of arguments must be at least 2'
    2955             :     !
    2956           2 :     RMSE_dp_1d = sqrt(MSE_dp_1d(x, y, mask = maske))
    2957             : 
    2958           2 :   END FUNCTION RMSE_dp_1d
    2959             : 
    2960           4 :   FUNCTION RMSE_sp_2d(x, y, mask)
    2961             : 
    2962             :     IMPLICIT NONE
    2963             : 
    2964             :     REAL(sp), DIMENSION(:, :), INTENT(IN) :: x, y
    2965             :     LOGICAL, DIMENSION(:, :), OPTIONAL, INTENT(IN) :: mask
    2966             :     REAL(sp) :: RMSE_sp_2d
    2967             : 
    2968             :     INTEGER(i4) :: n
    2969             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    2970           2 :     LOGICAL, DIMENSION(size(x, dim = 1), size(x, dim = 2)) :: maske
    2971             : 
    2972           2 :     if (present(mask)) then
    2973           3 :       shapemask = shape(mask)
    2974             :     else
    2975           3 :       shapemask = shape(x)
    2976             :     end if
    2977          12 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    2978           0 :             stop 'RMSE_sp_2d: shapes of inputs(x,y) or mask are not matching'
    2979             :     !
    2980           2 :     if (present(mask)) then
    2981          61 :       maske = mask
    2982          61 :       n = count(maske)
    2983             :     else
    2984          61 :       maske = .true.
    2985           1 :       n = size(x, dim = 1) * size(x, dim = 2)
    2986             :     end if
    2987           2 :     if (n .LE. 1_i4) stop 'RMSE_sp_2d: number of arguments must be at least 2'
    2988             :     !
    2989           2 :     RMSE_sp_2d = sqrt(MSE_sp_2d(x, y, mask = maske))
    2990             : 
    2991           2 :   END FUNCTION RMSE_sp_2d
    2992             : 
    2993           4 :   FUNCTION RMSE_dp_2d(x, y, mask)
    2994             : 
    2995             :     IMPLICIT NONE
    2996             : 
    2997             :     REAL(dp), DIMENSION(:, :), INTENT(IN) :: x, y
    2998             :     LOGICAL, DIMENSION(:, :), OPTIONAL, INTENT(IN) :: mask
    2999             :     REAL(dp) :: RMSE_dp_2d
    3000             : 
    3001             :     INTEGER(i4) :: n
    3002             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    3003           2 :     LOGICAL, DIMENSION(size(x, dim = 1), size(x, dim = 2)) :: maske
    3004             : 
    3005           2 :     if (present(mask)) then
    3006           3 :       shapemask = shape(mask)
    3007             :     else
    3008           3 :       shapemask = shape(x)
    3009             :     end if
    3010          12 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    3011           0 :             stop 'RMSE_dp_2d: shapes of inputs(x,y) or mask are not matching'
    3012             :     !
    3013           2 :     if (present(mask)) then
    3014          61 :       maske = mask
    3015          61 :       n = count(maske)
    3016             :     else
    3017          61 :       maske = .true.
    3018           1 :       n = size(x, dim = 1) * size(x, dim = 2)
    3019             :     end if
    3020           2 :     if (n .LE. 1_i4) stop 'RMSE_dp_2d: number of arguments must be at least 2'
    3021             :     !
    3022           2 :     RMSE_dp_2d = sqrt(MSE_dp_2d(x, y, mask = maske))
    3023             : 
    3024           2 :   END FUNCTION RMSE_dp_2d
    3025             : 
    3026           4 :   FUNCTION RMSE_sp_3d(x, y, mask)
    3027             : 
    3028             :     IMPLICIT NONE
    3029             : 
    3030             :     REAL(sp), DIMENSION(:, :, :), INTENT(IN) :: x, y
    3031             :     LOGICAL, DIMENSION(:, :, :), OPTIONAL, INTENT(IN) :: mask
    3032             :     REAL(sp) :: RMSE_sp_3d
    3033             : 
    3034             :     INTEGER(i4) :: n
    3035             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    3036             :     LOGICAL, DIMENSION(size(x, dim = 1), size(x, dim = 2), &
    3037           2 :             size(x, dim = 3)) :: maske
    3038             : 
    3039           2 :     if (present(mask)) then
    3040           4 :       shapemask = shape(mask)
    3041             :     else
    3042           4 :       shapemask = shape(x)
    3043             :     end if
    3044          16 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    3045           0 :             stop 'RMSE_sp_3d: shapes of inputs(x,y) or mask are not matching'
    3046             :     !
    3047           2 :     if (present(mask)) then
    3048         489 :       maske = mask
    3049         489 :       n = count(maske)
    3050             :     else
    3051         489 :       maske = .true.
    3052           1 :       n = size(x, dim = 1) * size(x, dim = 2) * size(x, dim = 3)
    3053             :     end if
    3054           2 :     if (n .LE. 1_i4) stop 'RMSE_sp_3d: number of arguments must be at least 2'
    3055             :     !
    3056           2 :     RMSE_sp_3d = sqrt(MSE_sp_3d(x, y, mask = maske))
    3057             : 
    3058           2 :   END FUNCTION RMSE_sp_3d
    3059             : 
    3060           4 :   FUNCTION RMSE_dp_3d(x, y, mask)
    3061             : 
    3062             :     IMPLICIT NONE
    3063             : 
    3064             :     REAL(dp), DIMENSION(:, :, :), INTENT(IN) :: x, y
    3065             :     LOGICAL, DIMENSION(:, :, :), OPTIONAL, INTENT(IN) :: mask
    3066             :     REAL(dp) :: RMSE_dp_3d
    3067             : 
    3068             :     INTEGER(i4) :: n
    3069             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    3070             :     LOGICAL, DIMENSION(size(x, dim = 1), size(x, dim = 2), &
    3071           2 :             size(x, dim = 3)) :: maske
    3072             : 
    3073           2 :     if (present(mask)) then
    3074           4 :       shapemask = shape(mask)
    3075             :     else
    3076           4 :       shapemask = shape(x)
    3077             :     end if
    3078          16 :     if ((any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask))) &
    3079           0 :             stop 'RMSE_dp_3d: shapes of inputs(x,y) or mask are not matching'
    3080             :     !
    3081           2 :     if (present(mask)) then
    3082         489 :       maske = mask
    3083         489 :       n = count(maske)
    3084             :     else
    3085         489 :       maske = .true.
    3086           1 :       n = size(x, dim = 1) * size(x, dim = 2) * size(x, dim = 3)
    3087             :     end if
    3088           2 :     if (n .LE. 1_i4) stop 'RMSE_dp_3d: number of arguments must be at least 2'
    3089             :     !
    3090           2 :     RMSE_dp_3d = sqrt(MSE_dp_3d(x, y, mask = maske))
    3091             : 
    3092           2 :   END FUNCTION RMSE_dp_3d
    3093             : 
    3094             :   ! ------------------------------------------------------------------
    3095             : 
    3096           4 :   FUNCTION wNSE_sp_1d(x, y, mask)
    3097             : 
    3098           2 :     USE mo_moment, ONLY: average
    3099             : 
    3100             :     IMPLICIT NONE
    3101             : 
    3102             :     REAL(sp), DIMENSION(:),           INTENT(IN)      :: x, y
    3103             :     LOGICAL,  DIMENSION(:), OPTIONAL, INTENT(IN)      :: mask
    3104             :     REAL(sp)                                          :: wNSE_sp_1d
    3105             : 
    3106             :     INTEGER(i4)                            :: n
    3107             :     INTEGER(i4), DIMENSION(size(shape(x))) :: shapemask
    3108         110 :     REAL(sp)                               :: xmean
    3109         332 :     REAL(sp), DIMENSION(size(x))           :: v1, v2, ww
    3110           2 :     LOGICAL,  DIMENSION(size(x))           :: maske
    3111             : 
    3112           2 :     if (present(mask)) then
    3113           2 :        shapemask = shape(mask)
    3114             :     else
    3115           2 :        shapemask =  shape(x)
    3116             :     end if
    3117           8 :     if ( (any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask)) ) &
    3118           0 :          stop 'wNSE_sp_1d: shapes of inputs(x,y) or mask are not matching'
    3119             :     !
    3120           2 :     if (present(mask)) then
    3121          56 :        maske = mask
    3122          55 :        n = count(maske)
    3123             :     else
    3124          55 :        maske = .true.
    3125           1 :        n = size(x)
    3126             :     end if
    3127           2 :     if (n .LE. 1_i4) stop 'wNSE_sp_1d: number of arguments must be at least 2'
    3128             :     ! mean of x
    3129           2 :     xmean = average(x, mask=maske)
    3130             :     !
    3131         112 :     v1 = merge(y - x    , 0.0_sp, maske)
    3132         112 :     v2 = merge(x - xmean, 0.0_sp, maske)
    3133         110 :     ww = merge(x        , 0.0_sp, maske)
    3134             :     !
    3135         220 :     wNSE_sp_1d = 1.0_sp - dot_product(ww * v1,v1) / dot_product(ww * v2,v2)
    3136             : 
    3137           2 :   END FUNCTION wNSE_sp_1d
    3138             : 
    3139           4 :   FUNCTION wNSE_dp_1d(x, y, mask)
    3140             : 
    3141           2 :     USE mo_moment, ONLY: average
    3142             : 
    3143             :     IMPLICIT NONE
    3144             : 
    3145             :     REAL(dp), DIMENSION(:),           INTENT(IN)      :: x, y
    3146             :     LOGICAL,  DIMENSION(:), OPTIONAL, INTENT(IN)      :: mask
    3147             :     REAL(dp)                                          :: wNSE_dp_1d
    3148             : 
    3149             :     INTEGER(i4)                                       :: n
    3150             :     INTEGER(i4), DIMENSION(size(shape(x)) )           :: shapemask
    3151         110 :     REAL(dp)                                          :: xmean
    3152         332 :     REAL(dp), DIMENSION(size(x))                      :: v1, v2, ww
    3153           2 :     LOGICAL,  DIMENSION(size(x))                      :: maske
    3154             : 
    3155           2 :     if (present(mask)) then
    3156           2 :        shapemask = shape(mask)
    3157             :     else
    3158           2 :        shapemask =  shape(x)
    3159             :     end if
    3160           8 :     if ( (any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask)) ) &
    3161           0 :          stop 'wNSE_dp_1d: shapes of inputs(x,y) or mask are not matching'
    3162             :     !
    3163           2 :     if (present(mask)) then
    3164          56 :        maske = mask
    3165          55 :        n = count(maske)
    3166             :     else
    3167          55 :        maske = .true.
    3168           1 :        n = size(x)
    3169             :     end if
    3170           2 :     if (n .LE. 1_i4) stop 'wNSE_dp_1d: number of arguments must be at least 2'
    3171             :     ! mean of x
    3172           2 :     xmean = average(x, mask=maske)
    3173             :     !
    3174         112 :     v1 = merge(y - x    , 0.0_dp, maske)
    3175         112 :     v2 = merge(x - xmean, 0.0_dp, maske)
    3176         110 :     ww = merge(x        , 0.0_dp, maske)
    3177             :     !
    3178         220 :     wNSE_dp_1d = 1.0_dp - dot_product(ww * v1,v1) / dot_product(ww * v2,v2)
    3179             : 
    3180           2 :   END FUNCTION wNSE_dp_1d
    3181             : 
    3182           4 :   FUNCTION wNSE_sp_2d(x, y, mask)
    3183             : 
    3184           2 :     USE mo_moment, ONLY: average
    3185             : 
    3186             :     IMPLICIT NONE
    3187             : 
    3188             :     REAL(sp), DIMENSION(:,:),           INTENT(IN)    :: x, y
    3189             :     LOGICAL,  DIMENSION(:,:), OPTIONAL, INTENT(IN)    :: mask
    3190             :     REAL(sp)                                          :: wNSE_sp_2d
    3191             : 
    3192             :     INTEGER(i4)                                       :: n
    3193             :     INTEGER(i4), DIMENSION(size(shape(x)) )           :: shapemask
    3194           2 :     REAL(sp)                                          :: xmean
    3195           2 :     LOGICAL, DIMENSION(size(x, dim=1), size(x, dim=2)):: maske
    3196             : 
    3197           2 :     if (present(mask)) then
    3198           3 :        shapemask = shape(mask)
    3199             :     else
    3200           3 :        shapemask =  shape(x)
    3201             :     end if
    3202          12 :     if ( (any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask)) ) &
    3203           0 :          stop 'wNSE_sp_2d: shapes of inputs(x,y) or mask are not matching'
    3204             :     !
    3205           2 :     if (present(mask)) then
    3206          62 :        maske = mask
    3207          61 :        n     = count(maske)
    3208             :     else
    3209          61 :        maske = .true.
    3210           1 :        n     = size(x, dim=1) * size(x, dim=2)
    3211             :     end if
    3212             :     !
    3213           2 :     if (n .LE. 1_i4) stop 'wNSE_sp_2d: number of arguments must be at least 2'
    3214             :     ! mean of x
    3215             :     xmean = average(reshape(x(:,:), (/size(x, dim=1)*size(x, dim=2)/)), &
    3216           6 :          mask=reshape(maske(:,:), (/size(x, dim=1)*size(x, dim=2)/)))
    3217             :     !
    3218         246 :     wNSE_sp_2d = 1.0_sp - sum(x * (y-x)*(y-x), mask=maske) / sum(x * (x-xmean)*(x-xmean), mask=maske)
    3219             :     !
    3220           2 :   END FUNCTION wNSE_sp_2d
    3221             : 
    3222           4 :   FUNCTION wNSE_dp_2d(x, y, mask)
    3223             : 
    3224           2 :     USE mo_moment, ONLY: average
    3225             : 
    3226             :     IMPLICIT NONE
    3227             : 
    3228             :     REAL(dp), DIMENSION(:,:),           INTENT(IN)    :: x, y
    3229             :     LOGICAL,  DIMENSION(:,:), OPTIONAL, INTENT(IN)    :: mask
    3230             :     REAL(dp)                                          :: wNSE_dp_2d
    3231             : 
    3232             :     INTEGER(i4)                                       :: n
    3233             :     INTEGER(i4), DIMENSION(size(shape(x)) )           :: shapemask
    3234           2 :     REAL(dp)                                          :: xmean
    3235           2 :     LOGICAL, DIMENSION(size(x, dim=1), size(x, dim=2)):: maske
    3236             : 
    3237           2 :     if (present(mask)) then
    3238           3 :        shapemask = shape(mask)
    3239             :     else
    3240           3 :        shapemask =  shape(x)
    3241             :     end if
    3242          12 :     if ( (any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask)) ) &
    3243           0 :          stop 'wNSE_dp_2d: shapes of inputs(x,y) or mask are not matching'
    3244             :     !
    3245           2 :     if (present(mask)) then
    3246          62 :        maske = mask
    3247          61 :        n     = count(maske)
    3248             :     else
    3249          61 :        maske = .true.
    3250           1 :        n     = size(x, dim=1) * size(x, dim=2)
    3251             :     end if
    3252             :     !
    3253           2 :     if (n .LE. 1_i4) stop 'wNSE_dp_2d: number of arguments must be at least 2'
    3254             :     ! mean of x
    3255             :     xmean = average(reshape(x(:,:), (/size(x, dim=1)*size(x, dim=2)/)), &
    3256           6 :          mask=reshape(maske(:,:), (/size(x, dim=1)*size(x, dim=2)/)))
    3257             :     !
    3258         246 :     wNSE_dp_2d = 1.0_dp - sum(x * (y-x)*(y-x), mask=maske) / sum(x * (x-xmean)*(x-xmean), mask=maske)
    3259             :     !
    3260           2 :   END FUNCTION wNSE_dp_2d
    3261             : 
    3262           4 :   FUNCTION wNSE_sp_3d(x, y, mask)
    3263             : 
    3264           2 :     USE mo_moment, ONLY: average
    3265             : 
    3266             :     IMPLICIT NONE
    3267             : 
    3268             :     REAL(sp), DIMENSION(:,:,:),           INTENT(IN)  :: x, y
    3269             :     LOGICAL,  DIMENSION(:,:,:), OPTIONAL, INTENT(IN)  :: mask
    3270             :     REAL(sp)                                          :: wNSE_sp_3d
    3271             : 
    3272             :     INTEGER(i4)                                       :: n
    3273             :     INTEGER(i4), DIMENSION(size(shape(x)) )           :: shapemask
    3274           2 :     REAL(sp)                                          :: xmean
    3275             :     LOGICAL,  DIMENSION(size(x, dim=1), &
    3276           2 :          size(x, dim=2), size(x, dim=3))    :: maske
    3277             : 
    3278           2 :     if (present(mask)) then
    3279           4 :        shapemask = shape(mask)
    3280             :     else
    3281           4 :        shapemask =  shape(x)
    3282             :     end if
    3283          16 :     if ( (any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask)) ) &
    3284           0 :          stop 'wNSE_sp_3d: shapes of inputs(x,y) or mask are not matching'
    3285             :     !
    3286           2 :     if (present(mask)) then
    3287         490 :        maske = mask
    3288         489 :        n     = count(maske)
    3289             :     else
    3290         489 :        maske = .true.
    3291           1 :        n     = size(x, dim=1) * size(x, dim=2) * size(x, dim=3)
    3292             :     end if
    3293             :     !
    3294           2 :     if (n .LE. 1_i4) stop 'wNSE_sp_3d: number of arguments must be at least 2'
    3295             :     ! mean of x
    3296             :     xmean = average(reshape(x(:,:,:), (/size(x, dim=1)*size(x, dim=2)*size(x, dim=3)/)), &
    3297           6 :          mask=reshape(maske(:,:,:), (/size(x, dim=1)*size(x, dim=2)*size(x, dim=3)/)))
    3298             :     !
    3299        1958 :     wNSE_sp_3d = 1.0_sp - sum(x * (y-x)*(y-x), mask=maske) / sum(x * (x-xmean)*(x-xmean), mask=maske)
    3300             :     !
    3301           2 :   END FUNCTION wNSE_sp_3d
    3302             : 
    3303           2 :   FUNCTION wNSE_dp_3d(x, y, mask)
    3304             : 
    3305           2 :     USE mo_moment, ONLY: average
    3306             : 
    3307             :     IMPLICIT NONE
    3308             : 
    3309             :     REAL(dp), DIMENSION(:,:,:),           INTENT(IN)  :: x, y
    3310             :     LOGICAL,  DIMENSION(:,:,:), OPTIONAL, INTENT(IN)  :: mask
    3311             :     REAL(dp)                                          :: wNSE_dp_3d
    3312             : 
    3313             :     INTEGER(i4)                                       :: n
    3314             :     INTEGER(i4), DIMENSION(size(shape(x)) )           :: shapemask
    3315           2 :     REAL(dp)                                          :: xmean
    3316             :     LOGICAL,  DIMENSION(size(x, dim=1), &
    3317           2 :          size(x, dim=2), size(x, dim=3))    :: maske
    3318             : 
    3319           2 :     if (present(mask)) then
    3320           4 :        shapemask = shape(mask)
    3321             :     else
    3322           4 :        shapemask =  shape(x)
    3323             :     end if
    3324          16 :     if ( (any(shape(x) .NE. shape(y))) .OR. (any(shape(x) .NE. shapemask)) ) &
    3325           0 :          stop 'wNSE_dp_3d: shapes of inputs(x,y) or mask are not matching'
    3326             :     !
    3327           2 :     if (present(mask)) then
    3328         490 :        maske = mask
    3329         489 :        n     = count(maske)
    3330             :     else
    3331         489 :        maske = .true.
    3332           1 :        n     = size(x, dim=1) * size(x, dim=2) * size(x, dim=3)
    3333             :     end if
    3334             :     !
    3335           2 :     if (n .LE. 1_i4) stop 'wNSE_dp_3d: number of arguments must be at least 2'
    3336             :     ! Average of x
    3337             :     xmean = average(reshape(x(:,:,:), (/size(x, dim=1)*size(x, dim=2)*size(x, dim=3)/)), &
    3338           6 :          mask=reshape(maske(:,:,:), (/size(x, dim=1)*size(x, dim=2)*size(x, dim=3)/)))
    3339             :     !
    3340        1958 :     wNSE_dp_3d = 1.0_dp - sum(x * (y-x)*(y-x), mask=maske) / sum(x * (x-xmean)*(x-xmean), mask=maske)
    3341             :     !
    3342           2 :   END FUNCTION wNSE_dp_3d
    3343             : 
    3344             : END MODULE mo_errormeasures

Generated by: LCOV version 1.16