libsim Versione 7.2.0

◆ vol7d_peeling()

subroutine, public vol7d_peeling ( type(vol7d), intent(inout)  this,
integer, dimension(:,:,:,:,:), intent(inout), optional, pointer  data_id,
character(len=*), dimension(:), intent(in), optional  keep_attr,
character(len=*), dimension(:), intent(in), optional  delete_attr,
logical, intent(in), optional  preserve,
logical, intent(in), optional  purgeana 
)

Remove data under the predefined grade of confidence.

If neither keep_attr nor delete_attr are passed, all the attributes will be deleted after peeling; if keep_attr is provided, only attributed listed in keep_attr will be kept in output, (delete_attr will be ignored); if delete_attr is provided, attributed listed in delete_attr will be deleted from output.

Parametri
[in,out]thisobject that has to be peeled
[in,out]data_iddata ID to use with dballe DB (for fast write of attributes)
[in]keep_attrBtable of attributes that should be kept after removing data
[in]delete_attrBtable of attributes that should be deleted after removing data
[in]preservepreserve all attributes if true (alternative to keep_attr and delete_attr)
[in]purgeanaif true remove ana with all data missing

Definizione alla linea 2865 del file modqc.F90.

2866! Copyright (C) 2010 ARPA-SIM <urpsim@smr.arpa.emr.it>
2867! authors:
2868! Davide Cesari <dcesari@arpa.emr.it>
2869! Paolo Patruno <ppatruno@arpa.emr.it>
2870
2871! This program is free software; you can redistribute it and/or
2872! modify it under the terms of the GNU General Public License as
2873! published by the Free Software Foundation; either version 2 of
2874! the License, or (at your option) any later version.
2875
2876! This program is distributed in the hope that it will be useful,
2877! but WITHOUT ANY WARRANTY; without even the implied warranty of
2878! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
2879! GNU General Public License for more details.
2880
2881! You should have received a copy of the GNU General Public License
2882! along with this program. If not, see <http://www.gnu.org/licenses/>.
2883#include "config.h"
2884
2887
3034module modqc
3035use kinds
3038use vol7d_class
3039
3040
3041implicit none
3042
3043
3045type :: qcpartype
3046 integer (kind=int_b):: att
3047 integer (kind=int_b):: gross_error ! special valuer for "*B33192" when gross error check failed
3048 integer (kind=int_b):: invalidated ! special valuer for "*B33196" when manual invalidation happen
3049end type qcpartype
3050
3052type(qcpartype) :: qcpar=qcpartype(10_int_b,0_int_b,1_int_b)
3053
3054integer, parameter :: nqcattrvars=4
3055CHARACTER(len=10),parameter :: qcattrvarsbtables(nqcattrvars)=(/"*B33196","*B33192","*B33193","*B33194"/)
3056
3057type :: qcattrvars
3058 TYPE(vol7d_var) :: vars(nqcattrvars)
3059 CHARACTER(len=10) :: btables(nqcattrvars)
3060end type qcattrvars
3061
3063interface init
3064 module procedure init_qcattrvars
3065end interface
3066
3068interface peeled
3069 module procedure peeledrb, peeleddb, peeledbb, peeledib, peeledcb &
3070 ,peeledri, peeleddi, peeledbi, peeledii, peeledci &
3071 ,peeledrr, peeleddr, peeledbr, peeledir, peeledcr &
3072 ,peeledrd, peeleddd, peeledbd, peeledid, peeledcd &
3073 ,peeledrc, peeleddc, peeledbc, peeledic, peeledcc
3074end interface
3075
3076
3078interface vd
3079 module procedure vdi,vdb,vdr,vdd,vdc
3080end interface
3081
3083interface vdge
3084 module procedure vdgei,vdgeb,vdger,vdged,vdgec
3085end interface
3086
3088interface invalidated
3089 module procedure invalidatedi,invalidatedb,invalidatedr,invalidatedd,invalidatedc
3090end interface
3091
3092private
3093
3094public vd, vdge, init, qcattrvars_new, invalidated, peeled, vol7d_peeling
3095public qcattrvars, nqcattrvars, qcattrvarsbtables
3096public qcpar, qcpartype, qcsummaryflagb ! ,qcsummaryflagi
3097
3098contains
3099
3100
3101! peeled routines
3102#undef VOL7D_POLY_SUBTYPE
3103#undef VOL7D_POLY_SUBTYPES
3104#undef VOL7D_POLY_ISC
3105#define VOL7D_POLY_SUBTYPE REAL
3106#define VOL7D_POLY_SUBTYPES r
3107
3108#undef VOL7D_POLY_TYPE
3109#undef VOL7D_POLY_TYPES
3110#undef VOL7D_POLY_ISC
3111#undef VOL7D_POLY_TYPES_SUBTYPES
3112#define VOL7D_POLY_TYPE REAL
3113#define VOL7D_POLY_TYPES r
3114#define VOL7D_POLY_TYPES_SUBTYPES rr
3115#include "modqc_peeled_include.F90"
3116#include "modqc_peel_util_include.F90"
3117#undef VOL7D_POLY_TYPE
3118#undef VOL7D_POLY_TYPES
3119#undef VOL7D_POLY_TYPES_SUBTYPES
3120#define VOL7D_POLY_TYPE DOUBLE PRECISION
3121#define VOL7D_POLY_TYPES d
3122#define VOL7D_POLY_TYPES_SUBTYPES dr
3123#include "modqc_peeled_include.F90"
3124#include "modqc_peel_util_include.F90"
3125#undef VOL7D_POLY_TYPE
3126#undef VOL7D_POLY_TYPES
3127#undef VOL7D_POLY_TYPES_SUBTYPES
3128#define VOL7D_POLY_TYPE INTEGER
3129#define VOL7D_POLY_TYPES i
3130#define VOL7D_POLY_TYPES_SUBTYPES ir
3131#include "modqc_peeled_include.F90"
3132#include "modqc_peel_util_include.F90"
3133#undef VOL7D_POLY_TYPE
3134#undef VOL7D_POLY_TYPES
3135#undef VOL7D_POLY_TYPES_SUBTYPES
3136#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
3137#define VOL7D_POLY_TYPES b
3138#define VOL7D_POLY_TYPES_SUBTYPES br
3139#include "modqc_peeled_include.F90"
3140#include "modqc_peel_util_include.F90"
3141#undef VOL7D_POLY_TYPE
3142#undef VOL7D_POLY_TYPES
3143#undef VOL7D_POLY_TYPES_SUBTYPES
3144#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
3145#define VOL7D_POLY_TYPES c
3146#define VOL7D_POLY_ISC = 1
3147#define VOL7D_POLY_TYPES_SUBTYPES cr
3148#include "modqc_peeled_include.F90"
3149#include "modqc_peel_util_include.F90"
3150
3151
3152#undef VOL7D_POLY_SUBTYPE
3153#undef VOL7D_POLY_SUBTYPES
3154#undef VOL7D_POLY_ISC
3155#define VOL7D_POLY_SUBTYPE DOUBLE PRECISION
3156#define VOL7D_POLY_SUBTYPES d
3157
3158#undef VOL7D_POLY_TYPE
3159#undef VOL7D_POLY_TYPES
3160#undef VOL7D_POLY_TYPES_SUBTYPES
3161#define VOL7D_POLY_TYPE REAL
3162#define VOL7D_POLY_TYPES r
3163#define VOL7D_POLY_TYPES_SUBTYPES rd
3164#include "modqc_peeled_include.F90"
3165#undef VOL7D_POLY_TYPE
3166#undef VOL7D_POLY_TYPES
3167#undef VOL7D_POLY_TYPES_SUBTYPES
3168#define VOL7D_POLY_TYPE DOUBLE PRECISION
3169#define VOL7D_POLY_TYPES d
3170#define VOL7D_POLY_TYPES_SUBTYPES dd
3171#include "modqc_peeled_include.F90"
3172#undef VOL7D_POLY_TYPE
3173#undef VOL7D_POLY_TYPES
3174#undef VOL7D_POLY_TYPES_SUBTYPES
3175#define VOL7D_POLY_TYPE INTEGER
3176#define VOL7D_POLY_TYPES i
3177#define VOL7D_POLY_TYPES_SUBTYPES id
3178#include "modqc_peeled_include.F90"
3179#undef VOL7D_POLY_TYPE
3180#undef VOL7D_POLY_TYPES
3181#undef VOL7D_POLY_TYPES_SUBTYPES
3182#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
3183#define VOL7D_POLY_TYPES b
3184#define VOL7D_POLY_TYPES_SUBTYPES bd
3185#include "modqc_peeled_include.F90"
3186#undef VOL7D_POLY_TYPE
3187#undef VOL7D_POLY_TYPES
3188#undef VOL7D_POLY_TYPES_SUBTYPES
3189#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
3190#define VOL7D_POLY_TYPES c
3191#define VOL7D_POLY_TYPES_SUBTYPES cd
3192#include "modqc_peeled_include.F90"
3193
3194
3195#undef VOL7D_POLY_SUBTYPE
3196#undef VOL7D_POLY_SUBTYPES
3197#undef VOL7D_POLY_ISC
3198#define VOL7D_POLY_SUBTYPE INTEGER
3199#define VOL7D_POLY_SUBTYPES i
3200
3201#undef VOL7D_POLY_TYPE
3202#undef VOL7D_POLY_TYPES
3203#undef VOL7D_POLY_TYPES_SUBTYPES
3204#define VOL7D_POLY_TYPE REAL
3205#define VOL7D_POLY_TYPES r
3206#define VOL7D_POLY_TYPES_SUBTYPES ri
3207#include "modqc_peeled_include.F90"
3208#undef VOL7D_POLY_TYPE
3209#undef VOL7D_POLY_TYPES
3210#undef VOL7D_POLY_TYPES_SUBTYPES
3211#define VOL7D_POLY_TYPE DOUBLE PRECISION
3212#define VOL7D_POLY_TYPES d
3213#define VOL7D_POLY_TYPES_SUBTYPES di
3214#include "modqc_peeled_include.F90"
3215#undef VOL7D_POLY_TYPE
3216#undef VOL7D_POLY_TYPES
3217#undef VOL7D_POLY_TYPES_SUBTYPES
3218#define VOL7D_POLY_TYPE INTEGER
3219#define VOL7D_POLY_TYPES i
3220#define VOL7D_POLY_TYPES_SUBTYPES ii
3221#include "modqc_peeled_include.F90"
3222#undef VOL7D_POLY_TYPE
3223#undef VOL7D_POLY_TYPES
3224#undef VOL7D_POLY_TYPES_SUBTYPES
3225#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
3226#define VOL7D_POLY_TYPES b
3227#define VOL7D_POLY_TYPES_SUBTYPES bi
3228#include "modqc_peeled_include.F90"
3229#undef VOL7D_POLY_TYPE
3230#undef VOL7D_POLY_TYPES
3231#undef VOL7D_POLY_TYPES_SUBTYPES
3232#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
3233#define VOL7D_POLY_TYPES c
3234#define VOL7D_POLY_ISC = 1
3235#define VOL7D_POLY_TYPES_SUBTYPES ci
3236#include "modqc_peeled_include.F90"
3237
3238
3239#undef VOL7D_POLY_SUBTYPE
3240#undef VOL7D_POLY_SUBTYPES
3241#undef VOL7D_POLY_ISC
3242#define VOL7D_POLY_SUBTYPE INTEGER(kind=int_b)
3243#define VOL7D_POLY_SUBTYPES b
3244
3245#undef VOL7D_POLY_TYPE
3246#undef VOL7D_POLY_TYPES
3247#undef VOL7D_POLY_TYPES_SUBTYPES
3248#define VOL7D_POLY_TYPE REAL
3249#define VOL7D_POLY_TYPES r
3250#define VOL7D_POLY_TYPES_SUBTYPES rb
3251#include "modqc_peeled_include.F90"
3252#undef VOL7D_POLY_TYPE
3253#undef VOL7D_POLY_TYPES
3254#undef VOL7D_POLY_TYPES_SUBTYPES
3255#define VOL7D_POLY_TYPE DOUBLE PRECISION
3256#define VOL7D_POLY_TYPES d
3257#define VOL7D_POLY_TYPES_SUBTYPES db
3258#include "modqc_peeled_include.F90"
3259#undef VOL7D_POLY_TYPE
3260#undef VOL7D_POLY_TYPES
3261#undef VOL7D_POLY_TYPES_SUBTYPES
3262#define VOL7D_POLY_TYPE INTEGER
3263#define VOL7D_POLY_TYPES i
3264#define VOL7D_POLY_TYPES_SUBTYPES ib
3265#include "modqc_peeled_include.F90"
3266#undef VOL7D_POLY_TYPE
3267#undef VOL7D_POLY_TYPES
3268#undef VOL7D_POLY_TYPES_SUBTYPES
3269#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
3270#define VOL7D_POLY_TYPES b
3271#define VOL7D_POLY_TYPES_SUBTYPES bb
3272#include "modqc_peeled_include.F90"
3273#undef VOL7D_POLY_TYPE
3274#undef VOL7D_POLY_TYPES
3275#undef VOL7D_POLY_TYPES_SUBTYPES
3276#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
3277#define VOL7D_POLY_TYPES c
3278#define VOL7D_POLY_ISC = 1
3279#define VOL7D_POLY_TYPES_SUBTYPES cb
3280#include "modqc_peeled_include.F90"
3281
3282
3283#undef VOL7D_POLY_SUBTYPE
3284#undef VOL7D_POLY_SUBTYPES
3285#undef VOL7D_POLY_ISC
3286#define VOL7D_POLY_SUBTYPE CHARACTER(len=vol7d_cdatalen)
3287#define VOL7D_POLY_SUBTYPES c
3288
3289#undef VOL7D_POLY_TYPE
3290#undef VOL7D_POLY_TYPES
3291#undef VOL7D_POLY_TYPES_SUBTYPES
3292#define VOL7D_POLY_TYPE REAL
3293#define VOL7D_POLY_TYPES r
3294#define VOL7D_POLY_TYPES_SUBTYPES rc
3295#include "modqc_peeled_include.F90"
3296#undef VOL7D_POLY_TYPE
3297#undef VOL7D_POLY_TYPES
3298#undef VOL7D_POLY_TYPES_SUBTYPES
3299#define VOL7D_POLY_TYPE DOUBLE PRECISION
3300#define VOL7D_POLY_TYPES d
3301#define VOL7D_POLY_TYPES_SUBTYPES dc
3302#include "modqc_peeled_include.F90"
3303#undef VOL7D_POLY_TYPE
3304#undef VOL7D_POLY_TYPES
3305#undef VOL7D_POLY_TYPES_SUBTYPES
3306#define VOL7D_POLY_TYPE INTEGER
3307#define VOL7D_POLY_TYPES i
3308#define VOL7D_POLY_TYPES_SUBTYPES ic
3309#include "modqc_peeled_include.F90"
3310#undef VOL7D_POLY_TYPE
3311#undef VOL7D_POLY_TYPES
3312#undef VOL7D_POLY_TYPES_SUBTYPES
3313#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
3314#define VOL7D_POLY_TYPES b
3315#define VOL7D_POLY_TYPES_SUBTYPES bc
3316#include "modqc_peeled_include.F90"
3317#undef VOL7D_POLY_TYPE
3318#undef VOL7D_POLY_TYPES
3319#undef VOL7D_POLY_TYPES_SUBTYPES
3320#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
3321#define VOL7D_POLY_TYPES c
3322#define VOL7D_POLY_ISC = 1
3323#define VOL7D_POLY_TYPES_SUBTYPES cc
3324#include "modqc_peeled_include.F90"
3325
3326
3327subroutine init_qcattrvars(this)
3328
3329type(qcattrvars),intent(inout) :: this
3330integer :: i
3331
3332this%btables(:) =qcattrvarsbtables
3333do i =1, nqcattrvars
3334 call init(this%vars(i),this%btables(i))
3335end do
3336
3337end subroutine init_qcattrvars
3338
3339
3340type(qcattrvars) function qcattrvars_new()
3341
3342call init(qcattrvars_new)
3343
3344end function qcattrvars_new
3345
3346
3354SUBROUTINE vol7d_peeling(this, data_id, keep_attr, delete_attr, preserve, purgeana)
3355TYPE(vol7d),INTENT(INOUT) :: this
3356integer,INTENT(inout),pointer,OPTIONAL :: data_id(:,:,:,:,:)
3357CHARACTER(len=*),INTENT(in),OPTIONAL :: keep_attr(:)
3358CHARACTER(len=*),INTENT(in),OPTIONAL :: delete_attr(:)
3359logical,intent(in),optional :: preserve
3360logical,intent(in),optional :: purgeana
3361
3362integer :: inddativar,inddatiattrinv,inddatiattrcli,inddatiattrtem,inddatiattrspa,inddativarattr
3363type(qcattrvars) :: attrvars
3364
3365INTEGER(kind=int_b),pointer :: invbb(:,:,:,:,:),clibb(:,:,:,:,:),tembb(:,:,:,:,:),spabb(:,:,:,:,:)
3366INTEGER,pointer :: invbi(:,:,:,:,:),clibi(:,:,:,:,:),tembi(:,:,:,:,:),spabi(:,:,:,:,:)
3367REAL,pointer :: invbr(:,:,:,:,:),clibr(:,:,:,:,:),tembr(:,:,:,:,:),spabr(:,:,:,:,:)
3368DOUBLE PRECISION,pointer :: invbd(:,:,:,:,:),clibd(:,:,:,:,:),tembd(:,:,:,:,:),spabd(:,:,:,:,:)
3369CHARACTER(len=vol7d_cdatalen),pointer :: invbc(:,:,:,:,:),clibc(:,:,:,:,:),tembc(:,:,:,:,:),spabc(:,:,:,:,:)
3370
3371call l4f_log(l4f_info,'starting peeling')
3372
3373call init(attrvars)
3374
3375! generate code per i vari tipi di dati di v7d
3376! tramite un template e il preprocessore
3377
3378
3379#undef VOL7D_POLY_SUBTYPE
3380#undef VOL7D_POLY_SUBTYPES
3381#define VOL7D_POLY_SUBTYPE REAL
3382#define VOL7D_POLY_SUBTYPES r
3383
3384#undef VOL7D_POLY_TYPE
3385#undef VOL7D_POLY_TYPES
3386#define VOL7D_POLY_TYPE REAL
3387#define VOL7D_POLY_TYPES r
3388#include "modqc_peeling_include.F90"
3389#undef VOL7D_POLY_TYPE
3390#undef VOL7D_POLY_TYPES
3391#define VOL7D_POLY_TYPE DOUBLE PRECISION
3392#define VOL7D_POLY_TYPES d
3393#include "modqc_peeling_include.F90"
3394#undef VOL7D_POLY_TYPE
3395#undef VOL7D_POLY_TYPES
3396#define VOL7D_POLY_TYPE INTEGER
3397#define VOL7D_POLY_TYPES i
3398#include "modqc_peeling_include.F90"
3399#undef VOL7D_POLY_TYPE
3400#undef VOL7D_POLY_TYPES
3401#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
3402#define VOL7D_POLY_TYPES b
3403#include "modqc_peeling_include.F90"
3404#undef VOL7D_POLY_TYPE
3405#undef VOL7D_POLY_TYPES
3406#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
3407#define VOL7D_POLY_TYPES c
3408#include "modqc_peeling_include.F90"
3409
3410
3411#undef VOL7D_POLY_SUBTYPE
3412#undef VOL7D_POLY_SUBTYPES
3413#define VOL7D_POLY_SUBTYPE DOUBLE PRECISION
3414#define VOL7D_POLY_SUBTYPES d
3415
3416#undef VOL7D_POLY_TYPE
3417#undef VOL7D_POLY_TYPES
3418#define VOL7D_POLY_TYPE REAL
3419#define VOL7D_POLY_TYPES r
3420#include "modqc_peeling_include.F90"
3421#undef VOL7D_POLY_TYPE
3422#undef VOL7D_POLY_TYPES
3423#define VOL7D_POLY_TYPE DOUBLE PRECISION
3424#define VOL7D_POLY_TYPES d
3425#include "modqc_peeling_include.F90"
3426#undef VOL7D_POLY_TYPE
3427#undef VOL7D_POLY_TYPES
3428#define VOL7D_POLY_TYPE INTEGER
3429#define VOL7D_POLY_TYPES i
3430#include "modqc_peeling_include.F90"
3431#undef VOL7D_POLY_TYPE
3432#undef VOL7D_POLY_TYPES
3433#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
3434#define VOL7D_POLY_TYPES b
3435#include "modqc_peeling_include.F90"
3436#undef VOL7D_POLY_TYPE
3437#undef VOL7D_POLY_TYPES
3438#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
3439#define VOL7D_POLY_TYPES c
3440#include "modqc_peeling_include.F90"
3441
3442
3443#undef VOL7D_POLY_SUBTYPE
3444#undef VOL7D_POLY_SUBTYPES
3445#define VOL7D_POLY_SUBTYPE INTEGER
3446#define VOL7D_POLY_SUBTYPES i
3447
3448#undef VOL7D_POLY_TYPE
3449#undef VOL7D_POLY_TYPES
3450#define VOL7D_POLY_TYPE REAL
3451#define VOL7D_POLY_TYPES r
3452#include "modqc_peeling_include.F90"
3453#undef VOL7D_POLY_TYPE
3454#undef VOL7D_POLY_TYPES
3455#define VOL7D_POLY_TYPE DOUBLE PRECISION
3456#define VOL7D_POLY_TYPES d
3457#include "modqc_peeling_include.F90"
3458#undef VOL7D_POLY_TYPE
3459#undef VOL7D_POLY_TYPES
3460#define VOL7D_POLY_TYPE INTEGER
3461#define VOL7D_POLY_TYPES i
3462#include "modqc_peeling_include.F90"
3463#undef VOL7D_POLY_TYPE
3464#undef VOL7D_POLY_TYPES
3465#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
3466#define VOL7D_POLY_TYPES b
3467#include "modqc_peeling_include.F90"
3468#undef VOL7D_POLY_TYPE
3469#undef VOL7D_POLY_TYPES
3470#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
3471#define VOL7D_POLY_TYPES c
3472#include "modqc_peeling_include.F90"
3473
3474
3475#undef VOL7D_POLY_SUBTYPE
3476#undef VOL7D_POLY_SUBTYPES
3477#define VOL7D_POLY_SUBTYPE INTEGER(kind=int_b)
3478#define VOL7D_POLY_SUBTYPES b
3479
3480#undef VOL7D_POLY_TYPE
3481#undef VOL7D_POLY_TYPES
3482#define VOL7D_POLY_TYPE REAL
3483#define VOL7D_POLY_TYPES r
3484#include "modqc_peeling_include.F90"
3485#undef VOL7D_POLY_TYPE
3486#undef VOL7D_POLY_TYPES
3487#define VOL7D_POLY_TYPE DOUBLE PRECISION
3488#define VOL7D_POLY_TYPES d
3489#include "modqc_peeling_include.F90"
3490#undef VOL7D_POLY_TYPE
3491#undef VOL7D_POLY_TYPES
3492#define VOL7D_POLY_TYPE INTEGER
3493#define VOL7D_POLY_TYPES i
3494#include "modqc_peeling_include.F90"
3495#undef VOL7D_POLY_TYPE
3496#undef VOL7D_POLY_TYPES
3497#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
3498#define VOL7D_POLY_TYPES b
3499#include "modqc_peeling_include.F90"
3500#undef VOL7D_POLY_TYPE
3501#undef VOL7D_POLY_TYPES
3502#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
3503#define VOL7D_POLY_TYPES c
3504#include "modqc_peeling_include.F90"
3505
3506
3507
3508#undef VOL7D_POLY_SUBTYPE
3509#undef VOL7D_POLY_SUBTYPES
3510#define VOL7D_POLY_SUBTYPE CHARACTER(len=vol7d_cdatalen)
3511#define VOL7D_POLY_SUBTYPES c
3512
3513#undef VOL7D_POLY_TYPE
3514#undef VOL7D_POLY_TYPES
3515#define VOL7D_POLY_TYPE REAL
3516#define VOL7D_POLY_TYPES r
3517#include "modqc_peeling_include.F90"
3518#undef VOL7D_POLY_TYPE
3519#undef VOL7D_POLY_TYPES
3520#define VOL7D_POLY_TYPE DOUBLE PRECISION
3521#define VOL7D_POLY_TYPES d
3522#include "modqc_peeling_include.F90"
3523#undef VOL7D_POLY_TYPE
3524#undef VOL7D_POLY_TYPES
3525#define VOL7D_POLY_TYPE INTEGER
3526#define VOL7D_POLY_TYPES i
3527#include "modqc_peeling_include.F90"
3528#undef VOL7D_POLY_TYPE
3529#undef VOL7D_POLY_TYPES
3530#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
3531#define VOL7D_POLY_TYPES b
3532#include "modqc_peeling_include.F90"
3533#undef VOL7D_POLY_TYPE
3534#undef VOL7D_POLY_TYPES
3535#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
3536#define VOL7D_POLY_TYPES c
3537#include "modqc_peeling_include.F90"
3538
3539
3540
3541IF (.NOT.PRESENT(keep_attr) .AND. .NOT.PRESENT(delete_attr) .and. .not. optio_log(preserve)) THEN ! destroy all attributes
3542 IF (ASSOCIATED(this%voldatiattrr)) DEALLOCATE(this%voldatiattrr)
3543 IF (ASSOCIATED(this%voldatiattrd)) DEALLOCATE(this%voldatiattrd)
3544 IF (ASSOCIATED(this%voldatiattri)) DEALLOCATE(this%voldatiattri)
3545 IF (ASSOCIATED(this%voldatiattrb)) DEALLOCATE(this%voldatiattrb)
3546 IF (ASSOCIATED(this%voldatiattrc)) DEALLOCATE(this%voldatiattrc)
3547
3548 CALL delete(this%datiattr)
3549 CALL delete(this%dativarattr)
3550END IF
3551
3552IF (PRESENT(keep_attr)) THEN ! set to missing non requested attributes and reform
3553
3554 if (optio_log(preserve)) call l4f_log(l4f_warn,"preserve parameter ignored: keep_attr passed")
3555 CALL keep_var(this%datiattr%r)
3556 CALL keep_var(this%datiattr%d)
3557 CALL keep_var(this%datiattr%i)
3558 CALL keep_var(this%datiattr%b)
3559 CALL keep_var(this%datiattr%c)
3560 CALL qc_reform(this,data_id, miss=.true., purgeana=purgeana)
3561
3562ELSE IF (PRESENT(delete_attr)) THEN ! set to missing requested attributes and reform
3563
3564 if (optio_log(preserve)) call l4f_log(l4f_warn,"preserve parameter ignored: delete_attr passed")
3565 CALL delete_var(this%datiattr%r)
3566 CALL delete_var(this%datiattr%d)
3567 CALL delete_var(this%datiattr%i)
3568 CALL delete_var(this%datiattr%b)
3569 CALL delete_var(this%datiattr%c)
3570 CALL qc_reform(this,data_id, miss=.true., purgeana=purgeana)
3571
3572ELSE IF (PRESENT(purgeana)) THEN
3573
3574 CALL qc_reform(this,data_id, purgeana=purgeana)
3575
3576ENDIF
3577
3578
3579CONTAINS
3580
3581
3583subroutine qc_reform(this,data_id,miss, purgeana)
3584TYPE(vol7d),INTENT(INOUT) :: this
3585integer,INTENT(inout),pointer,OPTIONAL :: data_id(:,:,:,:,:)
3586logical,intent(in),optional :: miss
3587logical,intent(in),optional :: purgeana
3588
3589integer,pointer :: data_idtmp(:,:,:,:,:)
3590logical,allocatable :: llana(:)
3591integer,allocatable :: anaind(:)
3592integer :: i,j,nana
3593
3594if (optio_log(purgeana)) then
3595 allocate(llana(size(this%ana)))
3596 llana =.false.
3597 do i =1,size(this%ana)
3598 if (associated(this%voldatii)) llana(i)= llana(i) .or. any(c_e(this%voldatii(i,:,:,:,:,:)))
3599 if (associated(this%voldatir)) llana(i)= llana(i) .or. any(c_e(this%voldatir(i,:,:,:,:,:)))
3600 if (associated(this%voldatid)) llana(i)= llana(i) .or. any(c_e(this%voldatid(i,:,:,:,:,:)))
3601 if (associated(this%voldatib)) llana(i)= llana(i) .or. any(c_e(this%voldatib(i,:,:,:,:,:)))
3602 if (associated(this%voldatic)) llana(i)= llana(i) .or. any(c_e(this%voldatic(i,:,:,:,:,:)))
3603
3604#ifdef DEBUG
3605 if (.not. llana(i)) call l4f_log(l4f_debug,"remove station"//t2c(i))
3606#endif
3607
3608 end do
3609
3610 nana=count(llana)
3611
3612
3613 allocate(anaind(nana))
3614
3615 j=0
3616 do i=1,size(this%ana)
3617 if (llana(i)) then
3618 j=j+1
3619 anaind(j)=i
3620 end if
3621 end do
3622
3623
3624 if(present(data_id)) then
3625 allocate(data_idtmp(nana,size(data_id,2),size(data_id,3),size(data_id,4),size(data_id,5)))
3626 data_idtmp=data_id(anaind,:,:,:,:)
3627 if (associated(data_id))deallocate(data_id)
3628 data_id=>data_idtmp
3629 end if
3630
3631 call vol7d_reform(this,miss=miss,lana=llana)
3632
3633 deallocate(llana,anaind)
3634
3635else
3636
3637 call vol7d_reform(this,miss=miss)
3638
3639end if
3640
3641end subroutine qc_reform
3642
3643
3644SUBROUTINE keep_var(var)
3645TYPE(vol7d_var),intent(inout),POINTER :: var(:)
3646
3647INTEGER :: i
3648
3649IF (ASSOCIATED(var)) THEN
3650 if (size(var) == 0) then
3651 var%btable = vol7d_var_miss%btable
3652 else
3653 DO i = 1, SIZE(var)
3654 IF (all(var(i)%btable /= keep_attr(:))) THEN ! n.b. ALL((//)) = .TRUE.
3655 var(i)%btable = vol7d_var_miss%btable
3656 ENDIF
3657 ENDDO
3658 end if
3659ENDIF
3660
3661END SUBROUTINE keep_var
3662
3663SUBROUTINE delete_var(var)
3664TYPE(vol7d_var),intent(inout),POINTER :: var(:)
3665
3666INTEGER :: i
3667
3668IF (ASSOCIATED(var)) THEN
3669 if (size(var) == 0) then
3670 var%btable = vol7d_var_miss%btable
3671 else
3672 DO i = 1, SIZE(var)
3673 IF (any(var(i)%btable == delete_attr(:))) THEN ! n.b. ANY((//)) = .FALSE.
3674 var(i) = vol7d_var_miss
3675 ENDIF
3676 ENDDO
3677 end if
3678ENDIF
3679
3680END SUBROUTINE delete_var
3681
3682END SUBROUTINE vol7d_peeling
3683
3684
3685end module modqc
Variables user in Quality Control.
Definition: modqc.F90:386
Test di dato invalidato.
Definition: modqc.F90:411
Remove data under a defined grade of confidence.
Definition: modqc.F90:391
Check data validity based on single confidence.
Definition: modqc.F90:401
Check data validity based on gross error check.
Definition: modqc.F90:406
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.
Utilities and defines for quality control.
Definition: modqc.F90:357
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:368

Generated with Doxygen.