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