libsim Versione 7.1.11

◆ conv_func_convert()

elemental real function conv_func_convert ( type(conv_func), intent(in)  this,
real, intent(in)  values 
)

Return a copy of values converted by applying the conversion function this.

The numerical conversion (only linear at the moment) defined by the conv_func object this is applied to the values argument and the converted result is returned; missing values remain missing; if the conversion function is undefined (conv_func_miss) the values are unchanged. The method is ELEMENTAL, thus values can be also an array of any shape.

Parametri
[in]thisobject defining the conversion function
[in]valuesinput value to be converted

Definizione alla linea 1481 del file volgrid6d_var_class.F90.

1482! Copyright (C) 2010 ARPA-SIM <urpsim@smr.arpa.emr.it>
1483! authors:
1484! Davide Cesari <dcesari@arpa.emr.it>
1485! Paolo Patruno <ppatruno@arpa.emr.it>
1486
1487! This program is free software; you can redistribute it and/or
1488! modify it under the terms of the GNU General Public License as
1489! published by the Free Software Foundation; either version 2 of
1490! the License, or (at your option) any later version.
1491
1492! This program is distributed in the hope that it will be useful,
1493! but WITHOUT ANY WARRANTY; without even the implied warranty of
1494! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
1495! GNU General Public License for more details.
1496
1497! You should have received a copy of the GNU General Public License
1498! along with this program. If not, see <http://www.gnu.org/licenses/>.
1499#include "config.h"
1500
1512USE kinds
1514USE err_handling
1517USE grid_id_class
1518
1519IMPLICIT NONE
1520
1525TYPE volgrid6d_var
1526 integer :: centre
1527 integer :: category
1528 integer :: number
1529 integer :: discipline
1530 CHARACTER(len=65) :: description
1531 CHARACTER(len=24) :: unit
1532END TYPE volgrid6d_var
1533
1534TYPE(volgrid6d_var),PARAMETER :: volgrid6d_var_miss= &
1535 volgrid6d_var(imiss,imiss,imiss,imiss,cmiss,cmiss)
1536
1537TYPE(vol7d_var),PARAMETER :: vol7d_var_horstag(2) = (/ &
1538 vol7d_var('B11003', '', '', 0, 0, 0, 0, 0, 0), &
1539 vol7d_var('B11004', '', '', 0, 0, 0, 0, 0, 0) &
1540 /)
1541
1542TYPE(vol7d_var),PARAMETER :: vol7d_var_horcomp(4) = (/ &! RESHAPE( (/ &
1543 vol7d_var('B11003', '', '', 0, 0, 0, 0, 0, 0), &
1544 vol7d_var('B11004', '', '', 0, 0, 0, 0, 0, 0), &
1545 vol7d_var('B11200', '', '', 0, 0, 0, 0, 0, 0), &
1546 vol7d_var('B11201', '', '', 0, 0, 0, 0, 0, 0) &
1547/)
1548!/), (/2,2/)) ! bug in gfortran
1549
1558TYPE conv_func
1559 PRIVATE
1560 REAL :: a, b
1561END TYPE conv_func
1562
1563TYPE(conv_func), PARAMETER :: conv_func_miss=conv_func(rmiss,rmiss)
1564TYPE(conv_func), PARAMETER :: conv_func_identity=conv_func(1.0,0.0)
1565
1566TYPE vg6d_v7d_var_conv
1567 TYPE(volgrid6d_var) :: vg6d_var
1568 TYPE(vol7d_var) :: v7d_var
1569 TYPE(conv_func) :: c_func
1570! aggiungere informazioni ad es. su rotazione del vento
1571END TYPE vg6d_v7d_var_conv
1572
1573TYPE(vg6d_v7d_var_conv), PARAMETER :: vg6d_v7d_var_conv_miss= &
1574 vg6d_v7d_var_conv(volgrid6d_var_miss, vol7d_var_miss, conv_func_miss)
1575
1576TYPE(vg6d_v7d_var_conv), ALLOCATABLE :: conv_fwd(:), conv_bwd(:)
1577
1591INTERFACE init
1592 MODULE PROCEDURE volgrid6d_var_init
1593END INTERFACE
1594
1597INTERFACE delete
1598 MODULE PROCEDURE volgrid6d_var_delete
1599END INTERFACE
1600
1601INTERFACE c_e
1602 MODULE PROCEDURE volgrid6d_var_c_e
1603END INTERFACE
1604
1605
1610INTERFACE OPERATOR (==)
1611 MODULE PROCEDURE volgrid6d_var_eq, conv_func_eq
1612END INTERFACE
1613
1618INTERFACE OPERATOR (/=)
1619 MODULE PROCEDURE volgrid6d_var_ne, conv_func_ne
1620END INTERFACE
1621
1622#define VOL7D_POLY_TYPE TYPE(volgrid6d_var)
1623#define VOL7D_POLY_TYPES _var6d
1624#include "array_utilities_pre.F90"
1625
1627INTERFACE display
1628 MODULE PROCEDURE display_volgrid6d_var
1629END INTERFACE
1630
1635INTERFACE OPERATOR (*)
1636 MODULE PROCEDURE conv_func_mult
1637END INTERFACE OPERATOR (*)
1638
1641INTERFACE compute
1642 MODULE PROCEDURE conv_func_compute
1643END INTERFACE
1644
1647INTERFACE convert
1648 MODULE PROCEDURE varbufr2vargrib_convert, vargrib2varbufr_convert, &
1649 conv_func_convert
1650END INTERFACE
1651
1652PRIVATE
1653PUBLIC volgrid6d_var, volgrid6d_var_miss, volgrid6d_var_new, init, delete, &
1654 c_e, volgrid6d_var_normalize, &
1655 OPERATOR(==), OPERATOR(/=), OPERATOR(*), &
1656 count_distinct, pack_distinct, count_and_pack_distinct, &
1657 map_distinct, map_inv_distinct, &
1658 index, display, &
1659 vargrib2varbufr, varbufr2vargrib, &
1660 conv_func, conv_func_miss, compute, convert, &
1661 volgrid6d_var_hor_comp_index, volgrid6d_var_is_hor_comp
1662
1663
1664CONTAINS
1665
1666
1667ELEMENTAL FUNCTION volgrid6d_var_new(centre, category, number, &
1668 discipline, description, unit) RESULT(this)
1669integer,INTENT(in),OPTIONAL :: centre
1670integer,INTENT(in),OPTIONAL :: category
1671integer,INTENT(in),OPTIONAL :: number
1672integer,INTENT(in),OPTIONAL :: discipline
1673CHARACTER(len=*),INTENT(in),OPTIONAL :: description
1674CHARACTER(len=*),INTENT(in),OPTIONAL :: unit
1675
1676TYPE(volgrid6d_var) :: this
1677
1678CALL init(this, centre, category, number, discipline, description, unit)
1679
1680END FUNCTION volgrid6d_var_new
1681
1682
1683! documented in the interface
1684ELEMENTAL SUBROUTINE volgrid6d_var_init(this, centre, category, number, discipline,description,unit)
1685TYPE(volgrid6d_var),INTENT(INOUT) :: this ! object to be initialized
1686INTEGER,INTENT(in),OPTIONAL :: centre ! centre
1687INTEGER,INTENT(in),OPTIONAL :: category ! grib2: category / grib1: grib table version number
1688INTEGER,INTENT(in),OPTIONAL :: number ! parameter number
1689INTEGER,INTENT(in),OPTIONAL :: discipline ! grib2: discipline / grib1: 255
1690CHARACTER(len=*),INTENT(in),OPTIONAL :: description ! optional textual description of the variable
1691CHARACTER(len=*),INTENT(in),OPTIONAL :: unit ! optional textual description of the variable's unit
1692
1693IF (PRESENT(centre)) THEN
1694 this%centre = centre
1695ELSE
1696 this%centre = imiss
1697 this%category = imiss
1698 this%number = imiss
1699 this%discipline = imiss
1700 RETURN
1701ENDIF
1702
1703IF (PRESENT(category)) THEN
1704 this%category = category
1705ELSE
1706 this%category = imiss
1707 this%number = imiss
1708 this%discipline = imiss
1709 RETURN
1710ENDIF
1711
1712
1713IF (PRESENT(number)) THEN
1714 this%number = number
1715ELSE
1716 this%number = imiss
1717 this%discipline = imiss
1718 RETURN
1719ENDIF
1720
1721! se sono arrivato fino a qui ho impostato centre, category e number
1722!per il grib 1 manca discipline e imposto 255 (missing del grib2)
1723
1724IF (PRESENT(discipline)) THEN
1725 this%discipline = discipline
1726ELSE
1727 this%discipline = 255
1728ENDIF
1729
1730IF (PRESENT(description)) THEN
1731 this%description = description
1732ELSE
1733 this%description = cmiss
1734ENDIF
1735
1736IF (PRESENT(unit)) THEN
1737 this%unit = unit
1738ELSE
1739 this%unit = cmiss
1740ENDIF
1741
1742
1743
1744END SUBROUTINE volgrid6d_var_init
1745
1746
1747! documented in the interface
1748SUBROUTINE volgrid6d_var_delete(this)
1749TYPE(volgrid6d_var),INTENT(INOUT) :: this
1750
1751this%centre = imiss
1752this%category = imiss
1753this%number = imiss
1754this%discipline = imiss
1755this%description = cmiss
1756this%unit = cmiss
1757
1758END SUBROUTINE volgrid6d_var_delete
1759
1760
1761ELEMENTAL FUNCTION volgrid6d_var_c_e(this) RESULT(c_e)
1762TYPE(volgrid6d_var),INTENT(IN) :: this
1763LOGICAL :: c_e
1764c_e = this /= volgrid6d_var_miss
1765END FUNCTION volgrid6d_var_c_e
1766
1767
1768ELEMENTAL FUNCTION volgrid6d_var_eq(this, that) RESULT(res)
1769TYPE(volgrid6d_var),INTENT(IN) :: this, that
1770LOGICAL :: res
1771
1772IF (this%discipline == that%discipline) THEN
1773
1774 IF (this%discipline == 255) THEN ! grib1, WMO tables are all equivalent
1775 res = ((this%category == that%category) .OR. &
1776 (this%category >= 1 .AND. this%category <=3 .AND. &
1777 that%category >= 1 .AND. that%category <=3)) .AND. &
1778 this%number == that%number
1779
1780 IF ((this%category >= 128 .AND. this%category <= 254) .OR. &
1781 (this%number >= 128 .AND. this%number <= 254)) THEN
1782 res = res .AND. this%centre == that%centre ! local definition, centre matters
1783 ENDIF
1784
1785 ELSE ! grib2
1786 res = this%category == that%category .AND. &
1787 this%number == that%number
1788
1789 IF ((this%discipline >= 192 .AND. this%discipline <= 254) .OR. &
1790 (this%category >= 192 .AND. this%category <= 254) .OR. &
1791 (this%number >= 192 .AND. this%number <= 254)) THEN
1792 res = res .AND. this%centre == that%centre ! local definition, centre matters
1793 ENDIF
1794 ENDIF
1795
1796ELSE ! different edition or different discipline
1797 res = .false.
1798ENDIF
1799
1800END FUNCTION volgrid6d_var_eq
1801
1802
1803ELEMENTAL FUNCTION volgrid6d_var_ne(this, that) RESULT(res)
1804TYPE(volgrid6d_var),INTENT(IN) :: this, that
1805LOGICAL :: res
1806
1807res = .NOT.(this == that)
1808
1809END FUNCTION volgrid6d_var_ne
1810
1811
1812#include "array_utilities_inc.F90"
1813
1814
1816SUBROUTINE display_volgrid6d_var(this)
1817TYPE(volgrid6d_var),INTENT(in) :: this
1818
1819print*,"GRIDVAR: ",this%centre,this%discipline,this%category,this%number
1820
1821END SUBROUTINE display_volgrid6d_var
1822
1823
1836SUBROUTINE vargrib2varbufr(vargrib, varbufr, c_func)
1837TYPE(volgrid6d_var),INTENT(in) :: vargrib(:)
1838TYPE(vol7d_var),INTENT(out) :: varbufr(:)
1839TYPE(conv_func),POINTER :: c_func(:)
1840
1841INTEGER :: i, n, stallo
1842
1843n = min(SIZE(varbufr), SIZE(vargrib))
1844ALLOCATE(c_func(n),stat=stallo)
1845IF (stallo /= 0) THEN
1846 call l4f_log(l4f_fatal,"allocating memory")
1847 call raise_fatal_error()
1848ENDIF
1849
1850DO i = 1, n
1851 varbufr(i) = convert(vargrib(i), c_func(i))
1852ENDDO
1853
1854END SUBROUTINE vargrib2varbufr
1855
1856
1867FUNCTION vargrib2varbufr_convert(vargrib, c_func) RESULT(convert)
1868TYPE(volgrid6d_var),INTENT(in) :: vargrib
1869TYPE(conv_func),INTENT(out),OPTIONAL :: c_func
1870TYPE(vol7d_var) :: convert
1871
1872INTEGER :: i
1873
1874IF (.NOT. ALLOCATED(conv_fwd)) CALL vg6d_v7d_var_conv_setup()
1875
1876DO i = 1, SIZE(conv_fwd)
1877 IF (vargrib == conv_fwd(i)%vg6d_var) THEN
1878 convert = conv_fwd(i)%v7d_var
1879 IF (PRESENT(c_func)) c_func = conv_fwd(i)%c_func
1880 RETURN
1881 ENDIF
1882ENDDO
1883! not found
1884convert = vol7d_var_miss
1885IF (PRESENT(c_func)) c_func = conv_func_miss
1886
1887! set hint for backwards conversion
1888convert%gribhint(:) = (/vargrib%centre, vargrib%category, vargrib%number, &
1889 vargrib%discipline/)
1890
1891CALL l4f_log(l4f_warn, 'vargrib2varbufr: variable '// &
1892 trim(to_char(vargrib%centre))//':'//trim(to_char(vargrib%category))//':'// &
1893 trim(to_char(vargrib%number))//':'//trim(to_char(vargrib%discipline))// &
1894 ' not found in table')
1895
1896END FUNCTION vargrib2varbufr_convert
1897
1898
1914SUBROUTINE varbufr2vargrib(varbufr, vargrib, c_func, grid_id_template)
1915TYPE(vol7d_var),INTENT(in) :: varbufr(:)
1916TYPE(volgrid6d_var),INTENT(out) :: vargrib(:)
1917TYPE(conv_func),POINTER :: c_func(:)
1918TYPE(grid_id),INTENT(in),OPTIONAL :: grid_id_template
1919
1920INTEGER :: i, n, stallo
1921
1922n = min(SIZE(varbufr), SIZE(vargrib))
1923ALLOCATE(c_func(n),stat=stallo)
1924IF (stallo /= 0) THEN
1925 CALL l4f_log(l4f_fatal,"allocating memory")
1926 CALL raise_fatal_error()
1927ENDIF
1928
1929DO i = 1, n
1930 vargrib(i) = convert(varbufr(i), c_func(i), grid_id_template)
1931ENDDO
1932
1933END SUBROUTINE varbufr2vargrib
1934
1935
1949FUNCTION varbufr2vargrib_convert(varbufr, c_func, grid_id_template) RESULT(convert)
1950TYPE(vol7d_var),INTENT(in) :: varbufr
1951TYPE(conv_func),INTENT(out),OPTIONAL :: c_func
1952TYPE(grid_id),INTENT(in),OPTIONAL :: grid_id_template
1953TYPE(volgrid6d_var) :: convert
1954
1955INTEGER :: i
1956#ifdef HAVE_LIBGRIBAPI
1957INTEGER :: gaid, editionnumber, category, centre
1958#endif
1959
1960IF (.NOT. ALLOCATED(conv_bwd)) CALL vg6d_v7d_var_conv_setup()
1961
1962#ifdef HAVE_LIBGRIBAPI
1963editionnumber=255; category=255; centre=255
1964#endif
1965IF (PRESENT(grid_id_template)) THEN
1966#ifdef HAVE_LIBGRIBAPI
1967 gaid = grid_id_get_gaid(grid_id_template)
1968 IF (c_e(gaid)) THEN
1969 CALL grib_get(gaid, 'GRIBEditionNumber', editionnumber)
1970 IF (editionnumber == 1) THEN
1971 CALL grib_get(gaid,'gribTablesVersionNo',category)
1972 ENDIF
1973 CALL grib_get(gaid,'centre',centre)
1974 ENDIF
1975#endif
1976ENDIF
1977
1978DO i = 1, SIZE(conv_bwd)
1979 IF (varbufr == conv_bwd(i)%v7d_var) THEN
1980#ifdef HAVE_LIBGRIBAPI
1981 IF (editionnumber /= 255) THEN ! further check required (gaid present)
1982 IF (editionnumber == 1) THEN
1983 IF (conv_bwd(i)%vg6d_var%discipline /= 255) cycle ! wrong edition
1984 ELSE IF (editionnumber == 2) THEN
1985 IF (conv_bwd(i)%vg6d_var%discipline == 255) cycle ! wrong edition
1986 ENDIF
1987 IF (conv_bwd(i)%vg6d_var%centre /= 255 .AND. &
1988 conv_bwd(i)%vg6d_var%centre /= centre) cycle ! wrong centre
1989 ENDIF
1990#endif
1991 convert = conv_bwd(i)%vg6d_var
1992 IF (PRESENT(c_func)) c_func = conv_bwd(i)%c_func
1993 RETURN
1994 ENDIF
1995ENDDO
1996! not found
1997convert = volgrid6d_var_miss
1998IF (PRESENT(c_func)) c_func = conv_func_miss
1999
2000! if hint available use it as a fallback
2001IF (any(varbufr%gribhint /= imiss)) THEN
2002 convert%centre = varbufr%gribhint(1)
2003 convert%category = varbufr%gribhint(2)
2004 convert%number = varbufr%gribhint(3)
2005 convert%discipline = varbufr%gribhint(4)
2006ENDIF
2007
2008CALL l4f_log(l4f_warn, 'varbufr2vargrib: variable '// &
2009 trim(varbufr%btable)//" : "//trim(varbufr%description)//" : "//trim(varbufr%unit)// &
2010 ' not found in table')
2011
2012END FUNCTION varbufr2vargrib_convert
2013
2014
2022SUBROUTINE volgrid6d_var_normalize(this, c_func, grid_id_template)
2023TYPE(volgrid6d_var),INTENT(inout) :: this
2024TYPE(conv_func),INTENT(out) :: c_func
2025TYPE(grid_id),INTENT(in) :: grid_id_template
2026
2027LOGICAL :: eqed, eqcentre
2028INTEGER :: gaid, editionnumber, centre
2029TYPE(volgrid6d_var) :: tmpgrib
2030TYPE(vol7d_var) :: tmpbufr
2031TYPE(conv_func) tmpc_func1, tmpc_func2
2032
2033eqed = .true.
2034eqcentre = .true.
2035c_func = conv_func_miss
2036
2037#ifdef HAVE_LIBGRIBAPI
2038gaid = grid_id_get_gaid(grid_id_template)
2039IF (c_e(gaid)) THEN
2040 CALL grib_get(gaid, 'GRIBEditionNumber', editionnumber)
2041 CALL grib_get(gaid, 'centre', centre)
2042 eqed = editionnumber == 1 .EQV. this%discipline == 255
2043 eqcentre = centre == this%centre
2044ENDIF
2045#endif
2046
2047IF (eqed .AND. eqcentre) RETURN ! nothing to do
2048
2049tmpbufr = convert(this, tmpc_func1)
2050tmpgrib = convert(tmpbufr, tmpc_func2, grid_id_template)
2051
2052IF (tmpgrib /= volgrid6d_var_miss) THEN
2053! conversion back and forth successful, set also conversion function
2054 this = tmpgrib
2055 c_func = tmpc_func1 * tmpc_func2
2056! set to missing in common case to avoid useless computation
2057 IF (c_func == conv_func_identity) c_func = conv_func_miss
2058ELSE IF (.NOT.eqed) THEN
2059! conversion back and forth unsuccessful and grib edition incompatible, set to miss
2060 this = tmpgrib
2061ENDIF
2062
2063END SUBROUTINE volgrid6d_var_normalize
2064
2065
2066! Private subroutine for reading forward and backward conversion tables
2067! todo: better error handling
2068SUBROUTINE vg6d_v7d_var_conv_setup()
2069INTEGER :: un, i, n, stallo
2070
2071! forward, grib to bufr
2072un = open_package_file('vargrib2bufr.csv', filetype_data)
2073n=0
2074DO WHILE(.true.)
2075 READ(un,*,END=100)
2076 n = n + 1
2077ENDDO
2078
2079100 CONTINUE
2080
2081rewind(un)
2082ALLOCATE(conv_fwd(n),stat=stallo)
2083IF (stallo /= 0) THEN
2084 CALL l4f_log(l4f_fatal,"allocating memory")
2085 CALL raise_fatal_error()
2086ENDIF
2087
2088conv_fwd(:) = vg6d_v7d_var_conv_miss
2089CALL import_var_conv(un, conv_fwd)
2090CLOSE(un)
2091
2092! backward, bufr to grib
2093un = open_package_file('vargrib2bufr.csv', filetype_data)
2094! use the same file for now
2095!un = open_package_file('varbufr2grib.csv', filetype_data)
2096n=0
2097DO WHILE(.true.)
2098 READ(un,*,END=300)
2099 n = n + 1
2100ENDDO
2101
2102300 CONTINUE
2103
2104rewind(un)
2105ALLOCATE(conv_bwd(n),stat=stallo)
2106IF (stallo /= 0) THEN
2107 CALL l4f_log(l4f_fatal,"allocating memory")
2108 CALL raise_fatal_error()
2109end if
2110
2111conv_bwd(:) = vg6d_v7d_var_conv_miss
2112CALL import_var_conv(un, conv_bwd)
2113DO i = 1, n
2114 conv_bwd(i)%c_func%a = 1./conv_bwd(i)%c_func%a
2115 conv_bwd(i)%c_func%b = - conv_bwd(i)%c_func%b
2116ENDDO
2117CLOSE(un)
2118
2119CONTAINS
2120
2121SUBROUTINE import_var_conv(un, conv_type)
2122INTEGER, INTENT(in) :: un
2123TYPE(vg6d_v7d_var_conv), INTENT(out) :: conv_type(:)
2124
2125INTEGER :: i
2126TYPE(csv_record) :: csv
2127CHARACTER(len=1024) :: line
2128CHARACTER(len=10) :: btable
2129INTEGER :: centre, category, number, discipline
2130
2131DO i = 1, SIZE(conv_type)
2132 READ(un,'(A)',END=200)line
2133 CALL init(csv, line)
2134 CALL csv_record_getfield(csv, btable)
2135 CALL csv_record_getfield(csv) ! skip fields for description and unit,
2136 CALL csv_record_getfield(csv) ! they correspond to grib information, not bufr Btable
2137 CALL init(conv_type(i)%v7d_var, btable=btable)
2138
2139 CALL csv_record_getfield(csv, centre)
2140 CALL csv_record_getfield(csv, category)
2141 CALL csv_record_getfield(csv, number)
2142 CALL csv_record_getfield(csv, discipline)
2143 CALL init(conv_type(i)%vg6d_var, centre=centre, category=category, &
2144 number=number, discipline=discipline) ! controllare l'ordine
2145
2146 CALL csv_record_getfield(csv, conv_type(i)%c_func%a)
2147 CALL csv_record_getfield(csv, conv_type(i)%c_func%b)
2148 CALL delete(csv)
2149ENDDO
2150
2151200 CONTINUE
2152
2153END SUBROUTINE import_var_conv
2154
2155END SUBROUTINE vg6d_v7d_var_conv_setup
2156
2157
2158ELEMENTAL FUNCTION conv_func_eq(this, that) RESULT(res)
2159TYPE(conv_func),INTENT(IN) :: this, that
2160LOGICAL :: res
2161
2162res = this%a == that%a .AND. this%b == that%b
2163
2164END FUNCTION conv_func_eq
2165
2166
2167ELEMENTAL FUNCTION conv_func_ne(this, that) RESULT(res)
2168TYPE(conv_func),INTENT(IN) :: this, that
2169LOGICAL :: res
2170
2171res = .NOT.(this == that)
2172
2173END FUNCTION conv_func_ne
2174
2175
2176FUNCTION conv_func_mult(this, that) RESULT(mult)
2177TYPE(conv_func),INTENT(in) :: this
2178TYPE(conv_func),INTENT(in) :: that
2179
2180TYPE(conv_func) :: mult
2181
2182IF (this == conv_func_miss .OR. that == conv_func_miss) THEN
2183 mult = conv_func_miss
2184ELSE
2185 mult%a = this%a*that%a
2186 mult%b = this%a*that%b+this%b
2187ENDIF
2188
2189END FUNCTION conv_func_mult
2190
2198ELEMENTAL SUBROUTINE conv_func_compute(this, values)
2199TYPE(conv_func),INTENT(in) :: this
2200REAL,INTENT(inout) :: values
2201
2202IF (this /= conv_func_miss) THEN
2203 IF (c_e(values)) values = values*this%a + this%b
2204ELSE
2205 values=rmiss
2206ENDIF
2207
2208END SUBROUTINE conv_func_compute
2209
2210
2218ELEMENTAL FUNCTION conv_func_convert(this, values) RESULT(convert)
2219TYPE(conv_func),intent(in) :: this
2220REAL,INTENT(in) :: values
2221REAL :: convert
2222
2223convert = values
2224CALL compute(this, convert)
2225
2226END FUNCTION conv_func_convert
2227
2228
2242SUBROUTINE volgrid6d_var_hor_comp_index(this, xind, yind)
2243TYPE(volgrid6d_var),INTENT(in) :: this(:)
2244INTEGER,POINTER :: xind(:), yind(:)
2245
2246TYPE(vol7d_var) :: varbufr(SIZE(this))
2247TYPE(conv_func),POINTER :: c_func(:)
2248INTEGER :: i, nv, counts(SIZE(vol7d_var_horcomp))
2249
2250NULLIFY(xind, yind)
2251counts(:) = 0
2252
2253CALL vargrib2varbufr(this, varbufr, c_func)
2254
2255DO i = 1, SIZE(vol7d_var_horcomp)
2256 counts(i) = count(varbufr(:) == vol7d_var_horcomp(i))
2257ENDDO
2258
2259IF (any(counts(1::2) > 1)) THEN
2260 CALL l4f_log(l4f_warn, '> 1 variable refer to x component of the same field, (un)rotation impossible')
2261 DEALLOCATE(c_func)
2262 RETURN
2263ENDIF
2264IF (any(counts(2::2) > 1)) THEN
2265 CALL l4f_log(l4f_warn, '> 1 variable refer to y component of the same field, (un)rotation impossible')
2266 DEALLOCATE(c_func)
2267 RETURN
2268ENDIF
2269
2270! check that variables are paired and count pairs
2271nv = 0
2272DO i = 1, SIZE(vol7d_var_horcomp), 2
2273 IF (counts(i) == 0 .AND. counts(i+1) > 0) THEN
2274 CALL l4f_log(l4f_warn, 'variable '//trim(vol7d_var_horcomp(i+1)%btable)// &
2275 ' present but the corresponding x-component '// &
2276 trim(vol7d_var_horcomp(i)%btable)//' is missing, (un)rotation impossible')
2277 RETURN
2278 ELSE IF (counts(i+1) == 0 .AND. counts(i) > 0) THEN
2279 CALL l4f_log(l4f_warn, 'variable '//trim(vol7d_var_horcomp(i)%btable)// &
2280 ' present but the corresponding y-component '// &
2281 trim(vol7d_var_horcomp(i+1)%btable)//' is missing, (un)rotation impossible')
2282 RETURN
2283 ENDIF
2284 IF (counts(i) == 1 .AND. counts(i+1) == 1) nv = nv + 1
2285ENDDO
2286
2287! repeat the loop storing indices
2288ALLOCATE(xind(nv), yind(nv))
2289nv = 0
2290DO i = 1, SIZE(vol7d_var_horcomp), 2
2291 IF (counts(i) == 1 .AND. counts(i+1) == 1) THEN
2292 nv = nv + 1
2293 xind(nv) = index(varbufr(:), vol7d_var_horcomp(i))
2294 yind(nv) = index(varbufr(:), vol7d_var_horcomp(i+1))
2295 ENDIF
2296ENDDO
2297DEALLOCATE(c_func)
2298
2299END SUBROUTINE volgrid6d_var_hor_comp_index
2300
2301
2306FUNCTION volgrid6d_var_is_hor_comp(this) RESULT(is_hor_comp)
2307TYPE(volgrid6d_var),INTENT(in) :: this
2308LOGICAL :: is_hor_comp
2309
2310TYPE(vol7d_var) :: varbufr
2311
2312varbufr = convert(this)
2313is_hor_comp = any(varbufr == vol7d_var_horcomp(:))
2314
2315END FUNCTION volgrid6d_var_is_hor_comp
2316
2317! before unstaggering??
2318
2319!IF (.NOT. ALLOCATED(conv_fwd)) CALL vg6d_v7d_var_conv_setup()
2320!
2321!call init(varu,btable="B11003")
2322!call init(varv,btable="B11004")
2323!
2324! test about presence of u and v in standard table
2325!if ( index(conv_fwd(:)%v7d_var,varu) == 0 .or. index(conv_fwd(:)%v7d_var,varv) == 0 )then
2326! call l4f_category_log(this%category,L4F_FATAL, &
2327! "variables B11003 and/or B11004 (wind components) not defined by vg6d_v7d_var_conv_setup")
2328! CALL raise_error()
2329! RETURN
2330!end if
2331!
2332!if (associated(this%var))then
2333! nvar=size(this%var)
2334! allocate(varbufr(nvar),stat=stallo)
2335! if (stallo /=0)then
2336! call l4f_log(L4F_FATAL,"allocating memory")
2337! call raise_fatal_error("allocating memory")
2338! end if
2339!
2340! CALL vargrib2varbufr(this%var, varbufr)
2341!ELSE
2342! CALL l4f_category_log(this%category, L4F_ERROR, &
2343! "trying to destagger an incomplete volgrid6d object")
2344! CALL raise_error()
2345! RETURN
2346!end if
2347!
2348!nvaru=COUNT(varbufr==varu)
2349!nvarv=COUNT(varbufr==varv)
2350!
2351!if (nvaru > 1 )then
2352! call l4f_category_log(this%category,L4F_WARN, &
2353! ">1 variables refer to u wind component, destaggering will not be done ")
2354! DEALLOCATE(varbufr)
2355! RETURN
2356!endif
2357!
2358!if (nvarv > 1 )then
2359! call l4f_category_log(this%category,L4F_WARN, &
2360! ">1 variables refer to v wind component, destaggering will not be done ")
2361! DEALLOCATE(varbufr)
2362! RETURN
2363!endif
2364!
2365!if (nvaru == 0 .and. nvarv == 0) then
2366! call l4f_category_log(this%category,L4F_WARN, &
2367! "no u or v wind component found in volume, nothing to do")
2368! DEALLOCATE(varbufr)
2369! RETURN
2370!endif
2371!
2372!if (COUNT(varbufr/=varu .and. varbufr/=varv) > 0) then
2373! call l4f_category_log(this%category,L4F_WARN, &
2374! "there are variables different from u and v wind component in C grid")
2375!endif
2376
2377
2378END MODULE volgrid6d_var_class
Index method.
Apply the conversion function this to values.
Apply the conversion function this to values.
Destructor for the corresponding object, it assigns it to a missing value.
Display on the screen a brief content of object.
Initialize a volgrid6d_var object with the optional arguments provided.
Gestione degli errori.
Utilities for managing files.
This module defines an abstract interface to different drivers for access to files containing gridded...
Definition of constants to be used for declaring variables of a desired type.
Definition: kinds.F90:251
Definitions of constants and functions for working with missing values.
Classe per la gestione delle variabili osservate da stazioni meteo e affini.
Class for managing physical variables in a grib 1/2 fashion.
Definisce una variabile meteorologica osservata o un suo attributo.
Class defining a real conversion function between units.
Definition of a physical variable in grib coding style.

Generated with Doxygen.