libsim Versione 7.2.4
|
◆ qcsummaryflagi()
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 931 del file modqc.F90. 932! Copyright (C) 2010 ARPA-SIM <urpsim@smr.arpa.emr.it>
933! authors:
934! Davide Cesari <dcesari@arpa.emr.it>
935! Paolo Patruno <ppatruno@arpa.emr.it>
936
937! This program is free software; you can redistribute it and/or
938! modify it under the terms of the GNU General Public License as
939! published by the Free Software Foundation; either version 2 of
940! the License, or (at your option) any later version.
941
942! This program is distributed in the hope that it will be useful,
943! but WITHOUT ANY WARRANTY; without even the implied warranty of
944! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
945! GNU General Public License for more details.
946
947! You should have received a copy of the GNU General Public License
948! along with this program. If not, see <http://www.gnu.org/licenses/>.
949#include "config.h"
950
953
1105
1106
1107implicit none
1108
1109
1112 integer (kind=int_b):: att
1113 integer (kind=int_b):: gross_error ! special valuer for "*B33192" when gross error check failed
1114 integer (kind=int_b):: invalidated ! special valuer for "*B33196" when manual invalidation happen
1116
1119
1120integer, parameter :: nqcattrvars=4
1121CHARACTER(len=10),parameter :: qcattrvarsbtables(nqcattrvars)=(/"*B33196","*B33192","*B33193","*B33194"/)
1122
1123type :: qcattrvars
1124 TYPE(vol7d_var) :: vars(nqcattrvars)
1125 CHARACTER(len=10) :: btables(nqcattrvars)
1126end type qcattrvars
1127
1130 module procedure init_qcattrvars
1131end interface
1132
1135 module procedure peeledrb, peeleddb, peeledbb, peeledib, peeledcb &
1136 ,peeledri, peeleddi, peeledbi, peeledii, peeledci &
1137 ,peeledrr, peeleddr, peeledbr, peeledir, peeledcr &
1138 ,peeledrd, peeleddd, peeledbd, peeledid, peeledcd &
1139 ,peeledrc, peeleddc, peeledbc, peeledic, peeledcc
1140end interface
1141
1142
1145 module procedure vdi,vdb,vdr,vdd,vdc
1146end interface
1147
1150 module procedure vdgei,vdgeb,vdger,vdged,vdgec
1151end interface
1152
1155 module procedure invalidatedi,invalidatedb,invalidatedr,invalidatedd,invalidatedc
1156end interface
1157
1158private
1159
1161public qcattrvars, nqcattrvars, qcattrvarsbtables
1163
1164contains
1165
1166
1167! peeled routines
1168#undef VOL7D_POLY_SUBTYPE
1169#undef VOL7D_POLY_SUBTYPES
1170#undef VOL7D_POLY_ISC
1171#define VOL7D_POLY_SUBTYPE REAL
1172#define VOL7D_POLY_SUBTYPES r
1173
1174#undef VOL7D_POLY_TYPE
1175#undef VOL7D_POLY_TYPES
1176#undef VOL7D_POLY_ISC
1177#undef VOL7D_POLY_TYPES_SUBTYPES
1178#define VOL7D_POLY_TYPE REAL
1179#define VOL7D_POLY_TYPES r
1180#define VOL7D_POLY_TYPES_SUBTYPES rr
1181#include "modqc_peeled_include.F90"
1182#include "modqc_peel_util_include.F90"
1183#undef VOL7D_POLY_TYPE
1184#undef VOL7D_POLY_TYPES
1185#undef VOL7D_POLY_TYPES_SUBTYPES
1186#define VOL7D_POLY_TYPE DOUBLE PRECISION
1187#define VOL7D_POLY_TYPES d
1188#define VOL7D_POLY_TYPES_SUBTYPES dr
1189#include "modqc_peeled_include.F90"
1190#include "modqc_peel_util_include.F90"
1191#undef VOL7D_POLY_TYPE
1192#undef VOL7D_POLY_TYPES
1193#undef VOL7D_POLY_TYPES_SUBTYPES
1194#define VOL7D_POLY_TYPE INTEGER
1195#define VOL7D_POLY_TYPES i
1196#define VOL7D_POLY_TYPES_SUBTYPES ir
1197#include "modqc_peeled_include.F90"
1198#include "modqc_peel_util_include.F90"
1199#undef VOL7D_POLY_TYPE
1200#undef VOL7D_POLY_TYPES
1201#undef VOL7D_POLY_TYPES_SUBTYPES
1202#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1203#define VOL7D_POLY_TYPES b
1204#define VOL7D_POLY_TYPES_SUBTYPES br
1205#include "modqc_peeled_include.F90"
1206#include "modqc_peel_util_include.F90"
1207#undef VOL7D_POLY_TYPE
1208#undef VOL7D_POLY_TYPES
1209#undef VOL7D_POLY_TYPES_SUBTYPES
1210#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1211#define VOL7D_POLY_TYPES c
1212#define VOL7D_POLY_ISC = 1
1213#define VOL7D_POLY_TYPES_SUBTYPES cr
1214#include "modqc_peeled_include.F90"
1215#include "modqc_peel_util_include.F90"
1216
1217
1218#undef VOL7D_POLY_SUBTYPE
1219#undef VOL7D_POLY_SUBTYPES
1220#undef VOL7D_POLY_ISC
1221#define VOL7D_POLY_SUBTYPE DOUBLE PRECISION
1222#define VOL7D_POLY_SUBTYPES d
1223
1224#undef VOL7D_POLY_TYPE
1225#undef VOL7D_POLY_TYPES
1226#undef VOL7D_POLY_TYPES_SUBTYPES
1227#define VOL7D_POLY_TYPE REAL
1228#define VOL7D_POLY_TYPES r
1229#define VOL7D_POLY_TYPES_SUBTYPES rd
1230#include "modqc_peeled_include.F90"
1231#undef VOL7D_POLY_TYPE
1232#undef VOL7D_POLY_TYPES
1233#undef VOL7D_POLY_TYPES_SUBTYPES
1234#define VOL7D_POLY_TYPE DOUBLE PRECISION
1235#define VOL7D_POLY_TYPES d
1236#define VOL7D_POLY_TYPES_SUBTYPES dd
1237#include "modqc_peeled_include.F90"
1238#undef VOL7D_POLY_TYPE
1239#undef VOL7D_POLY_TYPES
1240#undef VOL7D_POLY_TYPES_SUBTYPES
1241#define VOL7D_POLY_TYPE INTEGER
1242#define VOL7D_POLY_TYPES i
1243#define VOL7D_POLY_TYPES_SUBTYPES id
1244#include "modqc_peeled_include.F90"
1245#undef VOL7D_POLY_TYPE
1246#undef VOL7D_POLY_TYPES
1247#undef VOL7D_POLY_TYPES_SUBTYPES
1248#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1249#define VOL7D_POLY_TYPES b
1250#define VOL7D_POLY_TYPES_SUBTYPES bd
1251#include "modqc_peeled_include.F90"
1252#undef VOL7D_POLY_TYPE
1253#undef VOL7D_POLY_TYPES
1254#undef VOL7D_POLY_TYPES_SUBTYPES
1255#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1256#define VOL7D_POLY_TYPES c
1257#define VOL7D_POLY_TYPES_SUBTYPES cd
1258#include "modqc_peeled_include.F90"
1259
1260
1261#undef VOL7D_POLY_SUBTYPE
1262#undef VOL7D_POLY_SUBTYPES
1263#undef VOL7D_POLY_ISC
1264#define VOL7D_POLY_SUBTYPE INTEGER
1265#define VOL7D_POLY_SUBTYPES i
1266
1267#undef VOL7D_POLY_TYPE
1268#undef VOL7D_POLY_TYPES
1269#undef VOL7D_POLY_TYPES_SUBTYPES
1270#define VOL7D_POLY_TYPE REAL
1271#define VOL7D_POLY_TYPES r
1272#define VOL7D_POLY_TYPES_SUBTYPES ri
1273#include "modqc_peeled_include.F90"
1274#undef VOL7D_POLY_TYPE
1275#undef VOL7D_POLY_TYPES
1276#undef VOL7D_POLY_TYPES_SUBTYPES
1277#define VOL7D_POLY_TYPE DOUBLE PRECISION
1278#define VOL7D_POLY_TYPES d
1279#define VOL7D_POLY_TYPES_SUBTYPES di
1280#include "modqc_peeled_include.F90"
1281#undef VOL7D_POLY_TYPE
1282#undef VOL7D_POLY_TYPES
1283#undef VOL7D_POLY_TYPES_SUBTYPES
1284#define VOL7D_POLY_TYPE INTEGER
1285#define VOL7D_POLY_TYPES i
1286#define VOL7D_POLY_TYPES_SUBTYPES ii
1287#include "modqc_peeled_include.F90"
1288#undef VOL7D_POLY_TYPE
1289#undef VOL7D_POLY_TYPES
1290#undef VOL7D_POLY_TYPES_SUBTYPES
1291#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1292#define VOL7D_POLY_TYPES b
1293#define VOL7D_POLY_TYPES_SUBTYPES bi
1294#include "modqc_peeled_include.F90"
1295#undef VOL7D_POLY_TYPE
1296#undef VOL7D_POLY_TYPES
1297#undef VOL7D_POLY_TYPES_SUBTYPES
1298#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1299#define VOL7D_POLY_TYPES c
1300#define VOL7D_POLY_ISC = 1
1301#define VOL7D_POLY_TYPES_SUBTYPES ci
1302#include "modqc_peeled_include.F90"
1303
1304
1305#undef VOL7D_POLY_SUBTYPE
1306#undef VOL7D_POLY_SUBTYPES
1307#undef VOL7D_POLY_ISC
1308#define VOL7D_POLY_SUBTYPE INTEGER(kind=int_b)
1309#define VOL7D_POLY_SUBTYPES b
1310
1311#undef VOL7D_POLY_TYPE
1312#undef VOL7D_POLY_TYPES
1313#undef VOL7D_POLY_TYPES_SUBTYPES
1314#define VOL7D_POLY_TYPE REAL
1315#define VOL7D_POLY_TYPES r
1316#define VOL7D_POLY_TYPES_SUBTYPES rb
1317#include "modqc_peeled_include.F90"
1318#undef VOL7D_POLY_TYPE
1319#undef VOL7D_POLY_TYPES
1320#undef VOL7D_POLY_TYPES_SUBTYPES
1321#define VOL7D_POLY_TYPE DOUBLE PRECISION
1322#define VOL7D_POLY_TYPES d
1323#define VOL7D_POLY_TYPES_SUBTYPES db
1324#include "modqc_peeled_include.F90"
1325#undef VOL7D_POLY_TYPE
1326#undef VOL7D_POLY_TYPES
1327#undef VOL7D_POLY_TYPES_SUBTYPES
1328#define VOL7D_POLY_TYPE INTEGER
1329#define VOL7D_POLY_TYPES i
1330#define VOL7D_POLY_TYPES_SUBTYPES ib
1331#include "modqc_peeled_include.F90"
1332#undef VOL7D_POLY_TYPE
1333#undef VOL7D_POLY_TYPES
1334#undef VOL7D_POLY_TYPES_SUBTYPES
1335#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1336#define VOL7D_POLY_TYPES b
1337#define VOL7D_POLY_TYPES_SUBTYPES bb
1338#include "modqc_peeled_include.F90"
1339#undef VOL7D_POLY_TYPE
1340#undef VOL7D_POLY_TYPES
1341#undef VOL7D_POLY_TYPES_SUBTYPES
1342#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1343#define VOL7D_POLY_TYPES c
1344#define VOL7D_POLY_ISC = 1
1345#define VOL7D_POLY_TYPES_SUBTYPES cb
1346#include "modqc_peeled_include.F90"
1347
1348
1349#undef VOL7D_POLY_SUBTYPE
1350#undef VOL7D_POLY_SUBTYPES
1351#undef VOL7D_POLY_ISC
1352#define VOL7D_POLY_SUBTYPE CHARACTER(len=vol7d_cdatalen)
1353#define VOL7D_POLY_SUBTYPES c
1354
1355#undef VOL7D_POLY_TYPE
1356#undef VOL7D_POLY_TYPES
1357#undef VOL7D_POLY_TYPES_SUBTYPES
1358#define VOL7D_POLY_TYPE REAL
1359#define VOL7D_POLY_TYPES r
1360#define VOL7D_POLY_TYPES_SUBTYPES rc
1361#include "modqc_peeled_include.F90"
1362#undef VOL7D_POLY_TYPE
1363#undef VOL7D_POLY_TYPES
1364#undef VOL7D_POLY_TYPES_SUBTYPES
1365#define VOL7D_POLY_TYPE DOUBLE PRECISION
1366#define VOL7D_POLY_TYPES d
1367#define VOL7D_POLY_TYPES_SUBTYPES dc
1368#include "modqc_peeled_include.F90"
1369#undef VOL7D_POLY_TYPE
1370#undef VOL7D_POLY_TYPES
1371#undef VOL7D_POLY_TYPES_SUBTYPES
1372#define VOL7D_POLY_TYPE INTEGER
1373#define VOL7D_POLY_TYPES i
1374#define VOL7D_POLY_TYPES_SUBTYPES ic
1375#include "modqc_peeled_include.F90"
1376#undef VOL7D_POLY_TYPE
1377#undef VOL7D_POLY_TYPES
1378#undef VOL7D_POLY_TYPES_SUBTYPES
1379#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1380#define VOL7D_POLY_TYPES b
1381#define VOL7D_POLY_TYPES_SUBTYPES bc
1382#include "modqc_peeled_include.F90"
1383#undef VOL7D_POLY_TYPE
1384#undef VOL7D_POLY_TYPES
1385#undef VOL7D_POLY_TYPES_SUBTYPES
1386#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1387#define VOL7D_POLY_TYPES c
1388#define VOL7D_POLY_ISC = 1
1389#define VOL7D_POLY_TYPES_SUBTYPES cc
1390#include "modqc_peeled_include.F90"
1391
1392
1393subroutine init_qcattrvars(this)
1394
1395type(qcattrvars),intent(inout) :: this
1396integer :: i
1397
1398this%btables(:) =qcattrvarsbtables
1399do i =1, nqcattrvars
1401end do
1402
1403end subroutine init_qcattrvars
1404
1405
1406type(qcattrvars) function qcattrvars_new()
1407
1409
1410end function qcattrvars_new
1411
1412
1420SUBROUTINE vol7d_peeling(this, data_id, keep_attr, delete_attr, preserve, purgeana)
1421TYPE(vol7d),INTENT(INOUT) :: this
1422integer,INTENT(inout),pointer,OPTIONAL :: data_id(:,:,:,:,:)
1423CHARACTER(len=*),INTENT(in),OPTIONAL :: keep_attr(:)
1424CHARACTER(len=*),INTENT(in),OPTIONAL :: delete_attr(:)
1425logical,intent(in),optional :: preserve
1426logical,intent(in),optional :: purgeana
1427
1428integer :: inddativar,inddatiattrinv,inddatiattrcli,inddatiattrtem,inddatiattrspa,inddativarattr
1429type(qcattrvars) :: attrvars
1430
1431INTEGER(kind=int_b),pointer :: invbb(:,:,:,:,:),clibb(:,:,:,:,:),tembb(:,:,:,:,:),spabb(:,:,:,:,:)
1432INTEGER,pointer :: invbi(:,:,:,:,:),clibi(:,:,:,:,:),tembi(:,:,:,:,:),spabi(:,:,:,:,:)
1433REAL,pointer :: invbr(:,:,:,:,:),clibr(:,:,:,:,:),tembr(:,:,:,:,:),spabr(:,:,:,:,:)
1434DOUBLE PRECISION,pointer :: invbd(:,:,:,:,:),clibd(:,:,:,:,:),tembd(:,:,:,:,:),spabd(:,:,:,:,:)
1435CHARACTER(len=vol7d_cdatalen),pointer :: invbc(:,:,:,:,:),clibc(:,:,:,:,:),tembc(:,:,:,:,:),spabc(:,:,:,:,:)
1436
1437call l4f_log(l4f_info,'starting peeling')
1438
1440
1441! generate code per i vari tipi di dati di v7d
1442! tramite un template e il preprocessore
1443
1444
1445#undef VOL7D_POLY_SUBTYPE
1446#undef VOL7D_POLY_SUBTYPES
1447#define VOL7D_POLY_SUBTYPE REAL
1448#define VOL7D_POLY_SUBTYPES r
1449
1450#undef VOL7D_POLY_TYPE
1451#undef VOL7D_POLY_TYPES
1452#define VOL7D_POLY_TYPE REAL
1453#define VOL7D_POLY_TYPES r
1454#include "modqc_peeling_include.F90"
1455#undef VOL7D_POLY_TYPE
1456#undef VOL7D_POLY_TYPES
1457#define VOL7D_POLY_TYPE DOUBLE PRECISION
1458#define VOL7D_POLY_TYPES d
1459#include "modqc_peeling_include.F90"
1460#undef VOL7D_POLY_TYPE
1461#undef VOL7D_POLY_TYPES
1462#define VOL7D_POLY_TYPE INTEGER
1463#define VOL7D_POLY_TYPES i
1464#include "modqc_peeling_include.F90"
1465#undef VOL7D_POLY_TYPE
1466#undef VOL7D_POLY_TYPES
1467#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1468#define VOL7D_POLY_TYPES b
1469#include "modqc_peeling_include.F90"
1470#undef VOL7D_POLY_TYPE
1471#undef VOL7D_POLY_TYPES
1472#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1473#define VOL7D_POLY_TYPES c
1474#include "modqc_peeling_include.F90"
1475
1476
1477#undef VOL7D_POLY_SUBTYPE
1478#undef VOL7D_POLY_SUBTYPES
1479#define VOL7D_POLY_SUBTYPE DOUBLE PRECISION
1480#define VOL7D_POLY_SUBTYPES d
1481
1482#undef VOL7D_POLY_TYPE
1483#undef VOL7D_POLY_TYPES
1484#define VOL7D_POLY_TYPE REAL
1485#define VOL7D_POLY_TYPES r
1486#include "modqc_peeling_include.F90"
1487#undef VOL7D_POLY_TYPE
1488#undef VOL7D_POLY_TYPES
1489#define VOL7D_POLY_TYPE DOUBLE PRECISION
1490#define VOL7D_POLY_TYPES d
1491#include "modqc_peeling_include.F90"
1492#undef VOL7D_POLY_TYPE
1493#undef VOL7D_POLY_TYPES
1494#define VOL7D_POLY_TYPE INTEGER
1495#define VOL7D_POLY_TYPES i
1496#include "modqc_peeling_include.F90"
1497#undef VOL7D_POLY_TYPE
1498#undef VOL7D_POLY_TYPES
1499#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1500#define VOL7D_POLY_TYPES b
1501#include "modqc_peeling_include.F90"
1502#undef VOL7D_POLY_TYPE
1503#undef VOL7D_POLY_TYPES
1504#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1505#define VOL7D_POLY_TYPES c
1506#include "modqc_peeling_include.F90"
1507
1508
1509#undef VOL7D_POLY_SUBTYPE
1510#undef VOL7D_POLY_SUBTYPES
1511#define VOL7D_POLY_SUBTYPE INTEGER
1512#define VOL7D_POLY_SUBTYPES i
1513
1514#undef VOL7D_POLY_TYPE
1515#undef VOL7D_POLY_TYPES
1516#define VOL7D_POLY_TYPE REAL
1517#define VOL7D_POLY_TYPES r
1518#include "modqc_peeling_include.F90"
1519#undef VOL7D_POLY_TYPE
1520#undef VOL7D_POLY_TYPES
1521#define VOL7D_POLY_TYPE DOUBLE PRECISION
1522#define VOL7D_POLY_TYPES d
1523#include "modqc_peeling_include.F90"
1524#undef VOL7D_POLY_TYPE
1525#undef VOL7D_POLY_TYPES
1526#define VOL7D_POLY_TYPE INTEGER
1527#define VOL7D_POLY_TYPES i
1528#include "modqc_peeling_include.F90"
1529#undef VOL7D_POLY_TYPE
1530#undef VOL7D_POLY_TYPES
1531#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1532#define VOL7D_POLY_TYPES b
1533#include "modqc_peeling_include.F90"
1534#undef VOL7D_POLY_TYPE
1535#undef VOL7D_POLY_TYPES
1536#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1537#define VOL7D_POLY_TYPES c
1538#include "modqc_peeling_include.F90"
1539
1540
1541#undef VOL7D_POLY_SUBTYPE
1542#undef VOL7D_POLY_SUBTYPES
1543#define VOL7D_POLY_SUBTYPE INTEGER(kind=int_b)
1544#define VOL7D_POLY_SUBTYPES b
1545
1546#undef VOL7D_POLY_TYPE
1547#undef VOL7D_POLY_TYPES
1548#define VOL7D_POLY_TYPE REAL
1549#define VOL7D_POLY_TYPES r
1550#include "modqc_peeling_include.F90"
1551#undef VOL7D_POLY_TYPE
1552#undef VOL7D_POLY_TYPES
1553#define VOL7D_POLY_TYPE DOUBLE PRECISION
1554#define VOL7D_POLY_TYPES d
1555#include "modqc_peeling_include.F90"
1556#undef VOL7D_POLY_TYPE
1557#undef VOL7D_POLY_TYPES
1558#define VOL7D_POLY_TYPE INTEGER
1559#define VOL7D_POLY_TYPES i
1560#include "modqc_peeling_include.F90"
1561#undef VOL7D_POLY_TYPE
1562#undef VOL7D_POLY_TYPES
1563#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1564#define VOL7D_POLY_TYPES b
1565#include "modqc_peeling_include.F90"
1566#undef VOL7D_POLY_TYPE
1567#undef VOL7D_POLY_TYPES
1568#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1569#define VOL7D_POLY_TYPES c
1570#include "modqc_peeling_include.F90"
1571
1572
1573
1574#undef VOL7D_POLY_SUBTYPE
1575#undef VOL7D_POLY_SUBTYPES
1576#define VOL7D_POLY_SUBTYPE CHARACTER(len=vol7d_cdatalen)
1577#define VOL7D_POLY_SUBTYPES c
1578
1579#undef VOL7D_POLY_TYPE
1580#undef VOL7D_POLY_TYPES
1581#define VOL7D_POLY_TYPE REAL
1582#define VOL7D_POLY_TYPES r
1583#include "modqc_peeling_include.F90"
1584#undef VOL7D_POLY_TYPE
1585#undef VOL7D_POLY_TYPES
1586#define VOL7D_POLY_TYPE DOUBLE PRECISION
1587#define VOL7D_POLY_TYPES d
1588#include "modqc_peeling_include.F90"
1589#undef VOL7D_POLY_TYPE
1590#undef VOL7D_POLY_TYPES
1591#define VOL7D_POLY_TYPE INTEGER
1592#define VOL7D_POLY_TYPES i
1593#include "modqc_peeling_include.F90"
1594#undef VOL7D_POLY_TYPE
1595#undef VOL7D_POLY_TYPES
1596#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
1597#define VOL7D_POLY_TYPES b
1598#include "modqc_peeling_include.F90"
1599#undef VOL7D_POLY_TYPE
1600#undef VOL7D_POLY_TYPES
1601#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
1602#define VOL7D_POLY_TYPES c
1603#include "modqc_peeling_include.F90"
1604
1605
1606
1607IF (.NOT.PRESENT(keep_attr) .AND. .NOT.PRESENT(delete_attr) .and. .not. optio_log(preserve)) THEN ! destroy all attributes
1608 IF (ASSOCIATED(this%voldatiattrr)) DEALLOCATE(this%voldatiattrr)
1609 IF (ASSOCIATED(this%voldatiattrd)) DEALLOCATE(this%voldatiattrd)
1610 IF (ASSOCIATED(this%voldatiattri)) DEALLOCATE(this%voldatiattri)
1611 IF (ASSOCIATED(this%voldatiattrb)) DEALLOCATE(this%voldatiattrb)
1612 IF (ASSOCIATED(this%voldatiattrc)) DEALLOCATE(this%voldatiattrc)
1613
1614 CALL delete(this%datiattr)
1615 CALL delete(this%dativarattr)
1616END IF
1617
1618IF (PRESENT(keep_attr)) THEN ! set to missing non requested attributes and reform
1619
1620 if (optio_log(preserve)) call l4f_log(l4f_warn,"preserve parameter ignored: keep_attr passed")
1621 CALL keep_var(this%datiattr%r)
1622 CALL keep_var(this%datiattr%d)
1623 CALL keep_var(this%datiattr%i)
1624 CALL keep_var(this%datiattr%b)
1625 CALL keep_var(this%datiattr%c)
1626 CALL qc_reform(this,data_id, miss=.true., purgeana=purgeana)
1627
1628ELSE IF (PRESENT(delete_attr)) THEN ! set to missing requested attributes and reform
1629
1630 if (optio_log(preserve)) call l4f_log(l4f_warn,"preserve parameter ignored: delete_attr passed")
1631 CALL delete_var(this%datiattr%r)
1632 CALL delete_var(this%datiattr%d)
1633 CALL delete_var(this%datiattr%i)
1634 CALL delete_var(this%datiattr%b)
1635 CALL delete_var(this%datiattr%c)
1636 CALL qc_reform(this,data_id, miss=.true., purgeana=purgeana)
1637
1638ELSE IF (PRESENT(purgeana)) THEN
1639
1640 CALL qc_reform(this,data_id, purgeana=purgeana)
1641
1642ENDIF
1643
1644
1645CONTAINS
1646
1647
1649subroutine qc_reform(this,data_id,miss, purgeana)
1650TYPE(vol7d),INTENT(INOUT) :: this
1651integer,INTENT(inout),pointer,OPTIONAL :: data_id(:,:,:,:,:)
1652logical,intent(in),optional :: miss
1653logical,intent(in),optional :: purgeana
1654
1655integer,pointer :: data_idtmp(:,:,:,:,:)
1656logical,allocatable :: llana(:)
1657integer,allocatable :: anaind(:)
1658integer :: i,j,nana
1659
1660if (optio_log(purgeana)) then
1661 allocate(llana(size(this%ana)))
1662 llana =.false.
1663 do i =1,size(this%ana)
1664 if (associated(this%voldatii)) llana(i)= llana(i) .or. any(c_e(this%voldatii(i,:,:,:,:,:)))
1665 if (associated(this%voldatir)) llana(i)= llana(i) .or. any(c_e(this%voldatir(i,:,:,:,:,:)))
1666 if (associated(this%voldatid)) llana(i)= llana(i) .or. any(c_e(this%voldatid(i,:,:,:,:,:)))
1667 if (associated(this%voldatib)) llana(i)= llana(i) .or. any(c_e(this%voldatib(i,:,:,:,:,:)))
1668 if (associated(this%voldatic)) llana(i)= llana(i) .or. any(c_e(this%voldatic(i,:,:,:,:,:)))
1669
1670#ifdef DEBUG
1671 if (.not. llana(i)) call l4f_log(l4f_debug,"remove station"//t2c(i))
1672#endif
1673
1674 end do
1675
1676 nana=count(llana)
1677
1678
1679 allocate(anaind(nana))
1680
1681 j=0
1682 do i=1,size(this%ana)
1683 if (llana(i)) then
1684 j=j+1
1685 anaind(j)=i
1686 end if
1687 end do
1688
1689
1690 if(present(data_id)) then
1691 allocate(data_idtmp(nana,size(data_id,2),size(data_id,3),size(data_id,4),size(data_id,5)))
1692 data_idtmp=data_id(anaind,:,:,:,:)
1693 if (associated(data_id))deallocate(data_id)
1694 data_id=>data_idtmp
1695 end if
1696
1697 call vol7d_reform(this,miss=miss,lana=llana)
1698
1699 deallocate(llana,anaind)
1700
1701else
1702
1703 call vol7d_reform(this,miss=miss)
1704
1705end if
1706
1707end subroutine qc_reform
1708
1709
1710SUBROUTINE keep_var(var)
1711TYPE(vol7d_var),intent(inout),POINTER :: var(:)
1712
1713INTEGER :: i
1714
1715IF (ASSOCIATED(var)) THEN
1716 if (size(var) == 0) then
1717 var%btable = vol7d_var_miss%btable
1718 else
1719 DO i = 1, SIZE(var)
1720 IF (all(var(i)%btable /= keep_attr(:))) THEN ! n.b. ALL((//)) = .TRUE.
1721 var(i)%btable = vol7d_var_miss%btable
1722 ENDIF
1723 ENDDO
1724 end if
1725ENDIF
1726
1727END SUBROUTINE keep_var
1728
1729SUBROUTINE delete_var(var)
1730TYPE(vol7d_var),intent(inout),POINTER :: var(:)
1731
1732INTEGER :: i
1733
1734IF (ASSOCIATED(var)) THEN
1735 if (size(var) == 0) then
1736 var%btable = vol7d_var_miss%btable
1737 else
1738 DO i = 1, SIZE(var)
1739 IF (any(var(i)%btable == delete_attr(:))) THEN ! n.b. ANY((//)) = .FALSE.
1740 var(i) = vol7d_var_miss
1741 ENDIF
1742 ENDDO
1743 end if
1744ENDIF
1745
1746END SUBROUTINE delete_var
1747
1748END SUBROUTINE vol7d_peeling
1749
1750
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 |