libsim Versione 7.1.11

◆ qcsummaryflagc()

elemental logical function qcsummaryflagc ( character(len=vol7d_cdatalen), intent(in), optional  flag0,
character(len=vol7d_cdatalen), intent(in), optional  flag1,
character(len=vol7d_cdatalen), intent(in), optional  flag2,
character(len=vol7d_cdatalen), intent(in), optional  flag3 
)
private

Check data validity based on multiple confidences.

Compute final decision boolean flag Quality control is complete if one of 3 conditions is verified: a) invalidated data b) gross error check failed c) tot variable less -1 Controllo di validita' del dato basato su test multipli. Per il calcolo della validita' del dato (flag booleano B33007), si prendono in considerazione 3 test; il dato risulta invalidato (flag booleano posto a false) se e solo se uno dei test risulta soddisfatto: a) il dato e' stato invalidato a mano (flag0=B33196=0) b) il dato non ha passato il gross erro check (flag1=B33192=0) c) la variabile tot risulta minore a -1 La variabile tot e' il risultato del confronto tra controllo climatologico (flag1, B33192), controllo temporale (flag2, B33193) e controllo spaziale (flag3, B33194). Ad ognuno di tali controlli e' stato attribuito un punteggio a seconda che ciascuno dei valori relativi ai flag di qualita' risulti inferiore od uguale-maggiore di 10. Nel dettaglio: se B33192 < 10 tot=-1; se B33192>=10 tot=0 se B33193 < 10 tot=-1; se B33193>=10 tot=1 se B33194 < 10 tot=-1; se B33194>=10 tot=1 Ogni dato e' controllato nei 3 flag di qualita' presenti, e viene valutata la somma risultante di tot. Se tot risulta inferiore a -1, qcsummaryflag e' posto a false ed il dato e' invalitato (B33007=0). Se tot risulta maggiore od uguale a -1 qcsummaryflag e' true ed il dato e' valido.

Definizione alla linea 1301 del file modqc.F90.

