libsim Versione 7.2.6
|
◆ vol7d_normalize_vcoord()
Metodo per normalizzare la coordinata verticale. Per ora la normalizzazione effettuata riporta i valori di pressione nella sezione dati alla coordinata verticale sostituendo quella eventualmente presente. Classicamente serve per i dati con coordinata verticale model layer (105) Essendo che la pressione varia nello spazio orizzontale e nel tempo questo metodo restituisce un solo profilo verticale.
Definizione alla linea 1742 del file vol7d_class_compute.F90. 1743! Copyright (C) 2010 ARPA-SIM <urpsim@smr.arpa.emr.it>
1744! authors:
1745! Davide Cesari <dcesari@arpa.emr.it>
1746! Paolo Patruno <ppatruno@arpa.emr.it>
1747
1748! This program is free software; you can redistribute it and/or
1749! modify it under the terms of the GNU General Public License as
1750! published by the Free Software Foundation; either version 2 of
1751! the License, or (at your option) any later version.
1752
1753! This program is distributed in the hope that it will be useful,
1754! but WITHOUT ANY WARRANTY; without even the implied warranty of
1755! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
1756! GNU General Public License for more details.
1757
1758! You should have received a copy of the GNU General Public License
1759! along with this program. If not, see <http://www.gnu.org/licenses/>.
1760#include "config.h"
1761!> Extension of vol7d_class with methods for performing simple
1762!! statistical operations on entire volumes of data. This module
1763!! includes methods for performing operations such as (re)averaging
1764!! and (re)cumulation in time of entire vol7d_class::vol7d objects
1765!! containing real or double precision data.
1766!!
1767!! \ingroup vol7d
1773IMPLICIT NONE
1774
1775CONTAINS
1776
1777
1778!> General-purpose method for computing a statistical processing on
1779!! data in a vol7d object already processed with the same statistical
1780!! processing, on a different time interval specified by \a step and
1781!! \a start. This method tries to apply all the suitable specialized
1782!! statistical processing methods according to the input and output
1783!! statistical processing requested. The argument \a stat_proc_input
1784!! determines which data will be processed, while the \a stat_proc
1785!! argument determines the type of statistical process to be applied
1786!! and which will be owned by output data.
1787!!
1788!! The possible combinations are:
1789!!
1790!! - \a stat_proc_input = 254
1791!! - \a stat_proc = 0 average instantaneous observations
1792!! - \a stat_proc = 2 compute maximum of instantaneous observations
1793!! - \a stat_proc = 3 compute minimum of instantaneous observations
1794!! - \a stat_proc = 4 compute difference of instantaneous observations
1795!! - \a stat_proc = 6 compute standard deviation of instantaneous
1796!! observations
1797!! - \a stat_proc = 201 compute the prevailing direction (mode) on
1798!! specified sectors, suitable only for variables representing an
1799!! angle in degrees, e.g. wind direction
1800!!
1801!! processing is computed on longer intervals by aggregation, see the
1802!! description of vol7d_compute_stat_proc_agg()
1803!!
1804!! - \a stat_proc_input = *
1805!! - \a stat_proc = 254 consider statistically processed values as
1806!! instantaneous without any extra processing
1807!!
1808!! see the description of vol7d_decompute_stat_proc()
1809!!
1810!! - \a stat_proc_input = 0, 1, 2, 3, 4, 200
1811!! - \a stat_proc = \a stat_proc_input recompute input data on
1812!! different intervals
1813!!
1814!! the same statistical processing is applied to obtain data
1815!! processed on a different interval, either longer, by
1816!! aggregation, or shorter, by differences, see the description of
1817!! vol7d_recompute_stat_proc_agg() and
1818!! vol7d_recompute_stat_proc_diff() respectively; it is also
1819!! possible to provide \a stat_proc_input \a /= \a stat_proc, but
1820!! it has to be used with care.
1821!!
1822!! - \a stat_proc_input = 0
1823!! - \a stat_proc = 1
1824!!
1825!! a time-averaged rate or flux is transformed into a
1826!! time-integrated value (sometimes called accumulated) on the same
1827!! interval by multiplying the values by the length of the time
1828!! interval in seconds, keeping constant all the rest, including
1829!! the variable; the unit of the variable implicitly changes
1830!! accordingly, this is supported officially in grib2 standard, in
1831!! the other cases it is a forcing of the standards.
1832!!
1833!! - \a stat_proc_input = 1
1834!! - \a stat_proc = 0
1835!!
1836!! a time-integrated value (sometimes called accumulated) is
1837!! transformed into a time-averaged rate or flux on the same
1838!! interval by dividing the values by the length of the time
1839!! interval in seconds, see also the previous description of the
1840!! opposite computation.
1841!!
1842!! If a particular statistical processing cannot be performed on the
1843!! input data, the program continues with a warning and, if requested,
1844!! the input data is passed over to the volume specified by the \a
1845!! other argument, in order to allow continuation of processing. All
1846!! the other parameters are passed over to the specifical statistical
1847!! processing methods and are documented there.
1848SUBROUTINE vol7d_compute_stat_proc(this, that, stat_proc_input, stat_proc, &
1849 step, start, full_steps, frac_valid, max_step, weighted, other)
1850TYPE(vol7d),INTENT(inout) :: this !< volume providing data to be recomputed, it is not modified by the method, apart from performing a \a vol7d_alloc_vol on it
1851TYPE(vol7d),INTENT(out) :: that !< output volume which will contain the recomputed data
1852INTEGER,INTENT(in) :: stat_proc_input !< type of statistical processing of data that has to be processed (from grib2 table), only data having timerange of this type will be processed, the actual statistical processing performed and which will appear in the output volume, is however determined by \a stat_proc argument
1853INTEGER,INTENT(in) :: stat_proc !< type of statistical processing to be recomputed (from grib2 table), data in output volume \a that will have a timerange of this type
1854TYPE(timedelta),INTENT(in) :: step !< length of the step over which the statistical processing is performed
1855TYPE(datetime),INTENT(in),OPTIONAL :: start !< start of statistical processing interval
1856LOGICAL,INTENT(in),OPTIONAL :: full_steps !< if \a .TRUE. apply processing only on intervals starting at a forecast time or a reference time modulo \a step
1857REAL,INTENT(in),OPTIONAL :: frac_valid !< minimum fraction of valid data required for considering acceptable a recomputed value, default=1.
1858TYPE(timedelta),INTENT(in),OPTIONAL :: max_step !< maximum allowed distance in time between two contiguougs valid data within an interval, for the interval to be eligible for statistical processing
1859LOGICAL,INTENT(in),OPTIONAL :: weighted !< if provided and \c .TRUE., the statistical process is computed, if possible, by weighting every value with a weight proportional to its validity interval
1860TYPE(vol7d),INTENT(inout),OPTIONAL :: other !< optional volume that, on exit, is going to contain the data that did not contribute to the accumulation computation
1861
1862TYPE(vol7d) :: that1, that2, other1
1863INTEGER :: steps
1864
1865IF (stat_proc_input == 254) THEN
1866 CALL l4f_log(l4f_info, 'computing statistical processing by aggregation '//&
1868
1869 CALL vol7d_compute_stat_proc_agg(this, that, stat_proc, &
1870 step, start, full_steps, max_step, weighted, other)
1871
1872ELSE IF (stat_proc == 254) THEN
1873 CALL l4f_log(l4f_info, &
1874 'computing instantaneous data from statistically processed '//&
1876
1877! compute length of cumulation step in seconds
1879
1880 IF (any(this%timerange(:)%p2 == steps)) THEN ! data is ready
1881 CALL vol7d_decompute_stat_proc(this, that, step, other, stat_proc_input)
1882 ELSE
1883 IF (any(this%timerange(:)%p2 == steps/2)) THEN ! need to average
1884! average twice on step interval, with a shift of step/2, check full_steps
1885 CALL vol7d_recompute_stat_proc_agg(this, that1, stat_proc_input, &
1886 step, full_steps=.false., frac_valid=1.0)
1887 CALL vol7d_recompute_stat_proc_agg(this, that2, stat_proc_input, &
1888 step, start=that1%time(1)+step/2, full_steps=.false., frac_valid=1.0)
1889! merge the result
1891! and process it
1892 CALL vol7d_decompute_stat_proc(that1, that, step, other, stat_proc_input)
1895 ELSE
1896! do nothing
1897 ENDIF
1898 ENDIF
1899
1900ELSE IF (stat_proc_input == stat_proc .OR. &
1901 (stat_proc == 0 .OR. stat_proc == 2 .OR. stat_proc == 3)) THEN
1902! avg, min and max can be computed from any input, with care
1903 CALL l4f_log(l4f_info, &
1904 'recomputing statistically processed data by aggregation and difference '//&
1906
1907 IF (PRESENT(other)) THEN
1908 CALL vol7d_recompute_stat_proc_agg(this, that1, stat_proc, &
1909 step, start, full_steps, frac_valid, &
1910 other=other, stat_proc_input=stat_proc_input)
1911 CALL vol7d_recompute_stat_proc_diff(this, that2, stat_proc, &
1912 step, full_steps, start, other=other1)
1914 ELSE
1915 CALL vol7d_recompute_stat_proc_agg(this, that1, stat_proc, &
1916 step, start, full_steps, frac_valid, stat_proc_input=stat_proc_input)
1917 CALL vol7d_recompute_stat_proc_diff(this, that2, stat_proc, step, full_steps, &
1918 start)
1919 ENDIF
1920
1923 that = that1
1924ELSE ! IF (stat_proc_input /= stat_proc) THEN
1925 IF ((stat_proc_input == 0 .AND. stat_proc == 1) .OR. &
1926 (stat_proc_input == 1 .AND. stat_proc == 0)) THEN
1927 CALL l4f_log(l4f_info, &
1928 'computing statistically processed data by integration/differentiation '// &
1930 CALL vol7d_compute_stat_proc_metamorph(this, that, stat_proc_input, &
1931 stat_proc)
1932 ELSE
1933 CALL l4f_log(l4f_error, &
1935 ' not implemented or does not make sense')
1936 CALL raise_error()
1937 ENDIF
1938ENDIF
1939
1940END SUBROUTINE vol7d_compute_stat_proc
1941
1942
1943!> Specialized method for statistically processing a set of data
1944!! already processed with the same statistical processing, on a
1945!! different time interval. This method performs statistical
1946!! processing by aggregation of shorter intervals. Only floating
1947!! point single or double precision data are processed.
1948!!
1949!! The output \a that vol7d object contains elements from the original volume
1950!! \a this satisfying the conditions
1951!! - real single or double precision variables
1952!! - timerange (vol7d_timerange_class::vol7d_timerange::timerange)
1953!! of type \a stat_proc (or \a stat_proc_input if provided)
1954!! - any p1 (analysis/observation or forecast)
1955!! - p2 > 0 (processing interval non null, non instantaneous data)
1956!! and equal to a multiplier of \a step
1957!!
1958!! Output data will have timerange of type \a stat_proc and p2 = \a
1959!! step. The supported statistical processing methods (parameter \a
1960!! stat_proc) are:
1961!!
1962!! - 0 average
1963!! - 1 accumulation
1964!! - 2 maximum
1965!! - 3 minimum
1966!! - 4 difference
1967!! - 6 standard deviation
1968!! - 200 vectorial mean
1969!!
1970!! The start of processing period can be computed automatically from
1971!! the input intervals as the first possible interval modulo \a step,
1972!! or, for a better control, it can be specified explicitly by the
1973!! optional argument \a start. Notice that \a start indicates the
1974!! beginning of the processing interval, so in the final volume, the
1975!! first datum may have time equal to \a start \a + \a step, e.g. in
1976!! the case when time is the verification time, which is typical for
1977!! observed datasets.
1978!!
1979!! The purpose of the optional argument \a stat_proc_input is to allow
1980!! processing with a certain statistical processing operator a dataset
1981!! already processed with a different operator, by specifying the
1982!! latter as stat_proc_input; this is useful, for example, if one
1983!! wants to compute the monthly average of daily maximum temperatures;
1984!! however this has to be used with care since the resulting data
1985!! volume will not carry all the information about the processing
1986!! which has been done, in the previous case, for example, the
1987!! temperatures will simply look like monthly average temperatures.
1988SUBROUTINE vol7d_recompute_stat_proc_agg(this, that, stat_proc, &
1989 step, start, full_steps, frac_valid, other, stat_proc_input)
1990TYPE(vol7d),INTENT(inout) :: this !< volume providing data to be recomputed, it is not modified by the method, apart from performing a \a vol7d_alloc_vol on it
1991TYPE(vol7d),INTENT(out) :: that !< output volume which will contain the recomputed data
1992INTEGER,INTENT(in) :: stat_proc !< type of statistical processing to be recomputed (from grib2 table), only data having timerange of this type will be recomputed and will appear in the output volume
1993TYPE(timedelta),INTENT(in) :: step !< length of the step over which the statistical processing is performed
1994TYPE(datetime),INTENT(in),OPTIONAL :: start !< start of statistical processing interval
1995LOGICAL,INTENT(in),OPTIONAL :: full_steps !< if \a .TRUE. and \a start is not provided, apply processing only on intervals starting at a forecast time or a reference time modulo \a step
1996REAL,INTENT(in),OPTIONAL :: frac_valid !< minimum fraction of valid data required for considering acceptable a recomputed value, default=1.
1997TYPE(vol7d),INTENT(inout),OPTIONAL :: other !< optional volume that, on exit, is going to contain the data that did not contribute to the statistical processing
1998INTEGER,INTENT(in),OPTIONAL :: stat_proc_input !< to be used with care, type of statistical processing of data that has to be processed (from grib2 table), only data having timerange of this type will be recomputed, the actual statistical processing performed and which will appear in the output volume, is however determined by \a stat_proc argument
1999
2000INTEGER :: tri
2001INTEGER :: i, j, n, n1, ndtr, i1, i3, i5, i6
2002INTEGER :: linshape(1)
2003REAL :: lfrac_valid, frac_c, frac_m
2004LOGICAL,ALLOCATABLE :: ttr_mask(:,:)
2005TYPE(arrayof_ttr_mapper),POINTER :: map_ttr(:,:)
2006INTEGER,POINTER :: dtratio(:)
2007
2008
2009IF (PRESENT(stat_proc_input)) THEN
2010 tri = stat_proc_input
2011ELSE
2012 tri = stat_proc
2013ENDIF
2014IF (PRESENT(frac_valid)) THEN
2015 lfrac_valid = frac_valid
2016ELSE
2017 lfrac_valid = 1.0
2018ENDIF
2019
2020! be safe
2021CALL vol7d_alloc_vol(this)
2022! initial check
2023
2024! cleanup the original volume
2025CALL vol7d_smart_sort(this, lsort_time=.true.) ! time-ordered volume needed
2027
2029CALL vol7d_alloc(that, nana=SIZE(this%ana), nlevel=SIZE(this%level), &
2030 nnetwork=SIZE(this%network))
2031IF (ASSOCIATED(this%dativar%r)) THEN
2032 CALL vol7d_alloc(that, ndativarr=SIZE(this%dativar%r))
2033 that%dativar%r = this%dativar%r
2034ENDIF
2035IF (ASSOCIATED(this%dativar%d)) THEN
2036 CALL vol7d_alloc(that, ndativard=SIZE(this%dativar%d))
2037 that%dativar%d = this%dativar%d
2038ENDIF
2039that%ana = this%ana
2040that%level = this%level
2041that%network = this%network
2042
2043! compute the output time and timerange and all the required mappings
2044CALL recompute_stat_proc_agg_common(this%time, this%timerange, stat_proc, tri, &
2045 step, this%time_definition, that%time, that%timerange, map_ttr, dtratio, &
2046 start, full_steps)
2047CALL vol7d_alloc_vol(that)
2048
2049ALLOCATE(ttr_mask(SIZE(this%time), SIZE(this%timerange)))
2050linshape = (/SIZE(ttr_mask)/)
2051! finally perform computations
2052IF (ASSOCIATED(this%voldatir)) THEN
2053 DO j = 1, SIZE(that%timerange)
2054 DO i = 1, SIZE(that%time)
2055
2056 DO i1 = 1, SIZE(this%ana)
2057 DO i3 = 1, SIZE(this%level)
2058 DO i6 = 1, SIZE(this%network)
2059 DO i5 = 1, SIZE(this%dativar%r)
2060
2061 frac_m = 0.
2062 DO n1 = SIZE(dtratio), 1, -1 ! precedence to longer periods
2063 IF (dtratio(n1) <= 0) cycle ! safety check
2064 ttr_mask = .false.
2065 DO n = 1, map_ttr(i,j)%arraysize
2066 IF (map_ttr(i,j)%array(n)%extra_info == dtratio(n1)) THEN
2068 map_ttr(i,j)%array(n)%itr,i5,i6))) THEN
2069 ttr_mask(map_ttr(i,j)%array(n)%it, &
2070 map_ttr(i,j)%array(n)%itr) = .true.
2071 ENDIF
2072 ENDIF
2073 ENDDO
2074
2075 ndtr = count(ttr_mask)
2076 frac_c = real(ndtr)/real(dtratio(n1))
2077
2078 IF (ndtr > 0 .AND. frac_c >= max(lfrac_valid, frac_m)) THEN
2079 frac_m = frac_c
2080 SELECT CASE(stat_proc)
2081 CASE (0, 200) ! average, vectorial mean
2082 that%voldatir(i1,i,i3,j,i5,i6) = &
2083 sum(this%voldatir(i1,:,i3,:,i5,i6), &
2084 mask=ttr_mask)/ndtr
2085 CASE (1, 4) ! accumulation, difference
2086 that%voldatir(i1,i,i3,j,i5,i6) = &
2087 sum(this%voldatir(i1,:,i3,:,i5,i6), &
2088 mask=ttr_mask)
2089 CASE (2) ! maximum
2090 that%voldatir(i1,i,i3,j,i5,i6) = &
2091 maxval(this%voldatir(i1,:,i3,:,i5,i6), &
2092 mask=ttr_mask)
2093 CASE (3) ! minimum
2094 that%voldatir(i1,i,i3,j,i5,i6) = &
2095 minval(this%voldatir(i1,:,i3,:,i5,i6), &
2096 mask=ttr_mask)
2097 CASE (6) ! stddev
2098 that%voldatir(i1,i,i3,j,i5,i6) = &
2099 stat_stddev( &
2100 reshape(this%voldatir(i1,:,i3,:,i5,i6), shape=linshape), &
2101 mask=reshape(ttr_mask, shape=linshape))
2102 END SELECT
2103 ENDIF
2104
2105 ENDDO ! dtratio
2106 ENDDO ! var
2107 ENDDO ! network
2108 ENDDO ! level
2109 ENDDO ! ana
2111 ENDDO ! otime
2112 ENDDO ! otimerange
2113ENDIF
2114
2115IF (ASSOCIATED(this%voldatid)) THEN
2116 DO j = 1, SIZE(that%timerange)
2117 DO i = 1, SIZE(that%time)
2118
2119 DO i1 = 1, SIZE(this%ana)
2120 DO i3 = 1, SIZE(this%level)
2121 DO i6 = 1, SIZE(this%network)
2122 DO i5 = 1, SIZE(this%dativar%d)
2123
2124 frac_m = 0.
2125 DO n1 = SIZE(dtratio), 1, -1 ! precedence to longer periods
2126 IF (dtratio(n1) <= 0) cycle ! safety check
2127 ttr_mask = .false.
2128 DO n = 1, map_ttr(i,j)%arraysize
2129 IF (map_ttr(i,j)%array(n)%extra_info == dtratio(n1)) THEN
2131 map_ttr(i,j)%array(n)%itr,i5,i6))) THEN
2132 ttr_mask(map_ttr(i,j)%array(n)%it, &
2133 map_ttr(i,j)%array(n)%itr) = .true.
2134 ENDIF
2135 ENDIF
2136 ENDDO
2137
2138 ndtr = count(ttr_mask)
2139 frac_c = real(ndtr)/real(dtratio(n1))
2140
2141 IF (ndtr > 0 .AND. frac_c >= max(lfrac_valid, frac_m)) THEN
2142 frac_m = frac_c
2143 SELECT CASE(stat_proc)
2144 CASE (0) ! average
2145 that%voldatid(i1,i,i3,j,i5,i6) = &
2146 sum(this%voldatid(i1,:,i3,:,i5,i6), &
2147 mask=ttr_mask)/ndtr
2148 CASE (1, 4) ! accumulation, difference
2149 that%voldatid(i1,i,i3,j,i5,i6) = &
2150 sum(this%voldatid(i1,:,i3,:,i5,i6), &
2151 mask=ttr_mask)
2152 CASE (2) ! maximum
2153 that%voldatid(i1,i,i3,j,i5,i6) = &
2154 maxval(this%voldatid(i1,:,i3,:,i5,i6), &
2155 mask=ttr_mask)
2156 CASE (3) ! minimum
2157 that%voldatid(i1,i,i3,j,i5,i6) = &
2158 minval(this%voldatid(i1,:,i3,:,i5,i6), &
2159 mask=ttr_mask)
2160 CASE (6) ! stddev
2161 that%voldatid(i1,i,i3,j,i5,i6) = &
2162 stat_stddev( &
2163 reshape(this%voldatid(i1,:,i3,:,i5,i6), shape=linshape), &
2164 mask=reshape(ttr_mask, shape=linshape))
2165 END SELECT
2166 ENDIF
2167
2168 ENDDO ! dtratio
2169 ENDDO ! var
2170 ENDDO ! network
2171 ENDDO ! level
2172 ENDDO ! ana
2174 ENDDO ! otime
2175 ENDDO ! otimerange
2176ENDIF
2177
2178DEALLOCATE(ttr_mask)
2179
2180CALL makeother()
2181
2182CONTAINS
2183
2184SUBROUTINE makeother()
2185IF (PRESENT(other)) THEN ! create volume with the remaining data for further processing
2187 ltimerange=(this%timerange(:)%timerange /= tri .OR. this%timerange(:)%p2 == imiss &
2188 .OR. this%timerange(:)%p2 == 0)) ! or MOD(steps, this%timerange(:)%p2) == 0
2189ENDIF
2190END SUBROUTINE makeother
2191
2192END SUBROUTINE vol7d_recompute_stat_proc_agg
2193
2194
2195!> Method for statistically processing a set of instantaneous data.
2196!! This method performs statistical processing by aggregation of
2197!! instantaneous data. Only floating point single or double precision
2198!! data are processed.
2199!!
2200!! The output \a that vol7d object contains elements from the original volume
2201!! \a this satisfying the conditions
2202!! - real single or double precision variables
2203!! - timerange (vol7d_timerange_class::vol7d_timerange::timerange)
2204!! of type 254 (instantaeous)
2205!! - any p1 (analysis/observation or forecast)
2206!! - p2 = 0 (processing interval null, instantaneous data)
2207!!
2208!! Output data will have timerange of type \a stat_proc, and p2 = \a
2209!! step. The supported statistical processing methods (parameter \a
2210!! stat_proc) are:
2211!!
2212!! - 0 average
2213!! - 2 maximum
2214!! - 3 minimum
2215!! - 4 difference
2216!! - 6 standard deviation
2217!! - 201 mode (only for wind direction sectors)
2218!!
2219!! In the case of average, it is possible to weigh the data
2220!! proportionally to the length of the time interval for which every
2221!! single value is valid, i.e. halfway between the time level of the
2222!! value itself and the time of its nearest valid neighbours (argument
2223!! \a weighted). A maximum distance in time for input valid data can
2224!! be assigned with the optional argument \a max_step, in order to
2225!! filter datasets with too long "holes".
2226SUBROUTINE vol7d_compute_stat_proc_agg(this, that, stat_proc, &
2227 step, start, full_steps, max_step, weighted, other)
2228TYPE(vol7d),INTENT(inout) :: this !< volume providing data to be computed, it is not modified by the method, apart from performing a \a vol7d_alloc_vol on it
2229TYPE(vol7d),INTENT(out) :: that !< output volume which will contain the computed data
2230INTEGER,INTENT(in) :: stat_proc !< type of statistical processing to be computed (from grib2 table)
2231TYPE(timedelta),INTENT(in) :: step !< length of the step over which the statistical processing is performed
2232TYPE(datetime),INTENT(in),OPTIONAL :: start !< start of statistical processing interval
2233LOGICAL,INTENT(in),OPTIONAL :: full_steps !< if \a .TRUE. and \a start is not provided, apply processing only on intervals starting at a forecast time or a reference time modulo \a step
2234TYPE(timedelta),INTENT(in),OPTIONAL :: max_step !< maximum allowed distance in time between two contiguougs valid data within an interval, for the interval to be eligible for statistical processing
2235LOGICAL,INTENT(in),OPTIONAL :: weighted !< if provided and \c .TRUE., the statistical process is computed, if possible, by weighting every value with a weight proportional to its validity interval
2236TYPE(vol7d),INTENT(inout),OPTIONAL :: other !< optional volume that, on exit, is going to contain the data that did not contribute to the accumulation computation
2237!INTEGER,INTENT(in),OPTIONAL :: stat_proc_input !< to be used with care, type of statistical processing of data that has to be processed (from grib2 table), only data having timerange of this type will be recomputed, the actual statistical processing performed and which will appear in the output volume, is however determined by \a stat_proc argument
2238
2239TYPE(vol7d) :: v7dtmp
2240INTEGER :: tri
2241INTEGER :: i, j, n, ninp, ndtr, i1, i3, i5, i6, vartype, maxsize
2242TYPE(timedelta) :: lmax_step, act_max_step
2243TYPE(datetime) :: pstart, pend, reftime
2244TYPE(arrayof_ttr_mapper),POINTER :: map_ttr(:,:)
2245REAL,ALLOCATABLE :: tmpvolr(:)
2246DOUBLE PRECISION,ALLOCATABLE :: tmpvold(:), weights(:)
2247LOGICAL,ALLOCATABLE :: lin_mask(:)
2248LOGICAL :: lweighted
2249CHARACTER(len=8) :: env_var
2250
2251IF (PRESENT(max_step)) THEN
2252 lmax_step = max_step
2253ELSE
2254 lmax_step = timedelta_max
2255ENDIF
2256lweighted = optio_log(weighted)
2257tri = 254
2258! enable bad behavior for climat database
2259env_var = ''
2260CALL getenv('LIBSIM_CLIMAT_BEHAVIOR', env_var)
2261lweighted = lweighted .AND. len_trim(env_var) == 0
2262! only average can be weighted
2263lweighted = lweighted .AND. stat_proc == 0
2264
2265! be safe
2266CALL vol7d_alloc_vol(this)
2267! initial check
2268
2269! cleanup the original volume
2270CALL vol7d_smart_sort(this, lsort_time=.true.) ! time-ordered volume needed
2272! copy everything except time and timerange
2273CALL vol7d_copy(this, v7dtmp, ltime=(/.false./), ltimerange=(/.false./))
2274
2275! create new volume
2277! compute the output time and timerange and all the required mappings
2278CALL recompute_stat_proc_agg_common(this%time, this%timerange, stat_proc, tri, &
2279 step, this%time_definition, that%time, that%timerange, map_ttr, start=start, &
2280 full_steps=full_steps)
2281! merge with information from original volume
2282CALL vol7d_merge(that, v7dtmp)
2283
2284maxsize = maxval(map_ttr(:,:)%arraysize)
2285ALLOCATE(tmpvolr(maxsize), tmpvold(maxsize), lin_mask(maxsize), weights(maxsize))
2286do_otimerange: DO j = 1, SIZE(that%timerange)
2287 do_otime: DO i = 1, SIZE(that%time)
2288 ninp = map_ttr(i,j)%arraysize
2289 IF (ninp <= 0) cycle do_otime
2290! required for some computations
2291 CALL time_timerange_get_period(that%time(i), that%timerange(j), &
2292 that%time_definition, pstart, pend, reftime)
2293
2294 IF (ASSOCIATED(this%voldatir)) THEN
2295 DO i1 = 1, SIZE(this%ana)
2296 DO i3 = 1, SIZE(this%level)
2297 DO i6 = 1, SIZE(this%network)
2298 DO i5 = 1, SIZE(this%dativar%r)
2299! stat_proc difference treated separately here
2300 IF (stat_proc == 4) THEN
2301 IF (ninp >= 2) THEN
2302 IF (map_ttr(i,j)%array(1)%extra_info == 1 .AND. &
2303 map_ttr(i,j)%array(ninp)%extra_info == 2) THEN
2305 map_ttr(i,j)%array(1)%itr,i5,i6)) .AND. &
2306 c_e(this%voldatir(i1,map_ttr(i,j)%array(ninp)%it,i3, &
2307 map_ttr(i,j)%array(ninp)%itr,i5,i6))) THEN
2308 that%voldatir(i1,i,i3,j,i5,i6) = &
2309 this%voldatir(i1,map_ttr(i,j)%array(ninp)%it,i3, &
2310 map_ttr(i,j)%array(ninp)%itr,i5,i6) - &
2311 this%voldatir(i1,map_ttr(i,j)%array(1)%it,i3, &
2312 map_ttr(i,j)%array(1)%itr,i5,i6)
2313 ENDIF
2314 ENDIF
2315 ENDIF
2316 cycle
2317 ENDIF
2318! other stat_proc
2319 vartype = vol7d_vartype(this%dativar%r(i5))
2320 lin_mask = .false.
2321 ndtr = 0
2322 DO n = 1, ninp
2324 map_ttr(i,j)%array(n)%itr,i5,i6))) THEN
2325 ndtr = ndtr + 1
2326 tmpvolr(ndtr) = this%voldatir(i1,map_ttr(i,j)%array(n)%it,i3, &
2327 map_ttr(i,j)%array(n)%itr,i5,i6)
2328 lin_mask(n) = .true.
2329 ENDIF
2330 ENDDO
2331 IF (ndtr == 0) cycle
2332 IF (lweighted) THEN
2333 CALL compute_stat_proc_agg_sw(map_ttr(i,j)%array(1:ninp)%time, &
2334 pstart, pend, lin_mask(1:ninp), act_max_step, weights)
2335 ELSE
2336 CALL compute_stat_proc_agg_sw(map_ttr(i,j)%array(1:ninp)%time, &
2337 pstart, pend, lin_mask(1:ninp), act_max_step)
2338 ENDIF
2339 IF (act_max_step > lmax_step) cycle
2340
2341 SELECT CASE(stat_proc)
2342 CASE (0) ! average
2343 IF (lweighted) THEN
2344 that%voldatir(i1,i,i3,j,i5,i6) = &
2345 sum(real(weights(1:ndtr))*tmpvolr(1:ndtr))
2346 ELSE
2347 that%voldatir(i1,i,i3,j,i5,i6) = &
2348 sum(tmpvolr(1:ndtr))/ndtr
2349 ENDIF
2350 CASE (2) ! maximum
2351 that%voldatir(i1,i,i3,j,i5,i6) = &
2352 maxval(tmpvolr(1:ndtr))
2353 CASE (3) ! minimum
2354 that%voldatir(i1,i,i3,j,i5,i6) = &
2355 minval(tmpvolr(1:ndtr))
2356 CASE (6) ! stddev
2357 that%voldatir(i1,i,i3,j,i5,i6) = &
2358 stat_stddev(tmpvolr(1:ndtr))
2359 CASE (201) ! mode
2360! mode only for angles at the moment, with predefined histogram
2361 IF (vartype == var_dir360) THEN
2362! remove undefined wind direction (=0), improve check?
2363! and reduce to interval [22.5,382.5[
2364 WHERE (tmpvolr(1:ndtr) == 0.0)
2365 tmpvolr(1:ndtr) = rmiss
2366 ELSE WHERE (tmpvolr(1:ndtr) < 22.5 .AND. tmpvolr(1:ndtr) > 0.0)
2367 tmpvolr(1:ndtr) = tmpvolr(1:ndtr) + 360.
2368 END WHERE
2369 that%voldatir(i1,i,i3,j,i5,i6) = &
2370 stat_mode_histogram(tmpvolr(1:ndtr), &
2371 8, 22.5, 382.5)
2372 ENDIF
2373 END SELECT
2374 ENDDO
2375 ENDDO
2376 ENDDO
2377 ENDDO
2378 ENDIF
2379
2380 IF (ASSOCIATED(this%voldatid)) THEN
2381 DO i1 = 1, SIZE(this%ana)
2382 DO i3 = 1, SIZE(this%level)
2383 DO i6 = 1, SIZE(this%network)
2384 DO i5 = 1, SIZE(this%dativar%d)
2385! stat_proc difference treated separately here
2386 IF (stat_proc == 4) THEN
2387 IF (n >= 2) THEN
2388 IF (map_ttr(i,j)%array(1)%extra_info == 1 .AND. &
2389 map_ttr(i,j)%array(ninp)%extra_info == 2) THEN
2391 map_ttr(i,j)%array(1)%itr,i5,i6)) .AND. &
2392 c_e(this%voldatid(i1,map_ttr(i,j)%array(ninp)%it,i3, &
2393 map_ttr(i,j)%array(ninp)%itr,i5,i6))) THEN
2394 that%voldatid(i1,i,i3,j,i5,i6) = &
2395 this%voldatid(i1,map_ttr(i,j)%array(ninp)%it,i3, &
2396 map_ttr(i,j)%array(ninp)%itr,i5,i6) - &
2397 this%voldatid(i1,map_ttr(i,j)%array(1)%it,i3, &
2398 map_ttr(i,j)%array(1)%itr,i5,i6)
2399 ENDIF
2400 ENDIF
2401 ENDIF
2402 cycle
2403 ENDIF
2404! other stat_proc
2405 vartype = vol7d_vartype(this%dativar%d(i5))
2406 lin_mask = .false.
2407 ndtr = 0
2408 DO n = 1, ninp
2410 map_ttr(i,j)%array(n)%itr,i5,i6))) THEN
2411 ndtr = ndtr + 1
2412 tmpvold(ndtr) = this%voldatid(i1,map_ttr(i,j)%array(n)%it,i3, &
2413 map_ttr(i,j)%array(n)%itr,i5,i6)
2414 lin_mask(n) = .true.
2415 ENDIF
2416 ENDDO
2417 IF (ndtr == 0) cycle
2418 IF (lweighted) THEN
2419 CALL compute_stat_proc_agg_sw(map_ttr(i,j)%array(1:ninp)%time, &
2420 pstart, pend, lin_mask(1:ninp), act_max_step, weights)
2421 ELSE
2422 CALL compute_stat_proc_agg_sw(map_ttr(i,j)%array(1:ninp)%time, &
2423 pstart, pend, lin_mask(1:ninp), act_max_step)
2424 ENDIF
2425 IF (act_max_step > lmax_step) cycle
2426
2427 SELECT CASE(stat_proc)
2428 CASE (0) ! average
2429 IF (lweighted) THEN
2430 that%voldatid(i1,i,i3,j,i5,i6) = &
2431 sum(real(weights(1:ndtr))*tmpvold(1:ndtr))
2432 ELSE
2433 that%voldatid(i1,i,i3,j,i5,i6) = &
2434 sum(tmpvold(1:ndtr))/ndtr
2435 ENDIF
2436 CASE (2) ! maximum
2437 that%voldatid(i1,i,i3,j,i5,i6) = &
2438 maxval(tmpvold(1:ndtr))
2439 CASE (3) ! minimum
2440 that%voldatid(i1,i,i3,j,i5,i6) = &
2441 minval(tmpvold(1:ndtr))
2442 CASE (6) ! stddev
2443 that%voldatid(i1,i,i3,j,i5,i6) = &
2444 stat_stddev(tmpvold(1:ndtr))
2445 CASE (201) ! mode
2446! mode only for angles at the moment, with predefined histogram
2447 IF (vartype == var_dir360) THEN
2448! remove undefined wind direction (=0), improve check?
2449! and reduce to interval [22.5,382.5[
2450 WHERE (tmpvold(1:ndtr) == 0.0d0)
2451 tmpvold(1:ndtr) = dmiss
2452 ELSE WHERE (tmpvold(1:ndtr) < 22.5d0 .AND. tmpvold(1:ndtr) > 0.0d0)
2453 tmpvold(1:ndtr) = tmpvold(1:ndtr) + 360.0d0
2454 END WHERE
2455 that%voldatid(i1,i,i3,j,i5,i6) = &
2456 stat_mode_histogram(tmpvold(1:ndtr), &
2457 8, 22.5d0, 382.5d0)
2458 ENDIF
2459 END SELECT
2460 ENDDO
2461 ENDDO
2462 ENDDO
2463 ENDDO
2464 ENDIF
2465
2467 ENDDO do_otime
2468ENDDO do_otimerange
2469
2470DEALLOCATE(map_ttr)
2471DEALLOCATE(tmpvolr, tmpvold, lin_mask, weights)
2472
2473IF (PRESENT(other)) THEN ! create volume with the remaining data for further processing
2475 ltimerange=(this%timerange(:)%timerange /= tri))
2476ENDIF
2477
2478END SUBROUTINE vol7d_compute_stat_proc_agg
2479
2480
2481!> Method to transform the timerange of a set of data from
2482!! statistically processed to instantaneous. The data does not change,
2483!! only time and timerange descriptors change.
2484!!
2485!! The output \a that vol7d object contains elements from the original volume
2486!! \a this satisfying the conditions
2487!! - real single or double precision variables
2488!! - timerange (vol7d_timerange_class::vol7d_timerange::timerange)
2489!! of type 0 (average) or \a stat_proc_input if provided
2490!! - p1 = 0 (end of period == reference time, analysis/observation)
2491!! - p2 == \a step
2492!!
2493!! Output data will have timerange 254, p1 = 0 and p2 = 0; the time
2494!! dimension is shifted by half \a step so that it coincides with the
2495!! mid point of the input statistical processing interval.
2496SUBROUTINE vol7d_decompute_stat_proc(this, that, step, other, stat_proc_input)
2497TYPE(vol7d),INTENT(inout) :: this !< volume providing data to be recomputed, it is not modified by the method, apart from performing a \a vol7d_alloc_vol on it
2498TYPE(vol7d),INTENT(out) :: that !< output volume which will contain the recomputed data
2499TYPE(timedelta),INTENT(in) :: step !< length of the step over which the statistical processing is performed
2500TYPE(vol7d),INTENT(inout),OPTIONAL :: other !< optional volume that, on exit, is going to contain the data that did not contribute to the statistical processing
2501INTEGER,INTENT(in),OPTIONAL :: stat_proc_input !< type of statistical processing of data that has to be processed (from grib2 table), if not provided, averaged data (statistical processing = 0) is processed
2502
2503INTEGER :: i, tri, steps
2504
2505
2506IF (PRESENT(stat_proc_input)) THEN
2507 tri = stat_proc_input
2508ELSE
2509 tri = 0
2510ENDIF
2511! be safe
2512CALL vol7d_alloc_vol(this)
2513
2514! compute length of cumulation step in seconds
2516
2517! filter requested data
2519 ltimerange=(this%timerange(:)%timerange == tri .AND. &
2520 this%timerange(:)%p1 == 0 .AND. this%timerange(:)%p2 == steps))
2521
2522! convert timerange to instantaneous and go back half step in time
2523that%timerange(:)%timerange = 254
2524that%timerange(:)%p2 = 0
2525DO i = 1, SIZE(that%time(:))
2526 that%time(i) = that%time(i) - step/2
2527ENDDO
2528
2529IF (PRESENT(other)) THEN ! create volume with the remaining data for further processing
2531 ltimerange=(this%timerange(:)%timerange /= tri .OR. &
2532 this%timerange(:)%p1 /= 0 .OR. this%timerange(:)%p2 /= steps))
2533ENDIF
2534
2535END SUBROUTINE vol7d_decompute_stat_proc
2536
2537
2538!> Specialized method for statistically processing a set of data
2539!! already processed with the same statistical processing, on a
2540!! different time interval. This method performs statistical
2541!! processing by difference of different intervals. Only floating
2542!! point single or double precision data with analysis/observation
2543!! or forecast timerange are processed.
2544!!
2545!! The output \a that vol7d object contains elements from the original volume
2546!! \a this satisfying the conditions
2547!! - real single or double precision variables
2548!! - timerange (vol7d_timerange_class::vol7d_timerange::timerange)
2549!! of type \a stat_proc
2550!! - any p1 (analysis/observation or forecast)
2551!! - p2 > 0 (processing interval non null, non instantaneous data)
2552!! and equal to a multiplier of \a step if \a full_steps is \c .TRUE.
2553!!
2554!! Output data will have timerange of type \a stat_proc and p2 = \a
2555!! step. The supported statistical processing methods (parameter \a
2556!! stat_proc) are:
2557!!
2558!! - 0 average
2559!! - 1 cumulation
2560!! - 4 difference
2561!!
2562!! Input volume may have any value of \a this%time_definition, and
2563!! that value will be conserved in the output volume.
2564SUBROUTINE vol7d_recompute_stat_proc_diff(this, that, stat_proc, step, full_steps, start, other)
2565TYPE(vol7d),INTENT(inout) :: this !< volume providing data to be recomputed, it is not modified by the method, apart from performing a \a vol7d_alloc_vol on it
2566TYPE(vol7d),INTENT(out) :: that !< output volume which will contain the recomputed data
2567INTEGER,INTENT(in) :: stat_proc !< type of statistical processing to be recomputed (from grib2 table), only data having timerange of this type will be recomputed and will appear in the output volume
2568TYPE(timedelta),INTENT(in) :: step !< length of the step over which the statistical processing is performed
2569LOGICAL,INTENT(in),OPTIONAL :: full_steps !< if provided and \a .TRUE., process only data having processing interval (p2) equal to a multiplier of \a step
2570TYPE(datetime),INTENT(in),OPTIONAL :: start !< start of statistical processing interval
2571TYPE(vol7d),INTENT(out),OPTIONAL :: other !< optional volume that, on exit, is going to contain the data that did not contribute to the statistical processing
2572
2573INTEGER :: i1, i3, i5, i6, i, j, k, l, nitr, steps
2574INTEGER,ALLOCATABLE :: map_tr(:,:,:,:,:), f(:), keep_tr(:,:,:)
2575LOGICAL,ALLOCATABLE :: mask_timerange(:)
2576LOGICAL,ALLOCATABLE :: mask_time(:)
2577TYPE(vol7d) :: v7dtmp
2578
2579
2580! be safe
2581CALL vol7d_alloc_vol(this)
2582! initialise the template of the output volume
2584
2585! compute length of cumulation step in seconds
2587
2588! compute the statistical processing relations, output time and
2589! timerange are defined here
2590CALL recompute_stat_proc_diff_common(this%time, this%timerange, stat_proc, step, &
2591 that%time, that%timerange, map_tr, f, keep_tr, &
2592 this%time_definition, full_steps, start)
2593nitr = SIZE(f)
2594
2595! complete the definition of the empty output template
2596CALL vol7d_alloc(that, nana=0, nlevel=0, nnetwork=0)
2597CALL vol7d_alloc_vol(that)
2598! shouldn't we exit here with an empty volume if stat_proc/=0,1 ?
2599ALLOCATE(mask_time(SIZE(this%time)), mask_timerange(SIZE(this%timerange)))
2600DO l = 1, SIZE(this%time)
2601 mask_time(l) = any(this%time(l) == that%time(:))
2602ENDDO
2603DO l = 1, SIZE(this%timerange)
2604 mask_timerange(l) = any(this%timerange(l) == that%timerange(:))
2605ENDDO
2606! create template for the output volume, keep all ana, level, network
2607! and variables; copy only the timeranges already satisfying the
2608! requested step, if any and only the times already existing in the
2609! output
2611 ltimerange=mask_timerange(:), ltime=mask_time(:))
2612! merge output created so far with template
2613CALL vol7d_merge(that, v7dtmp, lanasimple=.true., llevelsimple=.true.)
2614
2615! compute statistical processing
2616IF (ASSOCIATED(this%voldatir)) THEN
2617 DO l = 1, SIZE(this%time)
2618 DO k = 1, nitr
2619 DO j = 1, SIZE(this%time)
2620 DO i = 1, nitr
2622 DO i6 = 1, SIZE(this%network)
2623 DO i5 = 1, SIZE(this%dativar%r)
2624 DO i3 = 1, SIZE(this%level)
2625 DO i1 = 1, SIZE(this%ana)
2628
2629 IF (stat_proc == 0) THEN ! average
2630 that%voldatir( &
2631 i1,map_tr(i,j,k,l,1),i3,map_tr(i,j,k,l,2),i5,i6) = &
2632 (this%voldatir(i1,l,i3,f(k),i5,i6)*this%timerange(f(k))%p2 - &
2633 this%voldatir(i1,j,i3,f(i),i5,i6)*this%timerange(f(i))%p2)/ &
2634 steps ! optimize avoiding conversions
2635 ELSE IF (stat_proc == 1 .OR. stat_proc == 4) THEN ! acc, diff
2636 that%voldatir( &
2637 i1,map_tr(i,j,k,l,1),i3,map_tr(i,j,k,l,2),i5,i6) = &
2638 this%voldatir(i1,l,i3,f(k),i5,i6) - &
2639 this%voldatir(i1,j,i3,f(i),i5,i6)
2640 ENDIF
2641
2642 ENDIF
2643 ENDDO
2644 ENDDO
2645 ENDDO
2646 ENDDO
2647 ENDIF
2648 ENDDO
2649 ENDDO
2650 ENDDO
2651 ENDDO
2652ENDIF
2653
2654IF (ASSOCIATED(this%voldatid)) THEN
2655 DO l = 1, SIZE(this%time)
2656 DO k = 1, nitr
2657 DO j = 1, SIZE(this%time)
2658 DO i = 1, nitr
2660 DO i6 = 1, SIZE(this%network)
2661 DO i5 = 1, SIZE(this%dativar%d)
2662 DO i3 = 1, SIZE(this%level)
2663 DO i1 = 1, SIZE(this%ana)
2666! IF (.NOT.c_e(that%voldatid( &
2667! i1,map_tr(i,j,k,l,1),i3,map_tr(i,j,k,l,2),i5,i6))) THEN
2668
2669 IF (stat_proc == 0) THEN ! average
2670 that%voldatid( &
2671 i1,map_tr(i,j,k,l,1),i3,map_tr(i,j,k,l,2),i5,i6) = &
2672 (this%voldatid(i1,l,i3,f(k),i5,i6)*this%timerange(f(k))%p2 - &
2673 this%voldatid(i1,j,i3,f(i),i5,i6)*this%timerange(f(i))%p2)/ &
2674 steps ! optimize avoiding conversions
2675 ELSE IF (stat_proc == 1 .OR. stat_proc == 4) THEN ! acc, diff
2676 that%voldatid( &
2677 i1,map_tr(i,j,k,l,1),i3,map_tr(i,j,k,l,2),i5,i6) = &
2678 this%voldatid(i1,l,i3,f(k),i5,i6) - &
2679 this%voldatid(i1,j,i3,f(i),i5,i6)
2680 ENDIF
2681
2682! ENDIF
2683 ENDIF
2684 ENDDO
2685 ENDDO
2686 ENDDO
2687 ENDDO
2688 ENDIF
2689 ENDDO
2690 ENDDO
2691 ENDDO
2692 ENDDO
2693ENDIF
2694
2695! this should be avoided by sorting descriptors upstream
2696! descriptors now are sorted upstream with a dirty and expensive trick
2697! but the order may be scrambled in the call to vol7d_merge above
2698CALL vol7d_smart_sort(that, lsort_time=.true., lsort_timerange=.true.)
2699
2700CALL makeother(.true.)
2701
2702CONTAINS
2703
2704SUBROUTINE makeother(filter)
2705LOGICAL,INTENT(in) :: filter
2706IF (PRESENT(other)) THEN
2707 IF (filter) THEN ! create volume with the remaining data for further processing
2709 ltimerange=(this%timerange(:)%timerange /= stat_proc))
2710 ELSE
2712 ENDIF
2713ENDIF
2714END SUBROUTINE makeother
2715
2716END SUBROUTINE vol7d_recompute_stat_proc_diff
2717
2718
2719!> Specialized method for statistically processing a set of data
2720!! by integration/differentiation.
2721!! This method performs statistical processing by integrating
2722!! (accumulating) in time values representing time-average rates or
2723!! fluxes, (\a stat_proc_input=0 \a stat_proc=1) or by transforming a
2724!! time-integrated (accumulated) value in a time-average rate or flux
2725!! (\a stat_proc_input=1 \a stat_proc=0). Analysis/observation or
2726!! forecast timeranges are processed. The only operation performed is
2727!! respectively multiplying or dividing the values by the length of
2728!! the time interval in seconds.
2729!!
2730!! The output \a that vol7d object contains elements from the original
2731!! volume \a this satisfying the conditions
2732!! - timerange (vol7d_timerange_class::vol7d_timerange::timerange)
2733!! of type \a stat_proc_input (0 or 1)
2734!! - any p1 (analysis/observation or forecast)
2735!! - p2 > 0 (processing interval non null, non instantaneous data)
2736!!
2737!! Output data will have timerange of type \a stat_proc (1 or 0) and
2738!! p1 and p2 equal to the corresponding input values. The supported
2739!! statistical processing methods (parameter \a stat_proc) are:
2740!!
2741!! - 0 average
2742!! - 1 cumulation
2743!!
2744!! Input volume may have any value of \a this%time_definition, and
2745!! that value will be conserved in the output volume.
2746SUBROUTINE vol7d_compute_stat_proc_metamorph(this, that, stat_proc_input, stat_proc)
2747TYPE(vol7d),INTENT(inout) :: this !< volume providing data to be recomputed, it is not modified by the method, apart from performing a \a vol7d_alloc_vol on it
2748TYPE(vol7d),INTENT(out) :: that !< output volume which will contain the recomputed data
2749INTEGER,INTENT(in) :: stat_proc_input !< type of statistical processing of data that has to be processed (from grib2 table), only data having timerange of this type will be processed, the actual statistical processing performed and which will appear in the output volume, is however determined by \a stat_proc argument
2750INTEGER,INTENT(in) :: stat_proc !< type of statistical processing to be recomputed (from grib2 table), data in output volume \a that will have a timerange of this type
2751
2752INTEGER :: j
2753LOGICAL,ALLOCATABLE :: tr_mask(:)
2754REAL,ALLOCATABLE :: int_ratio(:)
2755DOUBLE PRECISION,ALLOCATABLE :: int_ratiod(:)
2756
2757IF (.NOT.((stat_proc_input == 0 .AND. stat_proc == 1) .OR. &
2758 (stat_proc_input == 1 .AND. stat_proc == 0))) THEN
2759
2760 CALL l4f_log(l4f_warn, &
2761 'compute_stat_proc_metamorph, can only be applied to average->accumulated timerange and viceversa')
2762! return an empty volume, without signaling error
2764 CALL vol7d_alloc_vol(that)
2765 RETURN
2766ENDIF
2767
2768! be safe
2769CALL vol7d_alloc_vol(this)
2770
2771! useful timeranges
2772tr_mask = this%timerange(:)%timerange == stat_proc_input .AND. this%timerange(:)%p2 /= imiss &
2773 .AND. this%timerange(:)%p2 /= 0
2774
2775! initial check (necessary?)
2776IF (count(tr_mask) == 0) THEN
2777 CALL l4f_log(l4f_warn, &
2778 'vol7d_compute, no timeranges suitable for statistical processing by metamorphosis')
2780! CALL makeother()
2781 RETURN
2782ENDIF
2783
2784! copy required data and reset timerange
2785CALL vol7d_copy(this, that, ltimerange=tr_mask)
2786that%timerange(:)%timerange = stat_proc
2787! why next automatic f2003 allocation does not always work?
2788ALLOCATE(int_ratio(SIZE(that%timerange)), int_ratiod(SIZE(that%timerange)))
2789
2790IF (stat_proc == 0) THEN ! average -> integral
2791 int_ratio = 1./real(that%timerange(:)%p2)
2792 int_ratiod = 1./dble(that%timerange(:)%p2)
2793ELSE ! cumulation
2794 int_ratio = real(that%timerange(:)%p2)
2795 int_ratiod = dble(that%timerange(:)%p2)
2796ENDIF
2797
2798IF (ASSOCIATED(that%voldatir)) THEN
2799 DO j = 1, SIZE(that%timerange)
2801 that%voldatir(:,:,:,j,:,:) = that%voldatir(:,:,:,j,:,:)*int_ratio(j)
2802 ELSEWHERE
2803 that%voldatir(:,:,:,j,:,:) = rmiss
2804 END WHERE
2805 ENDDO
2806ENDIF
2807
2808IF (ASSOCIATED(that%voldatid)) THEN
2809 DO j = 1, SIZE(that%timerange)
2811 that%voldatid(:,:,:,j,:,:) = that%voldatid(:,:,:,j,:,:)*int_ratiod(j)
2812 ELSEWHERE
2813 that%voldatid(:,:,:,j,:,:) = rmiss
2814 END WHERE
2815 ENDDO
2816ENDIF
2817
2818
2819END SUBROUTINE vol7d_compute_stat_proc_metamorph
2820
2821
2822SUBROUTINE vol7d_recompute_stat_proc_agg_multiv(this, that, &
2823 step, start, frac_valid, multiv_proc)
2824TYPE(vol7d),INTENT(inout) :: this !< volume providing data to be recomputed, it is not modified by the method, apart from performing a \a vol7d_alloc_vol on it
2825TYPE(vol7d),INTENT(out) :: that !< output volume which will contain the recomputed data
2826!INTEGER,INTENT(in) :: stat_proc !< type of statistical processing to be recomputed (from grib2 table), only data having timerange of this type will be recomputed and will appear in the output volume
2827TYPE(timedelta),INTENT(in) :: step !< length of the step over which the statistical processing is performed
2828TYPE(datetime),INTENT(in),OPTIONAL :: start !< start of statistical processing interval
2829REAL,INTENT(in),OPTIONAL :: frac_valid !< minimum fraction of valid data required for considering acceptable a recomputed value, default=1.
2830!TYPE(vol7d),INTENT(inout),OPTIONAL :: other !< optional volume that, on exit, is going to contain the data that did not contribute to the statistical processing
2831!INTEGER,INTENT(in),OPTIONAL :: stat_proc_input !< to be used with care, type of statistical processing of data that has to be processed (from grib2 table), only data having timerange of this type will be recomputed, the actual statistical processing performed and which will appear in the output volume, is however determined by \a stat_proc argument
2832INTEGER,INTENT(in) :: multiv_proc !< index of multivariate specific operation
2833
2834INTEGER :: tri
2835INTEGER :: i, j, n, n1, ndtr, i1, i3, i5, i6
2836INTEGER :: linshape(1)
2837REAL :: lfrac_valid, frac_c, frac_m
2838LOGICAL,ALLOCATABLE :: ttr_mask(:,:)
2839TYPE(arrayof_ttr_mapper),POINTER :: map_ttr(:,:)
2840INTEGER,POINTER :: dtratio(:)
2841INTEGER :: stat_proc_input, stat_proc
2842
2843SELECT CASE(multiv_proc)
2844CASE (1) ! direction of maximum
2845 stat_proc_input = 205
2846 stat_proc=205
2847END SELECT
2848
2849tri = stat_proc_input
2850IF (PRESENT(frac_valid)) THEN
2851 lfrac_valid = frac_valid
2852ELSE
2853 lfrac_valid = 1.0
2854ENDIF
2855
2856! be safe
2857CALL vol7d_alloc_vol(this)
2858! initial check
2859
2860! cleanup the original volume
2861CALL vol7d_smart_sort(this, lsort_time=.true.) ! time-ordered volume needed
2863
2865CALL vol7d_alloc(that, nana=SIZE(this%ana), nlevel=SIZE(this%level), &
2866 nnetwork=SIZE(this%network))
2867IF (ASSOCIATED(this%dativar%r)) THEN
2868 CALL vol7d_alloc(that, ndativarr=SIZE(this%dativar%r))
2869 that%dativar%r = this%dativar%r
2870ENDIF
2871IF (ASSOCIATED(this%dativar%d)) THEN
2872 CALL vol7d_alloc(that, ndativard=SIZE(this%dativar%d))
2873 that%dativar%d = this%dativar%d
2874ENDIF
2875that%ana = this%ana
2876that%level = this%level
2877that%network = this%network
2878
2879! compute the output time and timerange and all the required mappings
2880CALL recompute_stat_proc_agg_common(this%time, this%timerange, stat_proc, tri, &
2881 step, this%time_definition, that%time, that%timerange, map_ttr, &
2882 dtratio=dtratio, start=start)
2883CALL vol7d_alloc_vol(that)
2884
2885ALLOCATE(ttr_mask(SIZE(this%time), SIZE(this%timerange)))
2886linshape = (/SIZE(ttr_mask)/)
2887! finally perform computations
2888IF (ASSOCIATED(this%voldatir)) THEN
2889 DO j = 1, SIZE(that%timerange)
2890 DO i = 1, SIZE(that%time)
2891
2892 DO i1 = 1, SIZE(this%ana)
2893 DO i3 = 1, SIZE(this%level)
2894 DO i6 = 1, SIZE(this%network)
2895 DO i5 = 1, SIZE(this%dativar%r)
2896
2897 frac_m = 0.
2898 DO n1 = SIZE(dtratio), 1, -1 ! precedence to longer periods
2899 IF (dtratio(n1) <= 0) cycle ! safety check
2900 ttr_mask = .false.
2901 DO n = 1, map_ttr(i,j)%arraysize
2902 IF (map_ttr(i,j)%array(n)%extra_info == dtratio(n1)) THEN
2904 map_ttr(i,j)%array(n)%itr,i5,i6))) THEN
2905 ttr_mask(map_ttr(i,j)%array(n)%it, &
2906 map_ttr(i,j)%array(n)%itr) = .true.
2907 ENDIF
2908 ENDIF
2909 ENDDO
2910
2911 ndtr = count(ttr_mask)
2912 frac_c = real(ndtr)/real(dtratio(n1))
2913
2914 IF (ndtr > 0 .AND. frac_c >= max(lfrac_valid, frac_m)) THEN
2915 frac_m = frac_c
2916 SELECT CASE(multiv_proc)
2917 CASE (1) ! average, vectorial mean
2918 that%voldatir(i1,i,i3,j,i5,i6) = &
2919 sum(this%voldatir(i1,:,i3,:,i5,i6), &
2920 mask=ttr_mask)/ndtr
2921 END SELECT
2922 ENDIF
2923
2924 ENDDO ! dtratio
2925 ENDDO ! var
2926 ENDDO ! network
2927 ENDDO ! level
2928 ENDDO ! ana
2930 ENDDO ! otime
2931 ENDDO ! otimerange
2932ENDIF
2933
2934IF (ASSOCIATED(this%voldatid)) THEN
2935 DO j = 1, SIZE(that%timerange)
2936 DO i = 1, SIZE(that%time)
2937
2938 DO i1 = 1, SIZE(this%ana)
2939 DO i3 = 1, SIZE(this%level)
2940 DO i6 = 1, SIZE(this%network)
2941 DO i5 = 1, SIZE(this%dativar%d)
2942
2943 frac_m = 0.
2944 DO n1 = SIZE(dtratio), 1, -1 ! precedence to longer periods
2945 IF (dtratio(n1) <= 0) cycle ! safety check
2946 ttr_mask = .false.
2947 DO n = 1, map_ttr(i,j)%arraysize
2948 IF (map_ttr(i,j)%array(n)%extra_info == dtratio(n1)) THEN
2950 map_ttr(i,j)%array(n)%itr,i5,i6))) THEN
2951 ttr_mask(map_ttr(i,j)%array(n)%it, &
2952 map_ttr(i,j)%array(n)%itr) = .true.
2953 ENDIF
2954 ENDIF
2955 ENDDO
2956
2957 ndtr = count(ttr_mask)
2958 frac_c = real(ndtr)/real(dtratio(n1))
2959
2960 IF (ndtr > 0 .AND. frac_c >= max(lfrac_valid, frac_m)) THEN
2961 frac_m = frac_c
2962 SELECT CASE(stat_proc)
2963 CASE (0) ! average
2964 that%voldatid(i1,i,i3,j,i5,i6) = &
2965 sum(this%voldatid(i1,:,i3,:,i5,i6), &
2966 mask=ttr_mask)/ndtr
2967 CASE (1, 4) ! accumulation, difference
2968 that%voldatid(i1,i,i3,j,i5,i6) = &
2969 sum(this%voldatid(i1,:,i3,:,i5,i6), &
2970 mask=ttr_mask)
2971 CASE (2) ! maximum
2972 that%voldatid(i1,i,i3,j,i5,i6) = &
2973 maxval(this%voldatid(i1,:,i3,:,i5,i6), &
2974 mask=ttr_mask)
2975 CASE (3) ! minimum
2976 that%voldatid(i1,i,i3,j,i5,i6) = &
2977 minval(this%voldatid(i1,:,i3,:,i5,i6), &
2978 mask=ttr_mask)
2979 CASE (6) ! stddev
2980 that%voldatid(i1,i,i3,j,i5,i6) = &
2981 stat_stddev( &
2982 reshape(this%voldatid(i1,:,i3,:,i5,i6), shape=linshape), &
2983 mask=reshape(ttr_mask, shape=linshape))
2984 END SELECT
2985 ENDIF
2986
2987 ENDDO ! dtratio
2988 ENDDO ! var
2989 ENDDO ! network
2990 ENDDO ! level
2991 ENDDO ! ana
2993 ENDDO ! otime
2994 ENDDO ! otimerange
2995ENDIF
2996
2997DEALLOCATE(ttr_mask)
2998
2999END SUBROUTINE vol7d_recompute_stat_proc_agg_multiv
3000
3001!> Riempimento dei buchi temporali in un volume.
3002!! Questo metodo crea, a partire da un volume originale, un nuovo
3003!! volume dati in cui la dimensione tempo contiene tutti gli istanti
3004!! tra \a start e \a stopp (o tra il primo e l'ultimo livello
3005!! temporale) ad intervalli \a step. Gli eventuali livelli mancanti
3006!! vengono aggiunti riempiendo le corrispondenti posizioni dei volumi
3007!! dati con valori mancanti. I livelli temporali che non sono ad
3008!! intervalli \a step interi a partire dall'inizio, oppure quelli che
3009!! giacciono fuori dall'intervallo \a start:stop non vengono toccati e
3010!! quindi rimangono immutati nel volume finale (si veda anche la
3011!! descrizione di ::vol7d_filter_time). Il volume originale non
3012!! viene modificato e quindi dovrà essere distrutto da parte del
3013!! programma chiamante se il suo contenuto non è più
3014!! richiesto. Attenzione, se necessario la dimensione tempo (vettore
3015!! \a this%time del volume \a this ) viene riordinata, come effetto
3016!! collaterale della chiamata.
3017SUBROUTINE vol7d_fill_time(this, that, step, start, stopp, cyclicdt)
3018TYPE(vol7d),INTENT(inout) :: this
3019TYPE(vol7d),INTENT(inout) :: that
3020TYPE(timedelta),INTENT(in) :: step
3021TYPE(datetime),INTENT(in),OPTIONAL :: start
3022TYPE(datetime),INTENT(in),OPTIONAL :: stopp
3023TYPE(cyclicdatetime),INTENT(in),OPTIONAL :: cyclicdt !< cyclic date and time
3024
3025TYPE(cyclicdatetime) :: lcyclicdt
3026TYPE(datetime) :: counter, lstart, lstop
3027INTEGER :: i, naddtime
3028
3029CALL safe_start_stop(this, lstart, lstop, start, stopp)
3031
3032lcyclicdt=cyclicdatetime_miss
3033if (present(cyclicdt)) then
3035end if
3036
3039
3040! Count the number of time levels required for completing the series
3041! valid also in the case (SIZE(this%time) == 0)
3042naddtime = 0
3043counter = lstart
3044i = 1
3045naddcount: DO WHILE(counter <= lstop)
3046 DO WHILE(i <= SIZE(this%time)) ! this%time(i) chases counter
3047 IF (counter < this%time(i)) THEN ! this%time(i) overtook counter
3048 i = max(i-1,1) ! go back if possible
3049 EXIT
3050 ELSE IF (counter == this%time(i) .OR. .NOT. counter == lcyclicdt) THEN ! found matching time
3051 counter = counter + step
3052 cycle naddcount
3053 ENDIF
3054 i = i + 1
3055 ENDDO
3056 naddtime = naddtime + 1
3057 counter = counter + step
3058ENDDO naddcount
3059
3060! old universal algorithm, not optimized, check that the new one is equivalent
3061!naddtime = 0
3062!counter = lstart
3063!DO WHILE(counter <= lstop)
3064! IF (.NOT.ANY(counter == this%time(:))) THEN
3065! naddtime = naddtime + 1
3066! ENDIF
3067! counter = counter + step
3068!ENDDO
3069
3070IF (naddtime > 0) THEN
3071
3073 CALL vol7d_alloc(that, ntime=naddtime)
3074 CALL vol7d_alloc_vol(that)
3075
3076 ! Repeat the count loop setting the time levels to be added
3077 naddtime = 0
3078 counter = lstart
3079 i = 1
3080 naddadd: DO WHILE(counter <= lstop)
3081 DO WHILE(i <= SIZE(this%time)) ! this%time(i) chases counter
3082 IF (counter < this%time(i)) THEN ! this%time(i) overtook counter
3083 i = max(i-1,1) ! go back if possible
3084 EXIT
3085 ELSE IF (counter == this%time(i) .OR. .NOT. counter == lcyclicdt) THEN ! found matching time
3086 counter = counter + step
3087 cycle naddadd
3088 ENDIF
3089 i = i + 1
3090 ENDDO
3091 naddtime = naddtime + 1
3092 that%time(naddtime) = counter ! only difference
3093 counter = counter + step
3094 ENDDO naddadd
3095
3097
3098ELSE
3099!! ? why sort all dimension ?
3100!! CALL vol7d_copy(this, that, lsort_time=.TRUE.)
3102ENDIF
3103
3104
3105END SUBROUTINE vol7d_fill_time
3106
3107
3108!> Filter time dimension inside a volume.
3109!! Questo metodo crea, a partire da un volume originale, un nuovo
3110!! volume dati in cui la dimensione tempo contiene solo gli
3111!! istanti tra \a start e \a stopp (o tra il primo e l'ultimo livello
3112!! temporale) ad intervalli \a step; se specificato cyclicdt solo i
3113!! corrispondenti istanti di tempo vengono ulteriormente selezionati.
3114!! Il volume originale non viene
3115!! modificato e quindi dovrà essere distrutto da parte del programma
3116!! chiamante se il suo contenuto non è più richiesto. Attenzione, se
3117!! necessario, la dimensione tempo (vettore \a this%time del volume \a
3118!! this ) viene riordinata, come effetto collaterale della chiamata.
3119SUBROUTINE vol7d_filter_time(this, that, step, start, stopp, cyclicdt)
3120TYPE(vol7d),INTENT(inout) :: this
3121TYPE(vol7d),INTENT(inout) :: that
3122TYPE(timedelta),INTENT(in),optional :: step !< missing value admitted
3123TYPE(datetime),INTENT(in),OPTIONAL :: start
3124TYPE(datetime),INTENT(in),OPTIONAL :: stopp
3125TYPE(cyclicdatetime),INTENT(in),OPTIONAL :: cyclicdt !< cyclic date and time
3126
3127TYPE(datetime) :: lstart, lstop
3128LOGICAL, ALLOCATABLE :: time_mask(:)
3129
3130CALL safe_start_stop(this, lstart, lstop, start, stopp)
3132
3135
3136ALLOCATE(time_mask(SIZE(this%time)))
3137
3138time_mask = this%time >= lstart .AND. this%time <= lstop
3139
3140IF (PRESENT(cyclicdt)) THEN
3142 time_mask = time_mask .AND. this%time == cyclicdt
3143 ENDIF
3144ENDIF
3145
3146IF (PRESENT(step)) THEN
3148 time_mask = time_mask .AND. mod(this%time - lstart, step) == timedelta_0
3149 ENDIF
3150ENDIF
3151
3152CALL vol7d_copy(this,that, ltime=time_mask)
3153
3154DEALLOCATE(time_mask)
3155
3156END SUBROUTINE vol7d_filter_time
3157
3158
3159!> Fill data volume
3160!! Nearest data in time is set in the time coordinate.
3161!! Take in account istantaneous values only.
3162SUBROUTINE vol7d_fill_data(this, step, start, stopp, tolerance)
3163TYPE(vol7d),INTENT(inout) :: this !< data volume to elaborate
3164TYPE(timedelta),INTENT(in) :: step !< interval in time where to fill data
3165TYPE(datetime),INTENT(in),OPTIONAL :: start !< start time where to fill
3166TYPE(datetime),INTENT(in),OPTIONAL :: stopp !< stop time where to fill
3167TYPE(timedelta),INTENT(in),optional :: tolerance !< tolerance in time to find data to fill (excluding extreme) (default to step)
3168
3169TYPE(datetime) :: lstart, lstop
3170integer :: indana , indtime ,indlevel ,indtimerange ,inddativarr, indnetwork, iindtime
3171type(timedelta) :: deltato,deltat, ltolerance
3172
3173CALL safe_start_stop(this, lstart, lstop, start, stopp)
3175
3178
3179
3180ltolerance=step/2
3181
3182if (present(tolerance))then
3184end if
3185
3186
3187do indtime=1,size(this%time)
3188
3189 IF (this%time(indtime) < lstart .OR. this%time(indtime) > lstop .OR. &
3190 mod(this%time(indtime) - lstart, step) /= timedelta_0) cycle
3191 do indtimerange=1,size(this%timerange)
3192 if (this%timerange(indtimerange)%timerange /= 254) cycle
3193 do indnetwork=1,size(this%network)
3194 do inddativarr=1,size(this%dativar%r)
3195 do indlevel=1,size(this%level)
3196 do indana=1,size(this%ana)
3197
3198 !find the nearest data in time if data is missing
3199 if (.not. c_e(this%voldatir(indana, indtime, indlevel, indtimerange, inddativarr, indnetwork)))then
3200 deltato=timedelta_miss
3201
3202 !do iindtime=max(indtime-20,1),min(indtime+20,size(this%time)) !check on a chunk: 20 should be enought
3203
3204 do iindtime=indtime+1,size(this%time) !check forward
3205
3206 if (c_e(this%voldatir (indana, iindtime, indlevel, indtimerange, inddativarr, indnetwork )))then
3207 deltat=this%time(iindtime)-this%time(indtime)
3208
3209 if (deltat >= ltolerance) exit
3210
3211 if (deltat < deltato) then
3212 this%voldatir(indana, indtime, indlevel, indtimerange, inddativarr, indnetwork) = &
3213 this%voldatir(indana, iindtime, indlevel, indtimerange, inddativarr, indnetwork)
3214 deltato=deltat
3215 end if
3216 end if
3217 end do
3218
3219 do iindtime=indtime-1,1,-1 !check backward
3220
3221 if (c_e(this%voldatir (indana, iindtime, indlevel, indtimerange, inddativarr, indnetwork )))then
3222 if (iindtime < indtime) then
3223 deltat=this%time(indtime)-this%time(iindtime)
3224 else if (iindtime > indtime) then
3225 deltat=this%time(iindtime)-this%time(indtime)
3226 else
3227 cycle
3228 end if
3229
3230 if (deltat >= ltolerance) exit
3231
3232 if (deltat < deltato) then
3233 this%voldatir(indana, indtime, indlevel, indtimerange, inddativarr, indnetwork) = &
3234 this%voldatir(indana, iindtime, indlevel, indtimerange, inddativarr, indnetwork)
3235 deltato=deltat
3236 end if
3237 end if
3238 end do
3239
3240 end if
3241 end do
3242 end do
3243 end do
3244 end do
3245 end do
3246end do
3247
3248END SUBROUTINE vol7d_fill_data
3249
3250
3251! private utility routine for checking interval and start-stop times
3252! in input missing start-stop values are treated as not present
3253! in output missing start-stop values mean "do nothing"
3254SUBROUTINE safe_start_stop(this, lstart, lstop, start, stopp)
3255TYPE(vol7d),INTENT(inout) :: this
3256TYPE(datetime),INTENT(out) :: lstart
3257TYPE(datetime),INTENT(out) :: lstop
3258TYPE(datetime),INTENT(in),OPTIONAL :: start
3259TYPE(datetime),INTENT(in),OPTIONAL :: stopp
3260
3261lstart = datetime_miss
3262lstop = datetime_miss
3263! initial safety operation
3264CALL vol7d_alloc_vol(this)
3265IF (SIZE(this%time) == 0) RETURN ! avoid segmentation fault in case of empty volume
3266CALL vol7d_smart_sort(this, lsort_time=.true.)
3267
3268IF (PRESENT(start)) THEN
3270 lstart = start
3271 ELSE
3272 lstart = this%time(1)
3273 ENDIF
3274ELSE
3275 lstart = this%time(1)
3276ENDIF
3277IF (PRESENT(stopp)) THEN
3279 lstop = stopp
3280 ELSE
3281 lstop = this%time(SIZE(this%time))
3282 ENDIF
3283ELSE
3284 lstop = this%time(SIZE(this%time))
3285ENDIF
3286
3287END SUBROUTINE safe_start_stop
3288
3289
3290!> Metodo per normalizzare la coordinata verticale.
3291!! Per ora la normalizzazione effettuata riporta i valori di pressione
3292!! nella sezione dati alla coordinata verticale sostituendo quella eventualmente presente.
3293!! Classicamente serve per i dati con coordinata verticale model layer (105)
3294!! Essendo che la pressione varia nello spazio orizzontale e nel tempo
3295!! questo metodo restituisce un solo profilo verticale.
3296SUBROUTINE vol7d_normalize_vcoord(this,that,ana,time,timerange,network)
3297TYPE(vol7d),INTENT(INOUT) :: this !< oggetto da normalizzare
3298TYPE(vol7d),INTENT(OUT) :: that !< oggetto normalizzato
3299integer,intent(in) :: time,ana,timerange,network !< indici dell'elemento da estrarre
3300
3301character(len=1) :: type
3302integer :: ind
3303TYPE(vol7d_var) :: var
3304LOGICAL,allocatable :: ltime(:),ltimerange(:),lana(:),lnetwork(:)
3305logical,allocatable :: maschera(:)
3306
3307
3308allocate(ltime(size(this%time)))
3309allocate(ltimerange(size(this%timerange)))
3310allocate(lana(size(this%ana)))
3311allocate(lnetwork(size(this%network)))
3312
3313ltime=.false.
3314ltimerange=.false.
3315lana=.false.
3316lnetwork=.false.
3317
3318ltime(time)=.true.
3319ltimerange(timerange)=.true.
3320lana(ana)=.true.
3321lnetwork(network)=.true.
3322
3323call vol7d_copy(this, that,unique=.true.,&
3324 ltime=ltime,ltimerange=ltimerange,lana=lana,lnetwork=lnetwork )
3325
3327type=cmiss
3328!type="i"
3329ind = index(that%dativar, var, type=type)
3330
3331allocate(maschera(size(that%level)))
3332
3333maschera = (&
3334 (that%level%level1 == 105.and.that%level%level2 == 105) .or. &
3335 (that%level%level1 == 103 .and. that%level%level2 == imiss ) .or. &
3336 (that%level%level1 == 102 .and. that%level%level2 == imiss )) &
3337 .and. c_e(that%voldatic(1,1,:,1,ind,1))
3338
3339
3340select case (type)
3341
3342case("d")
3343
3344 where (maschera)
3345 that%level%level1 = 100
3346 that%level%l1 = int(realdat(that%voldatid(1,1,:,1,ind,1),that%dativar%d(ind)))
3347 that%level%l1 = int(that%voldatid(1,1,:,1,ind,1))
3348 that%level%level2 = imiss
3349 that%level%l2 = imiss
3350 end where
3351
3352case("r")
3353
3354 where (maschera)
3355 that%level%level1 = 100
3356 that%level%l1 = int(realdat(that%voldatir(1,1,:,1,ind,1),that%dativar%r(ind)))
3357 that%level%level2 = imiss
3358 that%level%l2 = imiss
3359 end where
3360
3361case("i")
3362
3363 where (maschera)
3364 that%level%level1 = 100
3365 that%level%l1 = int(realdat(that%voldatii(1,1,:,1,ind,1),that%dativar%i(ind)))
3366 that%level%level2 = imiss
3367 that%level%l2 = imiss
3368 end where
3369
3370case("b")
3371
3372 where (maschera)
3373 that%level%level1 = 100
3374 that%level%l1 = int(realdat(that%voldatib(1,1,:,1,ind,1),that%dativar%b(ind)))
3375 that%level%level2 = imiss
3376 that%level%l2 = imiss
3377 end where
3378
3379case("c")
3380
3381 where (maschera)
3382 that%level%level1 = 100
3383 that%level%l1 = int(realdat(that%voldatic(1,1,:,1,ind,1),that%dativar%c(ind)))
3384 that%level%level2 = imiss
3385 that%level%l2 = imiss
3386 end where
3387
3388end select
3389
3390deallocate(ltime)
3391deallocate(ltimerange)
3392deallocate(lana)
3393deallocate(lnetwork)
3394
3395END SUBROUTINE vol7d_normalize_vcoord
3396
3397
3398!!$!> Metodo per calcolare variabili derivate.
3399!!$!! TO DO !!
3400!!$SUBROUTINE vol7d_compute_var(this,that,var)
3401!!$TYPE(vol7d),INTENT(INOUT) :: this !< oggetto da normalizzare
3402!!$TYPE(vol7d),INTENT(OUT) :: that !< oggetto normalizzato
3403!!$
3404!!$character(len=1) :: type
3405!!$TYPE(vol7d_var),intent(in) :: var
3406
3407
3408!!$call init(var, btable="B10004") ! Pressure
3409!!$type=cmiss
3410!!$call vol7d_varvect_index(that%dativar,var , type=type,index_v=ind)
3411!!$
3412!!$select case (type)
3413!!$
3414!!$case("d")
3415!!$
3416!!$ where (that%level%level1 == 105.and.that%level%level2 == 105)
3417!!$ that%level%level1 = 100
3418!!$ that%level%l1 = realdat(that%voldatid(1,1,:,1,ind,1),that%dativar%d(ind))
3419!!$ that%level%level2 = imiss
3420!!$ that%level%l2 = imiss
3421!!$ end where
3422!!$
3423!!$case("r")
3424!!$
3425!!$ where (that%level%level1 == 105.and.that%level%level2 == 105)
3426!!$ that%level%level1 = 100
3427!!$ that%level%l1 = realdat(that%voldatir(1,1,:,1,ind,1),that%dativar%r(ind))
3428!!$ that%level%level2 = imiss
3429!!$ that%level%l2 = imiss
3430!!$ end where
3431!!$
3432!!$case("i")
3433!!$
3434!!$ where (that%level%level1 == 105.and.that%level%level2 == 105)
3435!!$ that%level%level1 = 100
3436!!$ that%level%l1 = realdat(that%voldatii(1,1,:,1,ind,1),that%dativar%i(ind))
3437!!$ that%level%level2 = imiss
3438!!$ that%level%l2 = imiss
3439!!$ end where
3440!!$
3441!!$case("b")
3442!!$
3443!!$ where (that%level%level1 == 105.and.that%level%level2 == 105)
3444!!$ that%level%level1 = 100
3445!!$ that%level%l1 = realdat(that%voldatib(1,1,:,1,ind,1),that%dativar%b(ind))
3446!!$ that%level%level2 = imiss
3447!!$ that%level%l2 = imiss
3448!!$ end where
3449!!$
3450!!$case("c")
3451!!$
3452!!$ where (that%level%level1 == 105.and.that%level%level2 == 105)
3453!!$ that%level%level1 = 100
3454!!$ that%level%l1 = realdat(that%voldatic(1,1,:,1,ind,1),that%dativar%c(ind))
3455!!$ that%level%level2 = imiss
3456!!$ that%level%l2 = imiss
3457!!$ end where
3458!!$
3459!!$end select
3460
3461!!$
3462!!$END SUBROUTINE vol7d_compute_var
3463!!$
3464
3465
3466
Restituiscono il valore dell'oggetto nella forma desiderata. Definition datetime_class.F90:322 Costruttori per le classi datetime e timedelta. Definition datetime_class.F90:311 Functions that return a trimmed CHARACTER representation of the input variable. Definition datetime_class.F90:349 Restituiscono il valore dell'oggetto in forma di stringa stampabile. Definition datetime_class.F90:327 Compute the mode of the random variable provided taking into account missing data. Definition simple_stat.f90:160 Compute the standard deviation of the random variable provided, taking into account missing data. Definition simple_stat.f90:69 Module for basic statistical computations taking into account missing data. Definition simple_stat.f90:25 This module contains functions that are only for internal use of the library. Definition stat_proc_engine.F90:211 Extension of vol7d_class with methods for performing simple statistical operations on entire volumes ... Definition vol7d_class_compute.F90:214 Classe per la gestione di un volume completo di dati osservati. Definition vol7d_class.F90:273 |