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