1302! Copyright (C) 2010 ARPA-SIM <urpsim@smr.arpa.emr.it>
1303! authors:
1304! Davide Cesari <dcesari@arpa.emr.it>
1305! Paolo Patruno <ppatruno@arpa.emr.it>
1306
1307! This program is free software; you can redistribute it and/or
1308! modify it under the terms of the GNU General Public License as
1309! published by the Free Software Foundation; either version 2 of
1310! the License, or (at your option) any later version.
1311
1312! This program is distributed in the hope that it will be useful,
1313! but WITHOUT ANY WARRANTY; without even the implied warranty of
1314! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
1315! GNU General Public License for more details.
1316
1317! You should have received a copy of the GNU General Public License
1318! along with this program. If not, see <http://www.gnu.org/licenses/>.
1319#include "config.h"
1320
1323
1470module modqc
1471use kinds
1474use vol7d_class
1475
1476
1477implicit none
1478
1479
1481type :: qcpartype
1482 integer (kind=int_b):: att
1483 integer (kind=int_b):: gross_error ! special valuer for "*B33192" when gross error check failed
1484 integer (kind=int_b):: invalidated ! special valuer for "*B33196" when manual invalidation happen
1485end type qcpartype
1486
1488type(qcpartype) :: qcpar=qcpartype(10_int_b,0_int_b,1_int_b)
1489
1490integer, parameter :: nqcattrvars=4
1491CHARACTER(len=10),parameter :: qcattrvarsbtables(nqcattrvars)=(/"*B33196","*B33192","*B33193","*B33194"/)
1492
1493type :: qcattrvars
1494 TYPE(vol7d_var) :: vars(nqcattrvars)
1495 CHARACTER(len=10) :: btables(nqcattrvars)
1496end type qcattrvars
1497
1499interface init
1500 module procedure init_qcattrvars
1501end interface
1502
1504interface peeled
1505 module procedure peeledrb, peeleddb, peeledbb, peeledib, peeledcb &
1506 ,peeledri, peeleddi, peeledbi, peeledii, peeledci &
1507 ,peeledrr, peeleddr, peeledbr, peeledir, peeledcr &
1508 ,peeledrd, peeleddd, peeledbd, peeledid, peeledcd &
1509 ,peeledrc, peeleddc, peeledbc, peeledic, peeledcc
1510end interface
1511
1512
1514interface vd
1515 module procedure vdi,vdb,vdr,vdd,vdc
1516end interface
1517
1519interface vdge
1520 module procedure vdgei,vdgeb,vdger,vdged,vdgec
1521end interface
1522
1524interface invalidated
1525 module procedure invalidatedi,invalidatedb,invalidatedr,invalidatedd,invalidatedc
1526end interface
1527
1528private
1529
1530public vd, vdge, init, qcattrvars_new, invalidated, peeled, vol7d_peeling
1531public qcattrvars, nqcattrvars, qcattrvarsbtables
1532public qcpar, qcpartype, qcsummaryflagb ! ,qcsummaryflagi
1533
1534contains
1535
1536
1537! peeled routines
1538#undef VOL7D_POLY_SUBTYPE
1539#undef VOL7D_POLY_SUBTYPES
1540#undef VOL7D_POLY_ISC
1541#define VOL7D_POLY_SUBTYPE REAL
1542#define VOL7D_POLY_SUBTYPES r
1543
1544#undef VOL7D_POLY_TYPE
1545#undef VOL7D_POLY_TYPES
1546#undef VOL7D_POLY_ISC
1547#undef VOL7D_POLY_TYPES_SUBTYPES
1548#define VOL7D_POLY_TYPE REAL
1549#define VOL7D_POLY_TYPES r
1550#define VOL7D_POLY_TYPES_SUBTYPES rr
1551#include "modqc_peeled_include.F90"
1552#include "modqc_peel_util_include.F90"
1553#undef VOL7D_POLY_TYPE
1554#undef VOL7D_POLY_TYPES
1555#undef VOL7D_POLY_TYPES_SUBTYPES
1556#define VOL7D_POLY_TYPE DOUBLE PRECISION
1557#define VOL7D_POLY_TYPES d
1558#define VOL7D_POLY_TYPES_SUBTYPES dr
1559#include "modqc_peeled_include.F90"
1560#include "modqc_peel_util_include.F90"
1561#undef VOL7D_POLY_TYPE
1562#undef VOL7D_POLY_TYPES
1563#undef VOL7D_POLY_TYPES_SUBTYPES
1564#define VOL7D_POLY_TYPE INTEGER
1565#define VOL7D_POLY_TYPES i
1566#define VOL7D_POLY_TYPES_SUBTYPES ir
1567#include "modqc_peeled_include.F90"
1568#include "modqc_peel_util_include.F90"
1569#undef VOL7D_POLY_TYPE
1570#undef VOL7D_POLY_TYPES
1571#undef VOL7D_POLY_TYPES_SUBTYPES
1572#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1573#define VOL7D_POLY_TYPES b
1574#define VOL7D_POLY_TYPES_SUBTYPES br
1575#include "modqc_peeled_include.F90"
1576#include "modqc_peel_util_include.F90"
1577#undef VOL7D_POLY_TYPE
1578#undef VOL7D_POLY_TYPES
1579#undef VOL7D_POLY_TYPES_SUBTYPES
1580#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1581#define VOL7D_POLY_TYPES c
1582#define VOL7D_POLY_ISC = 1
1583#define VOL7D_POLY_TYPES_SUBTYPES cr
1584#include "modqc_peeled_include.F90"
1585#include "modqc_peel_util_include.F90"
1586
1587
1588#undef VOL7D_POLY_SUBTYPE
1589#undef VOL7D_POLY_SUBTYPES
1590#undef VOL7D_POLY_ISC
1591#define VOL7D_POLY_SUBTYPE DOUBLE PRECISION
1592#define VOL7D_POLY_SUBTYPES d
1593
1594#undef VOL7D_POLY_TYPE
1595#undef VOL7D_POLY_TYPES
1596#undef VOL7D_POLY_TYPES_SUBTYPES
1597#define VOL7D_POLY_TYPE REAL
1598#define VOL7D_POLY_TYPES r
1599#define VOL7D_POLY_TYPES_SUBTYPES rd
1600#include "modqc_peeled_include.F90"
1601#undef VOL7D_POLY_TYPE
1602#undef VOL7D_POLY_TYPES
1603#undef VOL7D_POLY_TYPES_SUBTYPES
1604#define VOL7D_POLY_TYPE DOUBLE PRECISION
1605#define VOL7D_POLY_TYPES d
1606#define VOL7D_POLY_TYPES_SUBTYPES dd
1607#include "modqc_peeled_include.F90"
1608#undef VOL7D_POLY_TYPE
1609#undef VOL7D_POLY_TYPES
1610#undef VOL7D_POLY_TYPES_SUBTYPES
1611#define VOL7D_POLY_TYPE INTEGER
1612#define VOL7D_POLY_TYPES i
1613#define VOL7D_POLY_TYPES_SUBTYPES id
1614#include "modqc_peeled_include.F90"
1615#undef VOL7D_POLY_TYPE
1616#undef VOL7D_POLY_TYPES
1617#undef VOL7D_POLY_TYPES_SUBTYPES
1618#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1619#define VOL7D_POLY_TYPES b
1620#define VOL7D_POLY_TYPES_SUBTYPES bd
1621#include "modqc_peeled_include.F90"
1622#undef VOL7D_POLY_TYPE
1623#undef VOL7D_POLY_TYPES
1624#undef VOL7D_POLY_TYPES_SUBTYPES
1625#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1626#define VOL7D_POLY_TYPES c
1627#define VOL7D_POLY_TYPES_SUBTYPES cd
1628#include "modqc_peeled_include.F90"
1629
1630
1631#undef VOL7D_POLY_SUBTYPE
1632#undef VOL7D_POLY_SUBTYPES
1633#undef VOL7D_POLY_ISC
1634#define VOL7D_POLY_SUBTYPE INTEGER
1635#define VOL7D_POLY_SUBTYPES i
1636
1637#undef VOL7D_POLY_TYPE
1638#undef VOL7D_POLY_TYPES
1639#undef VOL7D_POLY_TYPES_SUBTYPES
1640#define VOL7D_POLY_TYPE REAL
1641#define VOL7D_POLY_TYPES r
1642#define VOL7D_POLY_TYPES_SUBTYPES ri
1643#include "modqc_peeled_include.F90"
1644#undef VOL7D_POLY_TYPE
1645#undef VOL7D_POLY_TYPES
1646#undef VOL7D_POLY_TYPES_SUBTYPES
1647#define VOL7D_POLY_TYPE DOUBLE PRECISION
1648#define VOL7D_POLY_TYPES d
1649#define VOL7D_POLY_TYPES_SUBTYPES di
1650#include "modqc_peeled_include.F90"
1651#undef VOL7D_POLY_TYPE
1652#undef VOL7D_POLY_TYPES
1653#undef VOL7D_POLY_TYPES_SUBTYPES
1654#define VOL7D_POLY_TYPE INTEGER
1655#define VOL7D_POLY_TYPES i
1656#define VOL7D_POLY_TYPES_SUBTYPES ii
1657#include "modqc_peeled_include.F90"
1658#undef VOL7D_POLY_TYPE
1659#undef VOL7D_POLY_TYPES
1660#undef VOL7D_POLY_TYPES_SUBTYPES
1661#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1662#define VOL7D_POLY_TYPES b
1663#define VOL7D_POLY_TYPES_SUBTYPES bi
1664#include "modqc_peeled_include.F90"
1665#undef VOL7D_POLY_TYPE
1666#undef VOL7D_POLY_TYPES
1667#undef VOL7D_POLY_TYPES_SUBTYPES
1668#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1669#define VOL7D_POLY_TYPES c
1670#define VOL7D_POLY_ISC = 1
1671#define VOL7D_POLY_TYPES_SUBTYPES ci
1672#include "modqc_peeled_include.F90"
1673
1674
1675#undef VOL7D_POLY_SUBTYPE
1676#undef VOL7D_POLY_SUBTYPES
1677#undef VOL7D_POLY_ISC
1678#define VOL7D_POLY_SUBTYPE INTEGER(kind=int_b)
1679#define VOL7D_POLY_SUBTYPES b
1680
1681#undef VOL7D_POLY_TYPE
1682#undef VOL7D_POLY_TYPES
1683#undef VOL7D_POLY_TYPES_SUBTYPES
1684#define VOL7D_POLY_TYPE REAL
1685#define VOL7D_POLY_TYPES r
1686#define VOL7D_POLY_TYPES_SUBTYPES rb
1687#include "modqc_peeled_include.F90"
1688#undef VOL7D_POLY_TYPE
1689#undef VOL7D_POLY_TYPES
1690#undef VOL7D_POLY_TYPES_SUBTYPES
1691#define VOL7D_POLY_TYPE DOUBLE PRECISION
1692#define VOL7D_POLY_TYPES d
1693#define VOL7D_POLY_TYPES_SUBTYPES db
1694#include "modqc_peeled_include.F90"
1695#undef VOL7D_POLY_TYPE
1696#undef VOL7D_POLY_TYPES
1697#undef VOL7D_POLY_TYPES_SUBTYPES
1698#define VOL7D_POLY_TYPE INTEGER
1699#define VOL7D_POLY_TYPES i
1700#define VOL7D_POLY_TYPES_SUBTYPES ib
1701#include "modqc_peeled_include.F90"
1702#undef VOL7D_POLY_TYPE
1703#undef VOL7D_POLY_TYPES
1704#undef VOL7D_POLY_TYPES_SUBTYPES
1705#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1706#define VOL7D_POLY_TYPES b
1707#define VOL7D_POLY_TYPES_SUBTYPES bb
1708#include "modqc_peeled_include.F90"
1709#undef VOL7D_POLY_TYPE
1710#undef VOL7D_POLY_TYPES
1711#undef VOL7D_POLY_TYPES_SUBTYPES
1712#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1713#define VOL7D_POLY_TYPES c
1714#define VOL7D_POLY_ISC = 1
1715#define VOL7D_POLY_TYPES_SUBTYPES cb
1716#include "modqc_peeled_include.F90"
1717
1718
1719#undef VOL7D_POLY_SUBTYPE
1720#undef VOL7D_POLY_SUBTYPES
1721#undef VOL7D_POLY_ISC
1722#define VOL7D_POLY_SUBTYPE CHARACTER(len=vol7d_cdatalen)
1723#define VOL7D_POLY_SUBTYPES c
1724
1725#undef VOL7D_POLY_TYPE
1726#undef VOL7D_POLY_TYPES
1727#undef VOL7D_POLY_TYPES_SUBTYPES
1728#define VOL7D_POLY_TYPE REAL
1729#define VOL7D_POLY_TYPES r
1730#define VOL7D_POLY_TYPES_SUBTYPES rc
1731#include "modqc_peeled_include.F90"
1732#undef VOL7D_POLY_TYPE
1733#undef VOL7D_POLY_TYPES
1734#undef VOL7D_POLY_TYPES_SUBTYPES
1735#define VOL7D_POLY_TYPE DOUBLE PRECISION
1736#define VOL7D_POLY_TYPES d
1737#define VOL7D_POLY_TYPES_SUBTYPES dc
1738#include "modqc_peeled_include.F90"
1739#undef VOL7D_POLY_TYPE
1740#undef VOL7D_POLY_TYPES
1741#undef VOL7D_POLY_TYPES_SUBTYPES
1742#define VOL7D_POLY_TYPE INTEGER
1743#define VOL7D_POLY_TYPES i
1744#define VOL7D_POLY_TYPES_SUBTYPES ic
1745#include "modqc_peeled_include.F90"
1746#undef VOL7D_POLY_TYPE
1747#undef VOL7D_POLY_TYPES
1748#undef VOL7D_POLY_TYPES_SUBTYPES
1749#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1750#define VOL7D_POLY_TYPES b
1751#define VOL7D_POLY_TYPES_SUBTYPES bc
1752#include "modqc_peeled_include.F90"
1753#undef VOL7D_POLY_TYPE
1754#undef VOL7D_POLY_TYPES
1755#undef VOL7D_POLY_TYPES_SUBTYPES
1756#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1757#define VOL7D_POLY_TYPES c
1758#define VOL7D_POLY_ISC = 1
1759#define VOL7D_POLY_TYPES_SUBTYPES cc
1760#include "modqc_peeled_include.F90"
1761
1762
1763subroutine init_qcattrvars(this)
1764
1765type(qcattrvars),intent(inout) :: this
1766integer :: i
1767
1768this%btables(:) =qcattrvarsbtables
1769do i =1, nqcattrvars
1770 call init(this%vars(i),this%btables(i))
1771end do
1772
1773end subroutine init_qcattrvars
1774
1775
1776type(qcattrvars) function qcattrvars_new()
1777
1778call init(qcattrvars_new)
1779
1780end function qcattrvars_new
1781
1782
1790SUBROUTINE vol7d_peeling(this, data_id, keep_attr, delete_attr, preserve, purgeana)
1791TYPE(vol7d),INTENT(INOUT) :: this
1792integer,INTENT(inout),pointer,OPTIONAL :: data_id(:,:,:,:,:)
1793CHARACTER(len=*),INTENT(in),OPTIONAL :: keep_attr(:)
1794CHARACTER(len=*),INTENT(in),OPTIONAL :: delete_attr(:)
1795logical,intent(in),optional :: preserve
1796logical,intent(in),optional :: purgeana
1797
1798integer :: inddativar,inddatiattrinv,inddatiattrcli,inddatiattrtem,inddatiattrspa,inddativarattr
1799type(qcattrvars) :: attrvars
1800
1801INTEGER(kind=int_b),pointer :: invbb(:,:,:,:,:),clibb(:,:,:,:,:),tembb(:,:,:,:,:),spabb(:,:,:,:,:)
1802INTEGER,pointer :: invbi(:,:,:,:,:),clibi(:,:,:,:,:),tembi(:,:,:,:,:),spabi(:,:,:,:,:)
1803REAL,pointer :: invbr(:,:,:,:,:),clibr(:,:,:,:,:),tembr(:,:,:,:,:),spabr(:,:,:,:,:)
1804DOUBLE PRECISION,pointer :: invbd(:,:,:,:,:),clibd(:,:,:,:,:),tembd(:,:,:,:,:),spabd(:,:,:,:,:)
1805CHARACTER(len=vol7d_cdatalen),pointer :: invbc(:,:,:,:,:),clibc(:,:,:,:,:),tembc(:,:,:,:,:),spabc(:,:,:,:,:)
1806
1807call l4f_log(l4f_info,'starting peeling')
1808
1809call init(attrvars)
1810
1811! generate code per i vari tipi di dati di v7d
1812! tramite un template e il preprocessore
1813
1814
1815#undef VOL7D_POLY_SUBTYPE
1816#undef VOL7D_POLY_SUBTYPES
1817#define VOL7D_POLY_SUBTYPE REAL
1818#define VOL7D_POLY_SUBTYPES r
1819
1820#undef VOL7D_POLY_TYPE
1821#undef VOL7D_POLY_TYPES
1822#define VOL7D_POLY_TYPE REAL
1823#define VOL7D_POLY_TYPES r
1824#include "modqc_peeling_include.F90"
1825#undef VOL7D_POLY_TYPE
1826#undef VOL7D_POLY_TYPES
1827#define VOL7D_POLY_TYPE DOUBLE PRECISION
1828#define VOL7D_POLY_TYPES d
1829#include "modqc_peeling_include.F90"
1830#undef VOL7D_POLY_TYPE
1831#undef VOL7D_POLY_TYPES
1832#define VOL7D_POLY_TYPE INTEGER
1833#define VOL7D_POLY_TYPES i
1834#include "modqc_peeling_include.F90"
1835#undef VOL7D_POLY_TYPE
1836#undef VOL7D_POLY_TYPES
1837#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1838#define VOL7D_POLY_TYPES b
1839#include "modqc_peeling_include.F90"
1840#undef VOL7D_POLY_TYPE
1841#undef VOL7D_POLY_TYPES
1842#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1843#define VOL7D_POLY_TYPES c
1844#include "modqc_peeling_include.F90"
1845
1846
1847#undef VOL7D_POLY_SUBTYPE
1848#undef VOL7D_POLY_SUBTYPES
1849#define VOL7D_POLY_SUBTYPE DOUBLE PRECISION
1850#define VOL7D_POLY_SUBTYPES d
1851
1852#undef VOL7D_POLY_TYPE
1853#undef VOL7D_POLY_TYPES
1854#define VOL7D_POLY_TYPE REAL
1855#define VOL7D_POLY_TYPES r
1856#include "modqc_peeling_include.F90"
1857#undef VOL7D_POLY_TYPE
1858#undef VOL7D_POLY_TYPES
1859#define VOL7D_POLY_TYPE DOUBLE PRECISION
1860#define VOL7D_POLY_TYPES d
1861#include "modqc_peeling_include.F90"
1862#undef VOL7D_POLY_TYPE
1863#undef VOL7D_POLY_TYPES
1864#define VOL7D_POLY_TYPE INTEGER
1865#define VOL7D_POLY_TYPES i
1866#include "modqc_peeling_include.F90"
1867#undef VOL7D_POLY_TYPE
1868#undef VOL7D_POLY_TYPES
1869#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1870#define VOL7D_POLY_TYPES b
1871#include "modqc_peeling_include.F90"
1872#undef VOL7D_POLY_TYPE
1873#undef VOL7D_POLY_TYPES
1874#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1875#define VOL7D_POLY_TYPES c
1876#include "modqc_peeling_include.F90"
1877
1878
1879#undef VOL7D_POLY_SUBTYPE
1880#undef VOL7D_POLY_SUBTYPES
1881#define VOL7D_POLY_SUBTYPE INTEGER
1882#define VOL7D_POLY_SUBTYPES i
1883
1884#undef VOL7D_POLY_TYPE
1885#undef VOL7D_POLY_TYPES
1886#define VOL7D_POLY_TYPE REAL
1887#define VOL7D_POLY_TYPES r
1888#include "modqc_peeling_include.F90"
1889#undef VOL7D_POLY_TYPE
1890#undef VOL7D_POLY_TYPES
1891#define VOL7D_POLY_TYPE DOUBLE PRECISION
1892#define VOL7D_POLY_TYPES d
1893#include "modqc_peeling_include.F90"
1894#undef VOL7D_POLY_TYPE
1895#undef VOL7D_POLY_TYPES
1896#define VOL7D_POLY_TYPE INTEGER
1897#define VOL7D_POLY_TYPES i
1898#include "modqc_peeling_include.F90"
1899#undef VOL7D_POLY_TYPE
1900#undef VOL7D_POLY_TYPES
1901#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1902#define VOL7D_POLY_TYPES b
1903#include "modqc_peeling_include.F90"
1904#undef VOL7D_POLY_TYPE
1905#undef VOL7D_POLY_TYPES
1906#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1907#define VOL7D_POLY_TYPES c
1908#include "modqc_peeling_include.F90"
1909
1910
1911#undef VOL7D_POLY_SUBTYPE
1912#undef VOL7D_POLY_SUBTYPES
1913#define VOL7D_POLY_SUBTYPE INTEGER(kind=int_b)
1914#define VOL7D_POLY_SUBTYPES b
1915
1916#undef VOL7D_POLY_TYPE
1917#undef VOL7D_POLY_TYPES
1918#define VOL7D_POLY_TYPE REAL
1919#define VOL7D_POLY_TYPES r
1920#include "modqc_peeling_include.F90"
1921#undef VOL7D_POLY_TYPE
1922#undef VOL7D_POLY_TYPES
1923#define VOL7D_POLY_TYPE DOUBLE PRECISION
1924#define VOL7D_POLY_TYPES d
1925#include "modqc_peeling_include.F90"
1926#undef VOL7D_POLY_TYPE
1927#undef VOL7D_POLY_TYPES
1928#define VOL7D_POLY_TYPE INTEGER
1929#define VOL7D_POLY_TYPES i
1930#include "modqc_peeling_include.F90"
1931#undef VOL7D_POLY_TYPE
1932#undef VOL7D_POLY_TYPES
1933#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1934#define VOL7D_POLY_TYPES b
1935#include "modqc_peeling_include.F90"
1936#undef VOL7D_POLY_TYPE
1937#undef VOL7D_POLY_TYPES
1938#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1939#define VOL7D_POLY_TYPES c
1940#include "modqc_peeling_include.F90"
1941
1942
1943
1944#undef VOL7D_POLY_SUBTYPE
1945#undef VOL7D_POLY_SUBTYPES
1946#define VOL7D_POLY_SUBTYPE CHARACTER(len=vol7d_cdatalen)
1947#define VOL7D_POLY_SUBTYPES c
1948
1949#undef VOL7D_POLY_TYPE
1950#undef VOL7D_POLY_TYPES
1951#define VOL7D_POLY_TYPE REAL
1952#define VOL7D_POLY_TYPES r
1953#include "modqc_peeling_include.F90"
1954#undef VOL7D_POLY_TYPE
1955#undef VOL7D_POLY_TYPES
1956#define VOL7D_POLY_TYPE DOUBLE PRECISION
1957#define VOL7D_POLY_TYPES d
1958#include "modqc_peeling_include.F90"
1959#undef VOL7D_POLY_TYPE
1960#undef VOL7D_POLY_TYPES
1961#define VOL7D_POLY_TYPE INTEGER
1962#define VOL7D_POLY_TYPES i
1963#include "modqc_peeling_include.F90"
1964#undef VOL7D_POLY_TYPE
1965#undef VOL7D_POLY_TYPES
1966#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1967#define VOL7D_POLY_TYPES b
1968#include "modqc_peeling_include.F90"
1969#undef VOL7D_POLY_TYPE
1970#undef VOL7D_POLY_TYPES
1971#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1972#define VOL7D_POLY_TYPES c
1973#include "modqc_peeling_include.F90"
1974
1975
1976
1977IF (.NOT.PRESENT(keep_attr) .AND. .NOT.PRESENT(delete_attr) .and. .not. optio_log(preserve)) THEN ! destroy all attributes
1978 IF (ASSOCIATED(this%voldatiattrr)) DEALLOCATE(this%voldatiattrr)
1979 IF (ASSOCIATED(this%voldatiattrd)) DEALLOCATE(this%voldatiattrd)
1980 IF (ASSOCIATED(this%voldatiattri)) DEALLOCATE(this%voldatiattri)
1981 IF (ASSOCIATED(this%voldatiattrb)) DEALLOCATE(this%voldatiattrb)
1982 IF (ASSOCIATED(this%voldatiattrc)) DEALLOCATE(this%voldatiattrc)
1983
1984 CALL delete(this%datiattr)
1985 CALL delete(this%dativarattr)
1986END IF
1987
1988IF (PRESENT(keep_attr)) THEN ! set to missing non requested attributes and reform
1989
1990 if (optio_log(preserve)) call l4f_log(l4f_warn,"preserve parameter ignored: keep_attr passed")
1991 CALL keep_var(this%datiattr%r)
1992 CALL keep_var(this%datiattr%d)
1993 CALL keep_var(this%datiattr%i)
1994 CALL keep_var(this%datiattr%b)
1995 CALL keep_var(this%datiattr%c)
1996 CALL qc_reform(this,data_id, miss=.true., purgeana=purgeana)
1997
1998ELSE IF (PRESENT(delete_attr)) THEN ! set to missing requested attributes and reform
1999
2000 if (optio_log(preserve)) call l4f_log(l4f_warn,"preserve parameter ignored: delete_attr passed")
2001 CALL delete_var(this%datiattr%r)
2002 CALL delete_var(this%datiattr%d)
2003 CALL delete_var(this%datiattr%i)
2004 CALL delete_var(this%datiattr%b)
2005 CALL delete_var(this%datiattr%c)
2006 CALL qc_reform(this,data_id, miss=.true., purgeana=purgeana)
2007
2008ELSE IF (PRESENT(purgeana)) THEN
2009
2010 CALL qc_reform(this,data_id, purgeana=purgeana)
2011
2012ENDIF
2013
2014
2015CONTAINS
2016
2017
2019subroutine qc_reform(this,data_id,miss, purgeana)
2020TYPE(vol7d),INTENT(INOUT) :: this
2021integer,INTENT(inout),pointer,OPTIONAL :: data_id(:,:,:,:,:)
2022logical,intent(in),optional :: miss
2023logical,intent(in),optional :: purgeana
2024
2025integer,pointer :: data_idtmp(:,:,:,:,:)
2026logical,allocatable :: llana(:)
2027integer,allocatable :: anaind(:)
2028integer :: i,j,nana
2029
2030if (optio_log(purgeana)) then
2031 allocate(llana(size(this%ana)))
2032 llana =.false.
2033 do i =1,size(this%ana)
2034 if (associated(this%voldatii)) llana(i)= llana(i) .or. any(c_e(this%voldatii(i,:,:,:,:,:)))
2035 if (associated(this%voldatir)) llana(i)= llana(i) .or. any(c_e(this%voldatir(i,:,:,:,:,:)))
2036 if (associated(this%voldatid)) llana(i)= llana(i) .or. any(c_e(this%voldatid(i,:,:,:,:,:)))
2037 if (associated(this%voldatib)) llana(i)= llana(i) .or. any(c_e(this%voldatib(i,:,:,:,:,:)))
2038 if (associated(this%voldatic)) llana(i)= llana(i) .or. any(c_e(this%voldatic(i,:,:,:,:,:)))
2039
2040#ifdef DEBUG
2041 if (.not. llana(i)) call l4f_log(l4f_debug,"remove station"//t2c(i))
2042#endif
2043
2044 end do
2045
2046 nana=count(llana)
2047
2048
2049 allocate(anaind(nana))
2050
2051 j=0
2052 do i=1,size(this%ana)
2053 if (llana(i)) then
2054 j=j+1
2055 anaind(j)=i
2056 end if
2057 end do
2058
2059
2060 if(present(data_id)) then
2061 allocate(data_idtmp(nana,size(data_id,2),size(data_id,3),size(data_id,4),size(data_id,5)))
2062 data_idtmp=data_id(anaind,:,:,:,:)
2063 if (associated(data_id))deallocate(data_id)
2064 data_id=>data_idtmp
2065 end if
2066
2067 call vol7d_reform(this,miss=miss,lana=llana)
2068
2069 deallocate(llana,anaind)
2070
2071else
2072
2073 call vol7d_reform(this,miss=miss)
2074
2075end if
2076
2077end subroutine qc_reform
2078
2079
2080SUBROUTINE keep_var(var)
2081TYPE(vol7d_var),intent(inout),POINTER :: var(:)
2082
2083INTEGER :: i
2084
2085IF (ASSOCIATED(var)) THEN
2086 if (size(var) == 0) then
2087 var%btable = vol7d_var_miss%btable
2088 else
2089 DO i = 1, SIZE(var)
2090 IF (all(var(i)%btable /= keep_attr(:))) THEN ! n.b. ALL((//)) = .TRUE.
2091 var(i)%btable = vol7d_var_miss%btable
2092 ENDIF
2093 ENDDO
2094 end if
2095ENDIF
2096
2097END SUBROUTINE keep_var
2098
2099SUBROUTINE delete_var(var)
2100TYPE(vol7d_var),intent(inout),POINTER :: var(:)
2101
2102INTEGER :: i
2103
2104IF (ASSOCIATED(var)) THEN
2105 if (size(var) == 0) then
2106 var%btable = vol7d_var_miss%btable
2107 else
2108 DO i = 1, SIZE(var)
2109 IF (any(var(i)%btable == delete_attr(:))) THEN ! n.b. ANY((//)) = .FALSE.
2110 var(i) = vol7d_var_miss
2111 ENDIF
2112 ENDDO
2113 end if
2114ENDIF
2115
2116END SUBROUTINE delete_var
2117
2118END SUBROUTINE vol7d_peeling
2119
2120
2121end module modqc
Variables user in Quality Control.
Definition: modqc.F90:392
Test di dato invalidato.
Definition: modqc.F90:417
Remove data under a defined grade of confidence.
Definition: modqc.F90:397
Check data validity based on single confidence.
Definition: modqc.F90:407
Check data validity based on gross error check.
Definition: modqc.F90:412
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.
Utilities and defines for quality control.
Definition: modqc.F90:363
Module for quickly interpreting the OPTIONAL parameters passed to a subprogram.
Classe per la gestione di un volume completo di dati osservati.
Definisce il livello di attendibilità per i dati validi.
Definition: modqc.F90:374

Generated with Doxygen.