libsim Versione 7.2.6
|
◆ qcsummaryflagc()
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 1295 del file modqc.F90. 1296! Copyright (C) 2010 ARPA-SIM <urpsim@smr.arpa.emr.it>
1297! authors:
1298! Davide Cesari <dcesari@arpa.emr.it>
1299! Paolo Patruno <ppatruno@arpa.emr.it>
1300
1301! This program is free software; you can redistribute it and/or
1302! modify it under the terms of the GNU General Public License as
1303! published by the Free Software Foundation; either version 2 of
1304! the License, or (at your option) any later version.
1305
1306! This program is distributed in the hope that it will be useful,
1307! but WITHOUT ANY WARRANTY; without even the implied warranty of
1308! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
1309! GNU General Public License for more details.
1310
1311! You should have received a copy of the GNU General Public License
1312! along with this program. If not, see <http://www.gnu.org/licenses/>.
1313#include "config.h"
1314!> \defgroup qc Libsim package, qc library.
1315!! Procedures for data quality control.
1316!! At the moment only climatological quality control is implemented
1317
1318!> Utilities and defines for quality control.
1319!!
1320!! Concise, high-value definitions of Data Quality by expert users,
1321!! analysts, implementers and journalists. This is a great starting point
1322!! to learn about Data Quality.
1323!!
1324!! Data Quality: The Accuracy Dimension
1325!!
1326!! "Data quality is defined as follows: data has quality if it satisfies
1327!! the requirements of its intended use. It lacks quality to the extent
1328!! that it does not satisfy the requirement. In other words, data quality
1329!! depends as much on the intended use as it does on the data itself. To
1330!! satisfy the intended use, the data must be accurate, timely, relevant,
1331!! complete, understood, and trusted." Jack E. Olson
1332!!
1333!! No Data Left Behind: Federal Student Aid - A Case History
1334!!
1335!! "Data quality institutionalizes a set of repeatable processes to
1336!! continuously monitor data and improve data accuracy, completeness,
1337!! timeliness and relevance." Holly Hyland and Lisa Elliott, Federal
1338!! Student Aid
1339!!
1340!! Data Quality: It's a Family Affair
1341!!
1342!! Data Quality definition: "The state of completeness, consistency,
1343!! timeliness and accuracy that makes data appropriate for a specific
1344!! use." Wim Helmer, Dun & Bradstreet
1345!!
1346!! Data Quality and Quality Management - Examples of Quality Evaluation
1347!! Procedures and Quality Management in European National Mapping
1348!! Agencies
1349!!
1350!! "Quality is defined as the totality of characteristics of a product
1351!! that bear on its ability to satisfy stated and implied needs (ISO
1352!! 8402, 1994). In the new ISO/DIS 9000:2000 standard (2000) the
1353!! definition of quality is: 'Ability of a set of inherent
1354!! characteristics of a product, system or process to fulfill
1355!! requirements of customers and other interested parties.' This
1356!! indicates that data quality and quality management are very closely
1357!! related. Data quality is part of the organisation's total quality
1358!! management." Antti Jakobsson
1359!!
1360!! text below from Wikipedia
1361!! http://it.wikipedia.org/wiki/Test_di_verifica_d%27ipotesi
1362!! http://creativecommons.org/licenses/by-sa/3.0/deed.it
1363!! L'ambito statistico
1364!!
1365!! Nel secondo caso la situazione è modificata in quanto interviene un
1366!! elemento nuovo, ovvero il caso. Si supponga di avere una moneta
1367!! recante due facce contrassegnate con testa e croce. Volendo verificare
1368!! l'ipotesi di bilanciamento della moneta si eseguono 20 lanci e si
1369!! contano quelli che danno esito testa. La conseguenza del bilanciamento
1370!! consiste nell'osservare un valore di teste attorno a 10. Tuttavia
1371!! anche in ipotesi di bilanciamento non si può escludere di osservare 20
1372!! teste. D'altronde, l'ipotesi di bilanciamento è logicamente
1373!! compatibile con un numero di teste variante da 0 a 20. In tale
1374!! contesto una qualsiasi decisione in merito all'ipotesi da verificare
1375!! comporta un rischio di errore. Ad esempio rigettare l'ipotesi di
1376!! bilanciamento della moneta avendo osservato 20 teste su 20 lanci
1377!! comporta il rischio di prendere una decisione errata. Nel procedere
1378!! alla verifica dell'ipotesi di bilanciamento della moneta, si ricorre a
1379!! una variabile casuale X. Tale variabile casuale X è una variabile
1380!! aleatoria discreta con distribuzione binomiale B(20; 0,5), dove 20
1381!! indica il numero di lanci e 0,5 la probabilità che si verifichi
1382!! l'evento "testa".
1383!!
1384!! Il risultato sperimentale si deve quindi confrontare con tale
1385!! distribuzione: quanto è distante tale risultato dal valore medio della
1386!! distribuzione B(20; 0,5)? Per rispondere alla domanda si deve
1387!! individuare un valore caratteristico della distribuzione B(20;
1388!! 0,5). Nel nostro caso tale valore caratteristico è il valore medio
1389!! 20/2 = 10. Per valutare la distanza tra il valore sperimentale e
1390!! quello atteso si valuta la probabilità di ottenere un valore
1391!! sperimentale lontano dal valore medio di B(20; 0,5), ossia nel caso
1392!! che dal nostro esperimento risulti X=15 (15 teste dopo 20 lanci), si
1393!! calcola P{|X-10|>=15-10} quindi P{X<=5 oppure X>=15}=0,041.
1394!!
1395!! Quindi, usando una moneta ben bilanciata, la probabilità di ottenere
1396!! un numero di teste X >= 15 (oppure X <= 5) dopo 20 lanci è pari a
1397!! 0,041 ossia al 4,1%. Giudicando bassa tale probabilità si rifiuterà
1398!! l'ipotesi di bilanciamento della moneta in esame, accettando quindi il
1399!! rischio del 4,1% di compiere un errore nel rifiutarla. Di solito, il
1400!! valore della probabilità adottato per rifiutare l'ipotesi nulla è <
1401!! 0,05. Tale valore è detto livello di significatività ed è definibile
1402!! come segue: il livello di significatività sotto l'ipotesi nulla è la
1403!! probabilità di cadere nella zona di rifiuto quando l'ipotesi nulla è
1404!! vera. Tale livello di significatività si indica convenzionalmente con
1405!! α. Il livello di significatività osservato α del test per il quale si
1406!! rifiuterebbe l'ipotesi nulla è detto valore-p (p-value). Riprendendo
1407!! l'esempio sopra riportato il valore-p è pari a 0,041. Adottando
1408!! nell'esempio α = 0,05, si rifiuterà l'ipotesi se
1409!! P{|X-10|>=x}<0,05. Tale condizione si raggiunge appunto se X<6 oppure
1410!! X>14. Tale insieme di valori si definisce convenzionalmente come
1411!! regione di rifiuto. Viceversa l'insieme { 6,7...14} si definisce regione
1412!! di accettazione. In questo modo si è costruita una regola di
1413!! comportamento per verificare l'ipotesi di bilanciamento della
1414!! moneta. Tale regola definisce il test statistico.
1415!!
1416!! In termini tecnici l'ipotesi da verificare si chiama ipotesi nulla e
1417!! si indica con H0, mentre l'ipotesi alternativa con H1. Nel caso della
1418!! moneta, se p è la probabilità di ottenere testa in un lancio la
1419!! verifica di ipotesi si traduce nel seguente sistema:
1420!!
1421!! H_0: p = \frac{1}{2}
1422!! H_1: p \ne \frac{1}{2}
1423!!
1424!! Come già osservato, il modo di condurre un test statistico comporta un
1425!! rischio di errore. Nella pratica statistica si individuano due tipi di
1426!! errori:
1427!!
1428!! 1. rifiutare H0 quando è vera, errore di primo tipo (α) (o errore di prima specie);
1429!! 2. accettare H0 quando è falsa, errore di secondo tipo (β) (o errore di seconda specie).
1430!!
1431!! Tornando all'esempio della moneta in cui la regione di accettazione è
1432!! data dall'insieme di valori {6..14}, la probabilità di rifiutare H0
1433!! quando è vera è stato calcolato pari a 0,041.Tale probabilità
1434!! rappresenta il rischio di incorrere in un errore di primo tipo e si
1435!! indica con α. Per valutare la probabilità di un errore di secondo tipo
1436!! è necessario specificare un valore di p in caso di verità di H1. Si
1437!! supponga che p=0,80, in tal caso la distribuzione di X è una
1438!! B(20;0,80)
1439!!
1440!! Con tale distribuzione di probabilità, l'errore di tipo 2 si calcola
1441!! sommando le probabilità relative ai valori di X della zona di
1442!! accettazione. Si trova quindi che la probabilità cercata è pari a
1443!! circa 0,20. Tale probabilità quantifica il rischio di incorrere
1444!! nell'errore di tipo 2. e si indica convenzionalmente con β. La
1445!! quantità 1-β si chiama potenza del test ed esprime quindi la capacità
1446!! di un test statistico riconoscere la falsità di H0 quando questa è
1447!! effettivamente falsa. La potenza del test trova applicazione nella
1448!! pratica statistica in fase di pianificazione di un esperimento.
1449!!
1450!!Scope of quality checks on observation values
1451!!Checks applied to determine the quality of an observation can range from the very simple to the
1452!!very complex. In roughly increasing order of complexity they can include:
1453!! * Syntactic checks (e.g. an air temperature must be a number to at most 1 decimal
1454!! place);
1455!! * Numeric ranges (e.g. the temperature must fall in the range -90 to +70);
1456!! * Climate range checks (i.e. is the datum consistent with climatology?)
1457!! * Intra-record consistency (e.g. the air temperature must not be less than the dew
1458!! point);
1459!! * Time-series consistency (e.g. the difference between two successive temperatures at
1460!! a site must be 'plausible'); and
1461!! * Spatial consistency (e.g. the station-dependent limits of plausible difference between
1462!! the temperatures at a station and its neighbours must not be violated).
1463!!\ingroup qc
1469
1470
1471implicit none
1472
1473
1474!> Definisce il livello di attendibilità per i dati validi
1476 integer (kind=int_b):: att !< confidence for "*B33192" "*B33193" "*B33194"
1477 integer (kind=int_b):: gross_error ! special valuer for "*B33192" when gross error check failed
1478 integer (kind=int_b):: invalidated ! special valuer for "*B33196" when manual invalidation happen
1480
1481!> Default: data with confidence less or equal 10 are rejected
1483
1484integer, parameter :: nqcattrvars=4
1485CHARACTER(len=10),parameter :: qcattrvarsbtables(nqcattrvars)=(/"*B33196","*B33192","*B33193","*B33194"/)
1486
1487type :: qcattrvars
1488 TYPE(vol7d_var) :: vars(nqcattrvars)
1489 CHARACTER(len=10) :: btables(nqcattrvars)
1490end type qcattrvars
1491
1492!> Variables user in Quality Control
1494 module procedure init_qcattrvars
1495end interface
1496
1497!> Remove data under a defined grade of confidence.
1499 module procedure peeledrb, peeleddb, peeledbb, peeledib, peeledcb &
1500 ,peeledri, peeleddi, peeledbi, peeledii, peeledci &
1501 ,peeledrr, peeleddr, peeledbr, peeledir, peeledcr &
1502 ,peeledrd, peeleddd, peeledbd, peeledid, peeledcd &
1503 ,peeledrc, peeleddc, peeledbc, peeledic, peeledcc
1504end interface
1505
1506
1507!> Check data validity based on single confidence
1509 module procedure vdi,vdb,vdr,vdd,vdc
1510end interface
1511
1512!> Check data validity based on gross error check
1514 module procedure vdgei,vdgeb,vdger,vdged,vdgec
1515end interface
1516
1517!> Test di dato invalidato
1519 module procedure invalidatedi,invalidatedb,invalidatedr,invalidatedd,invalidatedc
1520end interface
1521
1522private
1523
1525public qcattrvars, nqcattrvars, qcattrvarsbtables
1527
1528contains
1529
1530
1531! peeled routines
1532#undef VOL7D_POLY_SUBTYPE
1533#undef VOL7D_POLY_SUBTYPES
1534#undef VOL7D_POLY_ISC
1535#define VOL7D_POLY_SUBTYPE REAL
1536#define VOL7D_POLY_SUBTYPES r
1537
1538#undef VOL7D_POLY_TYPE
1539#undef VOL7D_POLY_TYPES
1540#undef VOL7D_POLY_ISC
1541#undef VOL7D_POLY_TYPES_SUBTYPES
1542#define VOL7D_POLY_TYPE REAL
1543#define VOL7D_POLY_TYPES r
1544#define VOL7D_POLY_TYPES_SUBTYPES rr
1545#include "modqc_peeled_include.F90"
1546#include "modqc_peel_util_include.F90"
1547#undef VOL7D_POLY_TYPE
1548#undef VOL7D_POLY_TYPES
1549#undef VOL7D_POLY_TYPES_SUBTYPES
1550#define VOL7D_POLY_TYPE DOUBLE PRECISION
1551#define VOL7D_POLY_TYPES d
1552#define VOL7D_POLY_TYPES_SUBTYPES dr
1553#include "modqc_peeled_include.F90"
1554#include "modqc_peel_util_include.F90"
1555#undef VOL7D_POLY_TYPE
1556#undef VOL7D_POLY_TYPES
1557#undef VOL7D_POLY_TYPES_SUBTYPES
1558#define VOL7D_POLY_TYPE INTEGER
1559#define VOL7D_POLY_TYPES i
1560#define VOL7D_POLY_TYPES_SUBTYPES ir
1561#include "modqc_peeled_include.F90"
1562#include "modqc_peel_util_include.F90"
1563#undef VOL7D_POLY_TYPE
1564#undef VOL7D_POLY_TYPES
1565#undef VOL7D_POLY_TYPES_SUBTYPES
1566#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1567#define VOL7D_POLY_TYPES b
1568#define VOL7D_POLY_TYPES_SUBTYPES br
1569#include "modqc_peeled_include.F90"
1570#include "modqc_peel_util_include.F90"
1571#undef VOL7D_POLY_TYPE
1572#undef VOL7D_POLY_TYPES
1573#undef VOL7D_POLY_TYPES_SUBTYPES
1574#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1575#define VOL7D_POLY_TYPES c
1576#define VOL7D_POLY_ISC = 1
1577#define VOL7D_POLY_TYPES_SUBTYPES cr
1578#include "modqc_peeled_include.F90"
1579#include "modqc_peel_util_include.F90"
1580
1581
1582#undef VOL7D_POLY_SUBTYPE
1583#undef VOL7D_POLY_SUBTYPES
1584#undef VOL7D_POLY_ISC
1585#define VOL7D_POLY_SUBTYPE DOUBLE PRECISION
1586#define VOL7D_POLY_SUBTYPES d
1587
1588#undef VOL7D_POLY_TYPE
1589#undef VOL7D_POLY_TYPES
1590#undef VOL7D_POLY_TYPES_SUBTYPES
1591#define VOL7D_POLY_TYPE REAL
1592#define VOL7D_POLY_TYPES r
1593#define VOL7D_POLY_TYPES_SUBTYPES rd
1594#include "modqc_peeled_include.F90"
1595#undef VOL7D_POLY_TYPE
1596#undef VOL7D_POLY_TYPES
1597#undef VOL7D_POLY_TYPES_SUBTYPES
1598#define VOL7D_POLY_TYPE DOUBLE PRECISION
1599#define VOL7D_POLY_TYPES d
1600#define VOL7D_POLY_TYPES_SUBTYPES dd
1601#include "modqc_peeled_include.F90"
1602#undef VOL7D_POLY_TYPE
1603#undef VOL7D_POLY_TYPES
1604#undef VOL7D_POLY_TYPES_SUBTYPES
1605#define VOL7D_POLY_TYPE INTEGER
1606#define VOL7D_POLY_TYPES i
1607#define VOL7D_POLY_TYPES_SUBTYPES id
1608#include "modqc_peeled_include.F90"
1609#undef VOL7D_POLY_TYPE
1610#undef VOL7D_POLY_TYPES
1611#undef VOL7D_POLY_TYPES_SUBTYPES
1612#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1613#define VOL7D_POLY_TYPES b
1614#define VOL7D_POLY_TYPES_SUBTYPES bd
1615#include "modqc_peeled_include.F90"
1616#undef VOL7D_POLY_TYPE
1617#undef VOL7D_POLY_TYPES
1618#undef VOL7D_POLY_TYPES_SUBTYPES
1619#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1620#define VOL7D_POLY_TYPES c
1621#define VOL7D_POLY_TYPES_SUBTYPES cd
1622#include "modqc_peeled_include.F90"
1623
1624
1625#undef VOL7D_POLY_SUBTYPE
1626#undef VOL7D_POLY_SUBTYPES
1627#undef VOL7D_POLY_ISC
1628#define VOL7D_POLY_SUBTYPE INTEGER
1629#define VOL7D_POLY_SUBTYPES i
1630
1631#undef VOL7D_POLY_TYPE
1632#undef VOL7D_POLY_TYPES
1633#undef VOL7D_POLY_TYPES_SUBTYPES
1634#define VOL7D_POLY_TYPE REAL
1635#define VOL7D_POLY_TYPES r
1636#define VOL7D_POLY_TYPES_SUBTYPES ri
1637#include "modqc_peeled_include.F90"
1638#undef VOL7D_POLY_TYPE
1639#undef VOL7D_POLY_TYPES
1640#undef VOL7D_POLY_TYPES_SUBTYPES
1641#define VOL7D_POLY_TYPE DOUBLE PRECISION
1642#define VOL7D_POLY_TYPES d
1643#define VOL7D_POLY_TYPES_SUBTYPES di
1644#include "modqc_peeled_include.F90"
1645#undef VOL7D_POLY_TYPE
1646#undef VOL7D_POLY_TYPES
1647#undef VOL7D_POLY_TYPES_SUBTYPES
1648#define VOL7D_POLY_TYPE INTEGER
1649#define VOL7D_POLY_TYPES i
1650#define VOL7D_POLY_TYPES_SUBTYPES ii
1651#include "modqc_peeled_include.F90"
1652#undef VOL7D_POLY_TYPE
1653#undef VOL7D_POLY_TYPES
1654#undef VOL7D_POLY_TYPES_SUBTYPES
1655#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1656#define VOL7D_POLY_TYPES b
1657#define VOL7D_POLY_TYPES_SUBTYPES bi
1658#include "modqc_peeled_include.F90"
1659#undef VOL7D_POLY_TYPE
1660#undef VOL7D_POLY_TYPES
1661#undef VOL7D_POLY_TYPES_SUBTYPES
1662#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1663#define VOL7D_POLY_TYPES c
1664#define VOL7D_POLY_ISC = 1
1665#define VOL7D_POLY_TYPES_SUBTYPES ci
1666#include "modqc_peeled_include.F90"
1667
1668
1669#undef VOL7D_POLY_SUBTYPE
1670#undef VOL7D_POLY_SUBTYPES
1671#undef VOL7D_POLY_ISC
1672#define VOL7D_POLY_SUBTYPE INTEGER(kind=int_b)
1673#define VOL7D_POLY_SUBTYPES b
1674
1675#undef VOL7D_POLY_TYPE
1676#undef VOL7D_POLY_TYPES
1677#undef VOL7D_POLY_TYPES_SUBTYPES
1678#define VOL7D_POLY_TYPE REAL
1679#define VOL7D_POLY_TYPES r
1680#define VOL7D_POLY_TYPES_SUBTYPES rb
1681#include "modqc_peeled_include.F90"
1682#undef VOL7D_POLY_TYPE
1683#undef VOL7D_POLY_TYPES
1684#undef VOL7D_POLY_TYPES_SUBTYPES
1685#define VOL7D_POLY_TYPE DOUBLE PRECISION
1686#define VOL7D_POLY_TYPES d
1687#define VOL7D_POLY_TYPES_SUBTYPES db
1688#include "modqc_peeled_include.F90"
1689#undef VOL7D_POLY_TYPE
1690#undef VOL7D_POLY_TYPES
1691#undef VOL7D_POLY_TYPES_SUBTYPES
1692#define VOL7D_POLY_TYPE INTEGER
1693#define VOL7D_POLY_TYPES i
1694#define VOL7D_POLY_TYPES_SUBTYPES ib
1695#include "modqc_peeled_include.F90"
1696#undef VOL7D_POLY_TYPE
1697#undef VOL7D_POLY_TYPES
1698#undef VOL7D_POLY_TYPES_SUBTYPES
1699#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1700#define VOL7D_POLY_TYPES b
1701#define VOL7D_POLY_TYPES_SUBTYPES bb
1702#include "modqc_peeled_include.F90"
1703#undef VOL7D_POLY_TYPE
1704#undef VOL7D_POLY_TYPES
1705#undef VOL7D_POLY_TYPES_SUBTYPES
1706#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1707#define VOL7D_POLY_TYPES c
1708#define VOL7D_POLY_ISC = 1
1709#define VOL7D_POLY_TYPES_SUBTYPES cb
1710#include "modqc_peeled_include.F90"
1711
1712
1713#undef VOL7D_POLY_SUBTYPE
1714#undef VOL7D_POLY_SUBTYPES
1715#undef VOL7D_POLY_ISC
1716#define VOL7D_POLY_SUBTYPE CHARACTER(len=vol7d_cdatalen)
1717#define VOL7D_POLY_SUBTYPES c
1718
1719#undef VOL7D_POLY_TYPE
1720#undef VOL7D_POLY_TYPES
1721#undef VOL7D_POLY_TYPES_SUBTYPES
1722#define VOL7D_POLY_TYPE REAL
1723#define VOL7D_POLY_TYPES r
1724#define VOL7D_POLY_TYPES_SUBTYPES rc
1725#include "modqc_peeled_include.F90"
1726#undef VOL7D_POLY_TYPE
1727#undef VOL7D_POLY_TYPES
1728#undef VOL7D_POLY_TYPES_SUBTYPES
1729#define VOL7D_POLY_TYPE DOUBLE PRECISION
1730#define VOL7D_POLY_TYPES d
1731#define VOL7D_POLY_TYPES_SUBTYPES dc
1732#include "modqc_peeled_include.F90"
1733#undef VOL7D_POLY_TYPE
1734#undef VOL7D_POLY_TYPES
1735#undef VOL7D_POLY_TYPES_SUBTYPES
1736#define VOL7D_POLY_TYPE INTEGER
1737#define VOL7D_POLY_TYPES i
1738#define VOL7D_POLY_TYPES_SUBTYPES ic
1739#include "modqc_peeled_include.F90"
1740#undef VOL7D_POLY_TYPE
1741#undef VOL7D_POLY_TYPES
1742#undef VOL7D_POLY_TYPES_SUBTYPES
1743#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1744#define VOL7D_POLY_TYPES b
1745#define VOL7D_POLY_TYPES_SUBTYPES bc
1746#include "modqc_peeled_include.F90"
1747#undef VOL7D_POLY_TYPE
1748#undef VOL7D_POLY_TYPES
1749#undef VOL7D_POLY_TYPES_SUBTYPES
1750#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1751#define VOL7D_POLY_TYPES c
1752#define VOL7D_POLY_ISC = 1
1753#define VOL7D_POLY_TYPES_SUBTYPES cc
1754#include "modqc_peeled_include.F90"
1755
1756
1757subroutine init_qcattrvars(this)
1758
1759type(qcattrvars),intent(inout) :: this
1760integer :: i
1761
1762this%btables(:) =qcattrvarsbtables
1763do i =1, nqcattrvars
1765end do
1766
1767end subroutine init_qcattrvars
1768
1769
1770type(qcattrvars) function qcattrvars_new()
1771
1773
1774end function qcattrvars_new
1775
1776
1777!> Remove data under the predefined grade of confidence.
1778!! If neither \a keep_attr nor \a delete_attr are passed, all the
1779!! attributes will be deleted after peeling; if \a keep_attr is
1780!! provided, only attributed listed in \a keep_attr will be kept in
1781!! output, (\a delete_attr will be ignored); if \a delete_attr is
1782!! provided, attributed listed in \a delete_attr will be deleted from
1783!! output.
1784SUBROUTINE vol7d_peeling(this, data_id, keep_attr, delete_attr, preserve, purgeana)
1785TYPE(vol7d),INTENT(INOUT) :: this !< object that has to be peeled
1786integer,INTENT(inout),pointer,OPTIONAL :: data_id(:,:,:,:,:) !< data ID to use with dballe DB (for fast write of attributes)
1787CHARACTER(len=*),INTENT(in),OPTIONAL :: keep_attr(:) !< Btable of attributes that should be kept after removing data
1788CHARACTER(len=*),INTENT(in),OPTIONAL :: delete_attr(:) !< Btable of attributes that should be deleted after removing data
1789logical,intent(in),optional :: preserve !< preserve all attributes if true (alternative to keep_attr and delete_attr)
1790logical,intent(in),optional :: purgeana !< if true remove ana with all data missing
1791
1792integer :: inddativar,inddatiattrinv,inddatiattrcli,inddatiattrtem,inddatiattrspa,inddativarattr
1793type(qcattrvars) :: attrvars
1794
1795INTEGER(kind=int_b),pointer :: invbb(:,:,:,:,:),clibb(:,:,:,:,:),tembb(:,:,:,:,:),spabb(:,:,:,:,:)
1796INTEGER,pointer :: invbi(:,:,:,:,:),clibi(:,:,:,:,:),tembi(:,:,:,:,:),spabi(:,:,:,:,:)
1797REAL,pointer :: invbr(:,:,:,:,:),clibr(:,:,:,:,:),tembr(:,:,:,:,:),spabr(:,:,:,:,:)
1798DOUBLE PRECISION,pointer :: invbd(:,:,:,:,:),clibd(:,:,:,:,:),tembd(:,:,:,:,:),spabd(:,:,:,:,:)
1799CHARACTER(len=vol7d_cdatalen),pointer :: invbc(:,:,:,:,:),clibc(:,:,:,:,:),tembc(:,:,:,:,:),spabc(:,:,:,:,:)
1800
1801call l4f_log(l4f_info,'starting peeling')
1802
1804
1805! generate code per i vari tipi di dati di v7d
1806! tramite un template e il preprocessore
1807
1808
1809#undef VOL7D_POLY_SUBTYPE
1810#undef VOL7D_POLY_SUBTYPES
1811#define VOL7D_POLY_SUBTYPE REAL
1812#define VOL7D_POLY_SUBTYPES r
1813
1814#undef VOL7D_POLY_TYPE
1815#undef VOL7D_POLY_TYPES
1816#define VOL7D_POLY_TYPE REAL
1817#define VOL7D_POLY_TYPES r
1818#include "modqc_peeling_include.F90"
1819#undef VOL7D_POLY_TYPE
1820#undef VOL7D_POLY_TYPES
1821#define VOL7D_POLY_TYPE DOUBLE PRECISION
1822#define VOL7D_POLY_TYPES d
1823#include "modqc_peeling_include.F90"
1824#undef VOL7D_POLY_TYPE
1825#undef VOL7D_POLY_TYPES
1826#define VOL7D_POLY_TYPE INTEGER
1827#define VOL7D_POLY_TYPES i
1828#include "modqc_peeling_include.F90"
1829#undef VOL7D_POLY_TYPE
1830#undef VOL7D_POLY_TYPES
1831#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1832#define VOL7D_POLY_TYPES b
1833#include "modqc_peeling_include.F90"
1834#undef VOL7D_POLY_TYPE
1835#undef VOL7D_POLY_TYPES
1836#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1837#define VOL7D_POLY_TYPES c
1838#include "modqc_peeling_include.F90"
1839
1840
1841#undef VOL7D_POLY_SUBTYPE
1842#undef VOL7D_POLY_SUBTYPES
1843#define VOL7D_POLY_SUBTYPE DOUBLE PRECISION
1844#define VOL7D_POLY_SUBTYPES d
1845
1846#undef VOL7D_POLY_TYPE
1847#undef VOL7D_POLY_TYPES
1848#define VOL7D_POLY_TYPE REAL
1849#define VOL7D_POLY_TYPES r
1850#include "modqc_peeling_include.F90"
1851#undef VOL7D_POLY_TYPE
1852#undef VOL7D_POLY_TYPES
1853#define VOL7D_POLY_TYPE DOUBLE PRECISION
1854#define VOL7D_POLY_TYPES d
1855#include "modqc_peeling_include.F90"
1856#undef VOL7D_POLY_TYPE
1857#undef VOL7D_POLY_TYPES
1858#define VOL7D_POLY_TYPE INTEGER
1859#define VOL7D_POLY_TYPES i
1860#include "modqc_peeling_include.F90"
1861#undef VOL7D_POLY_TYPE
1862#undef VOL7D_POLY_TYPES
1863#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1864#define VOL7D_POLY_TYPES b
1865#include "modqc_peeling_include.F90"
1866#undef VOL7D_POLY_TYPE
1867#undef VOL7D_POLY_TYPES
1868#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1869#define VOL7D_POLY_TYPES c
1870#include "modqc_peeling_include.F90"
1871
1872
1873#undef VOL7D_POLY_SUBTYPE
1874#undef VOL7D_POLY_SUBTYPES
1875#define VOL7D_POLY_SUBTYPE INTEGER
1876#define VOL7D_POLY_SUBTYPES i
1877
1878#undef VOL7D_POLY_TYPE
1879#undef VOL7D_POLY_TYPES
1880#define VOL7D_POLY_TYPE REAL
1881#define VOL7D_POLY_TYPES r
1882#include "modqc_peeling_include.F90"
1883#undef VOL7D_POLY_TYPE
1884#undef VOL7D_POLY_TYPES
1885#define VOL7D_POLY_TYPE DOUBLE PRECISION
1886#define VOL7D_POLY_TYPES d
1887#include "modqc_peeling_include.F90"
1888#undef VOL7D_POLY_TYPE
1889#undef VOL7D_POLY_TYPES
1890#define VOL7D_POLY_TYPE INTEGER
1891#define VOL7D_POLY_TYPES i
1892#include "modqc_peeling_include.F90"
1893#undef VOL7D_POLY_TYPE
1894#undef VOL7D_POLY_TYPES
1895#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1896#define VOL7D_POLY_TYPES b
1897#include "modqc_peeling_include.F90"
1898#undef VOL7D_POLY_TYPE
1899#undef VOL7D_POLY_TYPES
1900#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1901#define VOL7D_POLY_TYPES c
1902#include "modqc_peeling_include.F90"
1903
1904
1905#undef VOL7D_POLY_SUBTYPE
1906#undef VOL7D_POLY_SUBTYPES
1907#define VOL7D_POLY_SUBTYPE INTEGER(kind=int_b)
1908#define VOL7D_POLY_SUBTYPES b
1909
1910#undef VOL7D_POLY_TYPE
1911#undef VOL7D_POLY_TYPES
1912#define VOL7D_POLY_TYPE REAL
1913#define VOL7D_POLY_TYPES r
1914#include "modqc_peeling_include.F90"
1915#undef VOL7D_POLY_TYPE
1916#undef VOL7D_POLY_TYPES
1917#define VOL7D_POLY_TYPE DOUBLE PRECISION
1918#define VOL7D_POLY_TYPES d
1919#include "modqc_peeling_include.F90"
1920#undef VOL7D_POLY_TYPE
1921#undef VOL7D_POLY_TYPES
1922#define VOL7D_POLY_TYPE INTEGER
1923#define VOL7D_POLY_TYPES i
1924#include "modqc_peeling_include.F90"
1925#undef VOL7D_POLY_TYPE
1926#undef VOL7D_POLY_TYPES
1927#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1928#define VOL7D_POLY_TYPES b
1929#include "modqc_peeling_include.F90"
1930#undef VOL7D_POLY_TYPE
1931#undef VOL7D_POLY_TYPES
1932#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1933#define VOL7D_POLY_TYPES c
1934#include "modqc_peeling_include.F90"
1935
1936
1937
1938#undef VOL7D_POLY_SUBTYPE
1939#undef VOL7D_POLY_SUBTYPES
1940#define VOL7D_POLY_SUBTYPE CHARACTER(len=vol7d_cdatalen)
1941#define VOL7D_POLY_SUBTYPES c
1942
1943#undef VOL7D_POLY_TYPE
1944#undef VOL7D_POLY_TYPES
1945#define VOL7D_POLY_TYPE REAL
1946#define VOL7D_POLY_TYPES r
1947#include "modqc_peeling_include.F90"
1948#undef VOL7D_POLY_TYPE
1949#undef VOL7D_POLY_TYPES
1950#define VOL7D_POLY_TYPE DOUBLE PRECISION
1951#define VOL7D_POLY_TYPES d
1952#include "modqc_peeling_include.F90"
1953#undef VOL7D_POLY_TYPE
1954#undef VOL7D_POLY_TYPES
1955#define VOL7D_POLY_TYPE INTEGER
1956#define VOL7D_POLY_TYPES i
1957#include "modqc_peeling_include.F90"
1958#undef VOL7D_POLY_TYPE
1959#undef VOL7D_POLY_TYPES
1960#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1961#define VOL7D_POLY_TYPES b
1962#include "modqc_peeling_include.F90"
1963#undef VOL7D_POLY_TYPE
1964#undef VOL7D_POLY_TYPES
1965#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1966#define VOL7D_POLY_TYPES c
1967#include "modqc_peeling_include.F90"
1968
1969
1970
1971IF (.NOT.PRESENT(keep_attr) .AND. .NOT.PRESENT(delete_attr) .and. .not. optio_log(preserve)) THEN ! destroy all attributes
1972 IF (ASSOCIATED(this%voldatiattrr)) DEALLOCATE(this%voldatiattrr)
1973 IF (ASSOCIATED(this%voldatiattrd)) DEALLOCATE(this%voldatiattrd)
1974 IF (ASSOCIATED(this%voldatiattri)) DEALLOCATE(this%voldatiattri)
1975 IF (ASSOCIATED(this%voldatiattrb)) DEALLOCATE(this%voldatiattrb)
1976 IF (ASSOCIATED(this%voldatiattrc)) DEALLOCATE(this%voldatiattrc)
1977
1978 CALL delete(this%datiattr)
1979 CALL delete(this%dativarattr)
1980END IF
1981
1982IF (PRESENT(keep_attr)) THEN ! set to missing non requested attributes and reform
1983
1984 if (optio_log(preserve)) call l4f_log(l4f_warn,"preserve parameter ignored: keep_attr passed")
1985 CALL keep_var(this%datiattr%r)
1986 CALL keep_var(this%datiattr%d)
1987 CALL keep_var(this%datiattr%i)
1988 CALL keep_var(this%datiattr%b)
1989 CALL keep_var(this%datiattr%c)
1990 CALL qc_reform(this,data_id, miss=.true., purgeana=purgeana)
1991
1992ELSE IF (PRESENT(delete_attr)) THEN ! set to missing requested attributes and reform
1993
1994 if (optio_log(preserve)) call l4f_log(l4f_warn,"preserve parameter ignored: delete_attr passed")
1995 CALL delete_var(this%datiattr%r)
1996 CALL delete_var(this%datiattr%d)
1997 CALL delete_var(this%datiattr%i)
1998 CALL delete_var(this%datiattr%b)
1999 CALL delete_var(this%datiattr%c)
2000 CALL qc_reform(this,data_id, miss=.true., purgeana=purgeana)
2001
2002ELSE IF (PRESENT(purgeana)) THEN
2003
2004 CALL qc_reform(this,data_id, purgeana=purgeana)
2005
2006ENDIF
2007
2008
2009CONTAINS
2010
2011
2012!> Like vol7d_reform but manage data_id and have less options
2013subroutine qc_reform(this,data_id,miss, purgeana)
2014TYPE(vol7d),INTENT(INOUT) :: this !< object that has to be reformed
2015integer,INTENT(inout),pointer,OPTIONAL :: data_id(:,:,:,:,:) !< data ID to use with dballe DB (for fast write of attributes)
2016logical,intent(in),optional :: miss !< remove everithing related with missing position in description vector
2017logical,intent(in),optional :: purgeana !< if true remove ana with all data missing
2018
2019integer,pointer :: data_idtmp(:,:,:,:,:)
2020logical,allocatable :: llana(:)
2021integer,allocatable :: anaind(:)
2022integer :: i,j,nana
2023
2024if (optio_log(purgeana)) then
2025 allocate(llana(size(this%ana)))
2026 llana =.false.
2027 do i =1,size(this%ana)
2028 if (associated(this%voldatii)) llana(i)= llana(i) .or. any(c_e(this%voldatii(i,:,:,:,:,:)))
2029 if (associated(this%voldatir)) llana(i)= llana(i) .or. any(c_e(this%voldatir(i,:,:,:,:,:)))
2030 if (associated(this%voldatid)) llana(i)= llana(i) .or. any(c_e(this%voldatid(i,:,:,:,:,:)))
2031 if (associated(this%voldatib)) llana(i)= llana(i) .or. any(c_e(this%voldatib(i,:,:,:,:,:)))
2032 if (associated(this%voldatic)) llana(i)= llana(i) .or. any(c_e(this%voldatic(i,:,:,:,:,:)))
2033
2034#ifdef DEBUG
2035 if (.not. llana(i)) call l4f_log(l4f_debug,"remove station"//t2c(i))
2036#endif
2037
2038 end do
2039
2040 nana=count(llana)
2041
2042
2043 allocate(anaind(nana))
2044
2045 j=0
2046 do i=1,size(this%ana)
2047 if (llana(i)) then
2048 j=j+1
2049 anaind(j)=i
2050 end if
2051 end do
2052
2053
2054 if(present(data_id)) then
2055 allocate(data_idtmp(nana,size(data_id,2),size(data_id,3),size(data_id,4),size(data_id,5)))
2056 data_idtmp=data_id(anaind,:,:,:,:)
2057 if (associated(data_id))deallocate(data_id)
2058 data_id=>data_idtmp
2059 end if
2060
2061 call vol7d_reform(this,miss=miss,lana=llana)
2062
2063 deallocate(llana,anaind)
2064
2065else
2066
2067 call vol7d_reform(this,miss=miss)
2068
2069end if
2070
2071end subroutine qc_reform
2072
2073
2074SUBROUTINE keep_var(var)
2075TYPE(vol7d_var),intent(inout),POINTER :: var(:)
2076
2077INTEGER :: i
2078
2079IF (ASSOCIATED(var)) THEN
2080 if (size(var) == 0) then
2081 var%btable = vol7d_var_miss%btable
2082 else
2083 DO i = 1, SIZE(var)
2084 IF (all(var(i)%btable /= keep_attr(:))) THEN ! n.b. ALL((//)) = .TRUE.
2085 var(i)%btable = vol7d_var_miss%btable
2086 ENDIF
2087 ENDDO
2088 end if
2089ENDIF
2090
2091END SUBROUTINE keep_var
2092
2093SUBROUTINE delete_var(var)
2094TYPE(vol7d_var),intent(inout),POINTER :: var(:)
2095
2096INTEGER :: i
2097
2098IF (ASSOCIATED(var)) THEN
2099 if (size(var) == 0) then
2100 var%btable = vol7d_var_miss%btable
2101 else
2102 DO i = 1, SIZE(var)
2103 IF (any(var(i)%btable == delete_attr(:))) THEN ! n.b. ANY((//)) = .FALSE.
2104 var(i) = vol7d_var_miss
2105 ENDIF
2106 ENDDO
2107 end if
2108ENDIF
2109
2110END SUBROUTINE delete_var
2111
2112END SUBROUTINE vol7d_peeling
2113
2114
Definition of constants to be used for declaring variables of a desired type. Definition kinds.F90:245 Definitions of constants and functions for working with missing values. Definition missing_values.f90:50 Module for quickly interpreting the OPTIONAL parameters passed to a subprogram. Definition optional_values.f90:28 Classe per la gestione di un volume completo di dati osservati. Definition vol7d_class.F90:273 |