libsim Versione 7.2.1
modqc.F90
1! Copyright (C) 2010 ARPA-SIM <urpsim@smr.arpa.emr.it>
2! authors:
3! Davide Cesari <dcesari@arpa.emr.it>
4! Paolo Patruno <ppatruno@arpa.emr.it>
5
6! This program is free software; you can redistribute it and/or
7! modify it under the terms of the GNU General Public License as
8! published by the Free Software Foundation; either version 2 of
9! the License, or (at your option) any later version.
10
11! This program is distributed in the hope that it will be useful,
12! but WITHOUT ANY WARRANTY; without even the implied warranty of
13! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14! GNU General Public License for more details.
15
16! You should have received a copy of the GNU General Public License
17! along with this program. If not, see <http://www.gnu.org/licenses/>.
18#include "config.h"
19
22
169module modqc
170use kinds
173use vol7d_class
174
175
176implicit none
177
178
180type :: qcpartype
181 integer (kind=int_b):: att
182 integer (kind=int_b):: gross_error ! special valuer for "*B33192" when gross error check failed
183 integer (kind=int_b):: invalidated ! special valuer for "*B33196" when manual invalidation happen
184end type qcpartype
185
187type(qcpartype) :: qcpar=qcpartype(10_int_b,0_int_b,1_int_b)
188
189integer, parameter :: nqcattrvars=4
190CHARACTER(len=10),parameter :: qcattrvarsbtables(nqcattrvars)=(/"*B33196","*B33192","*B33193","*B33194"/)
191
192type :: qcattrvars
193 TYPE(vol7d_var) :: vars(nqcattrvars)
194 CHARACTER(len=10) :: btables(nqcattrvars)
195end type qcattrvars
196
198interface init
199 module procedure init_qcattrvars
200end interface
201
203interface peeled
204 module procedure peeledrb, peeleddb, peeledbb, peeledib, peeledcb &
205 ,peeledri, peeleddi, peeledbi, peeledii, peeledci &
206 ,peeledrr, peeleddr, peeledbr, peeledir, peeledcr &
207 ,peeledrd, peeleddd, peeledbd, peeledid, peeledcd &
208 ,peeledrc, peeleddc, peeledbc, peeledic, peeledcc
209end interface
210
211
213interface vd
214 module procedure vdi,vdb,vdr,vdd,vdc
215end interface
216
218interface vdge
219 module procedure vdgei,vdgeb,vdger,vdged,vdgec
220end interface
221
223interface invalidated
224 module procedure invalidatedi,invalidatedb,invalidatedr,invalidatedd,invalidatedc
225end interface
226
227private
228
229public vd, vdge, init, qcattrvars_new, invalidated, peeled, vol7d_peeling
230public qcattrvars, nqcattrvars, qcattrvarsbtables
231public qcpar, qcpartype, qcsummaryflagb ! ,qcsummaryflagi
232
233contains
234
235
236! peeled routines
237#undef VOL7D_POLY_SUBTYPE
238#undef VOL7D_POLY_SUBTYPES
239#undef VOL7D_POLY_ISC
240#define VOL7D_POLY_SUBTYPE REAL
241#define VOL7D_POLY_SUBTYPES r
242
243#undef VOL7D_POLY_TYPE
244#undef VOL7D_POLY_TYPES
245#undef VOL7D_POLY_ISC
246#undef VOL7D_POLY_TYPES_SUBTYPES
247#define VOL7D_POLY_TYPE REAL
248#define VOL7D_POLY_TYPES r
249#define VOL7D_POLY_TYPES_SUBTYPES rr
250#include "modqc_peeled_include.F90"
251#include "modqc_peel_util_include.F90"
252#undef VOL7D_POLY_TYPE
253#undef VOL7D_POLY_TYPES
254#undef VOL7D_POLY_TYPES_SUBTYPES
255#define VOL7D_POLY_TYPE DOUBLE PRECISION
256#define VOL7D_POLY_TYPES d
257#define VOL7D_POLY_TYPES_SUBTYPES dr
258#include "modqc_peeled_include.F90"
259#include "modqc_peel_util_include.F90"
260#undef VOL7D_POLY_TYPE
261#undef VOL7D_POLY_TYPES
262#undef VOL7D_POLY_TYPES_SUBTYPES
263#define VOL7D_POLY_TYPE INTEGER
264#define VOL7D_POLY_TYPES i
265#define VOL7D_POLY_TYPES_SUBTYPES ir
266#include "modqc_peeled_include.F90"
267#include "modqc_peel_util_include.F90"
268#undef VOL7D_POLY_TYPE
269#undef VOL7D_POLY_TYPES
270#undef VOL7D_POLY_TYPES_SUBTYPES
271#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
272#define VOL7D_POLY_TYPES b
273#define VOL7D_POLY_TYPES_SUBTYPES br
274#include "modqc_peeled_include.F90"
275#include "modqc_peel_util_include.F90"
276#undef VOL7D_POLY_TYPE
277#undef VOL7D_POLY_TYPES
278#undef VOL7D_POLY_TYPES_SUBTYPES
279#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
280#define VOL7D_POLY_TYPES c
281#define VOL7D_POLY_ISC = 1
282#define VOL7D_POLY_TYPES_SUBTYPES cr
283#include "modqc_peeled_include.F90"
284#include "modqc_peel_util_include.F90"
285
286
287#undef VOL7D_POLY_SUBTYPE
288#undef VOL7D_POLY_SUBTYPES
289#undef VOL7D_POLY_ISC
290#define VOL7D_POLY_SUBTYPE DOUBLE PRECISION
291#define VOL7D_POLY_SUBTYPES d
292
293#undef VOL7D_POLY_TYPE
294#undef VOL7D_POLY_TYPES
295#undef VOL7D_POLY_TYPES_SUBTYPES
296#define VOL7D_POLY_TYPE REAL
297#define VOL7D_POLY_TYPES r
298#define VOL7D_POLY_TYPES_SUBTYPES rd
299#include "modqc_peeled_include.F90"
300#undef VOL7D_POLY_TYPE
301#undef VOL7D_POLY_TYPES
302#undef VOL7D_POLY_TYPES_SUBTYPES
303#define VOL7D_POLY_TYPE DOUBLE PRECISION
304#define VOL7D_POLY_TYPES d
305#define VOL7D_POLY_TYPES_SUBTYPES dd
306#include "modqc_peeled_include.F90"
307#undef VOL7D_POLY_TYPE
308#undef VOL7D_POLY_TYPES
309#undef VOL7D_POLY_TYPES_SUBTYPES
310#define VOL7D_POLY_TYPE INTEGER
311#define VOL7D_POLY_TYPES i
312#define VOL7D_POLY_TYPES_SUBTYPES id
313#include "modqc_peeled_include.F90"
314#undef VOL7D_POLY_TYPE
315#undef VOL7D_POLY_TYPES
316#undef VOL7D_POLY_TYPES_SUBTYPES
317#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
318#define VOL7D_POLY_TYPES b
319#define VOL7D_POLY_TYPES_SUBTYPES bd
320#include "modqc_peeled_include.F90"
321#undef VOL7D_POLY_TYPE
322#undef VOL7D_POLY_TYPES
323#undef VOL7D_POLY_TYPES_SUBTYPES
324#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
325#define VOL7D_POLY_TYPES c
326#define VOL7D_POLY_TYPES_SUBTYPES cd
327#include "modqc_peeled_include.F90"
328
329
330#undef VOL7D_POLY_SUBTYPE
331#undef VOL7D_POLY_SUBTYPES
332#undef VOL7D_POLY_ISC
333#define VOL7D_POLY_SUBTYPE INTEGER
334#define VOL7D_POLY_SUBTYPES i
335
336#undef VOL7D_POLY_TYPE
337#undef VOL7D_POLY_TYPES
338#undef VOL7D_POLY_TYPES_SUBTYPES
339#define VOL7D_POLY_TYPE REAL
340#define VOL7D_POLY_TYPES r
341#define VOL7D_POLY_TYPES_SUBTYPES ri
342#include "modqc_peeled_include.F90"
343#undef VOL7D_POLY_TYPE
344#undef VOL7D_POLY_TYPES
345#undef VOL7D_POLY_TYPES_SUBTYPES
346#define VOL7D_POLY_TYPE DOUBLE PRECISION
347#define VOL7D_POLY_TYPES d
348#define VOL7D_POLY_TYPES_SUBTYPES di
349#include "modqc_peeled_include.F90"
350#undef VOL7D_POLY_TYPE
351#undef VOL7D_POLY_TYPES
352#undef VOL7D_POLY_TYPES_SUBTYPES
353#define VOL7D_POLY_TYPE INTEGER
354#define VOL7D_POLY_TYPES i
355#define VOL7D_POLY_TYPES_SUBTYPES ii
356#include "modqc_peeled_include.F90"
357#undef VOL7D_POLY_TYPE
358#undef VOL7D_POLY_TYPES
359#undef VOL7D_POLY_TYPES_SUBTYPES
360#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
361#define VOL7D_POLY_TYPES b
362#define VOL7D_POLY_TYPES_SUBTYPES bi
363#include "modqc_peeled_include.F90"
364#undef VOL7D_POLY_TYPE
365#undef VOL7D_POLY_TYPES
366#undef VOL7D_POLY_TYPES_SUBTYPES
367#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
368#define VOL7D_POLY_TYPES c
369#define VOL7D_POLY_ISC = 1
370#define VOL7D_POLY_TYPES_SUBTYPES ci
371#include "modqc_peeled_include.F90"
372
373
374#undef VOL7D_POLY_SUBTYPE
375#undef VOL7D_POLY_SUBTYPES
376#undef VOL7D_POLY_ISC
377#define VOL7D_POLY_SUBTYPE INTEGER(kind=int_b)
378#define VOL7D_POLY_SUBTYPES b
379
380#undef VOL7D_POLY_TYPE
381#undef VOL7D_POLY_TYPES
382#undef VOL7D_POLY_TYPES_SUBTYPES
383#define VOL7D_POLY_TYPE REAL
384#define VOL7D_POLY_TYPES r
385#define VOL7D_POLY_TYPES_SUBTYPES rb
386#include "modqc_peeled_include.F90"
387#undef VOL7D_POLY_TYPE
388#undef VOL7D_POLY_TYPES
389#undef VOL7D_POLY_TYPES_SUBTYPES
390#define VOL7D_POLY_TYPE DOUBLE PRECISION
391#define VOL7D_POLY_TYPES d
392#define VOL7D_POLY_TYPES_SUBTYPES db
393#include "modqc_peeled_include.F90"
394#undef VOL7D_POLY_TYPE
395#undef VOL7D_POLY_TYPES
396#undef VOL7D_POLY_TYPES_SUBTYPES
397#define VOL7D_POLY_TYPE INTEGER
398#define VOL7D_POLY_TYPES i
399#define VOL7D_POLY_TYPES_SUBTYPES ib
400#include "modqc_peeled_include.F90"
401#undef VOL7D_POLY_TYPE
402#undef VOL7D_POLY_TYPES
403#undef VOL7D_POLY_TYPES_SUBTYPES
404#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
405#define VOL7D_POLY_TYPES b
406#define VOL7D_POLY_TYPES_SUBTYPES bb
407#include "modqc_peeled_include.F90"
408#undef VOL7D_POLY_TYPE
409#undef VOL7D_POLY_TYPES
410#undef VOL7D_POLY_TYPES_SUBTYPES
411#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
412#define VOL7D_POLY_TYPES c
413#define VOL7D_POLY_ISC = 1
414#define VOL7D_POLY_TYPES_SUBTYPES cb
415#include "modqc_peeled_include.F90"
416
417
418#undef VOL7D_POLY_SUBTYPE
419#undef VOL7D_POLY_SUBTYPES
420#undef VOL7D_POLY_ISC
421#define VOL7D_POLY_SUBTYPE CHARACTER(len=vol7d_cdatalen)
422#define VOL7D_POLY_SUBTYPES c
423
424#undef VOL7D_POLY_TYPE
425#undef VOL7D_POLY_TYPES
426#undef VOL7D_POLY_TYPES_SUBTYPES
427#define VOL7D_POLY_TYPE REAL
428#define VOL7D_POLY_TYPES r
429#define VOL7D_POLY_TYPES_SUBTYPES rc
430#include "modqc_peeled_include.F90"
431#undef VOL7D_POLY_TYPE
432#undef VOL7D_POLY_TYPES
433#undef VOL7D_POLY_TYPES_SUBTYPES
434#define VOL7D_POLY_TYPE DOUBLE PRECISION
435#define VOL7D_POLY_TYPES d
436#define VOL7D_POLY_TYPES_SUBTYPES dc
437#include "modqc_peeled_include.F90"
438#undef VOL7D_POLY_TYPE
439#undef VOL7D_POLY_TYPES
440#undef VOL7D_POLY_TYPES_SUBTYPES
441#define VOL7D_POLY_TYPE INTEGER
442#define VOL7D_POLY_TYPES i
443#define VOL7D_POLY_TYPES_SUBTYPES ic
444#include "modqc_peeled_include.F90"
445#undef VOL7D_POLY_TYPE
446#undef VOL7D_POLY_TYPES
447#undef VOL7D_POLY_TYPES_SUBTYPES
448#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
449#define VOL7D_POLY_TYPES b
450#define VOL7D_POLY_TYPES_SUBTYPES bc
451#include "modqc_peeled_include.F90"
452#undef VOL7D_POLY_TYPE
453#undef VOL7D_POLY_TYPES
454#undef VOL7D_POLY_TYPES_SUBTYPES
455#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
456#define VOL7D_POLY_TYPES c
457#define VOL7D_POLY_ISC = 1
458#define VOL7D_POLY_TYPES_SUBTYPES cc
459#include "modqc_peeled_include.F90"
460
461
462subroutine init_qcattrvars(this)
463
464type(qcattrvars),intent(inout) :: this
465integer :: i
466
467this%btables(:) =qcattrvarsbtables
468do i =1, nqcattrvars
469 call init(this%vars(i),this%btables(i))
470end do
471
472end subroutine init_qcattrvars
473
474
475type(qcattrvars) function qcattrvars_new()
476
477call init(qcattrvars_new)
478
479end function qcattrvars_new
480
481
489SUBROUTINE vol7d_peeling(this, data_id, keep_attr, delete_attr, preserve, purgeana)
490TYPE(vol7d),INTENT(INOUT) :: this
491integer,INTENT(inout),pointer,OPTIONAL :: data_id(:,:,:,:,:)
492CHARACTER(len=*),INTENT(in),OPTIONAL :: keep_attr(:)
493CHARACTER(len=*),INTENT(in),OPTIONAL :: delete_attr(:)
494logical,intent(in),optional :: preserve
495logical,intent(in),optional :: purgeana
496
497integer :: inddativar,inddatiattrinv,inddatiattrcli,inddatiattrtem,inddatiattrspa,inddativarattr
498type(qcattrvars) :: attrvars
499
500INTEGER(kind=int_b),pointer :: invbb(:,:,:,:,:),clibb(:,:,:,:,:),tembb(:,:,:,:,:),spabb(:,:,:,:,:)
501INTEGER,pointer :: invbi(:,:,:,:,:),clibi(:,:,:,:,:),tembi(:,:,:,:,:),spabi(:,:,:,:,:)
502REAL,pointer :: invbr(:,:,:,:,:),clibr(:,:,:,:,:),tembr(:,:,:,:,:),spabr(:,:,:,:,:)
503DOUBLE PRECISION,pointer :: invbd(:,:,:,:,:),clibd(:,:,:,:,:),tembd(:,:,:,:,:),spabd(:,:,:,:,:)
504CHARACTER(len=vol7d_cdatalen),pointer :: invbc(:,:,:,:,:),clibc(:,:,:,:,:),tembc(:,:,:,:,:),spabc(:,:,:,:,:)
505
506call l4f_log(l4f_info,'starting peeling')
507
508call init(attrvars)
509
510! generate code per i vari tipi di dati di v7d
511! tramite un template e il preprocessore
512
513
514#undef VOL7D_POLY_SUBTYPE
515#undef VOL7D_POLY_SUBTYPES
516#define VOL7D_POLY_SUBTYPE REAL
517#define VOL7D_POLY_SUBTYPES r
518
519#undef VOL7D_POLY_TYPE
520#undef VOL7D_POLY_TYPES
521#define VOL7D_POLY_TYPE REAL
522#define VOL7D_POLY_TYPES r
523#include "modqc_peeling_include.F90"
524#undef VOL7D_POLY_TYPE
525#undef VOL7D_POLY_TYPES
526#define VOL7D_POLY_TYPE DOUBLE PRECISION
527#define VOL7D_POLY_TYPES d
528#include "modqc_peeling_include.F90"
529#undef VOL7D_POLY_TYPE
530#undef VOL7D_POLY_TYPES
531#define VOL7D_POLY_TYPE INTEGER
532#define VOL7D_POLY_TYPES i
533#include "modqc_peeling_include.F90"
534#undef VOL7D_POLY_TYPE
535#undef VOL7D_POLY_TYPES
536#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
537#define VOL7D_POLY_TYPES b
538#include "modqc_peeling_include.F90"
539#undef VOL7D_POLY_TYPE
540#undef VOL7D_POLY_TYPES
541#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
542#define VOL7D_POLY_TYPES c
543#include "modqc_peeling_include.F90"
544
545
546#undef VOL7D_POLY_SUBTYPE
547#undef VOL7D_POLY_SUBTYPES
548#define VOL7D_POLY_SUBTYPE DOUBLE PRECISION
549#define VOL7D_POLY_SUBTYPES d
550
551#undef VOL7D_POLY_TYPE
552#undef VOL7D_POLY_TYPES
553#define VOL7D_POLY_TYPE REAL
554#define VOL7D_POLY_TYPES r
555#include "modqc_peeling_include.F90"
556#undef VOL7D_POLY_TYPE
557#undef VOL7D_POLY_TYPES
558#define VOL7D_POLY_TYPE DOUBLE PRECISION
559#define VOL7D_POLY_TYPES d
560#include "modqc_peeling_include.F90"
561#undef VOL7D_POLY_TYPE
562#undef VOL7D_POLY_TYPES
563#define VOL7D_POLY_TYPE INTEGER
564#define VOL7D_POLY_TYPES i
565#include "modqc_peeling_include.F90"
566#undef VOL7D_POLY_TYPE
567#undef VOL7D_POLY_TYPES
568#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
569#define VOL7D_POLY_TYPES b
570#include "modqc_peeling_include.F90"
571#undef VOL7D_POLY_TYPE
572#undef VOL7D_POLY_TYPES
573#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
574#define VOL7D_POLY_TYPES c
575#include "modqc_peeling_include.F90"
576
577
578#undef VOL7D_POLY_SUBTYPE
579#undef VOL7D_POLY_SUBTYPES
580#define VOL7D_POLY_SUBTYPE INTEGER
581#define VOL7D_POLY_SUBTYPES i
582
583#undef VOL7D_POLY_TYPE
584#undef VOL7D_POLY_TYPES
585#define VOL7D_POLY_TYPE REAL
586#define VOL7D_POLY_TYPES r
587#include "modqc_peeling_include.F90"
588#undef VOL7D_POLY_TYPE
589#undef VOL7D_POLY_TYPES
590#define VOL7D_POLY_TYPE DOUBLE PRECISION
591#define VOL7D_POLY_TYPES d
592#include "modqc_peeling_include.F90"
593#undef VOL7D_POLY_TYPE
594#undef VOL7D_POLY_TYPES
595#define VOL7D_POLY_TYPE INTEGER
596#define VOL7D_POLY_TYPES i
597#include "modqc_peeling_include.F90"
598#undef VOL7D_POLY_TYPE
599#undef VOL7D_POLY_TYPES
600#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
601#define VOL7D_POLY_TYPES b
602#include "modqc_peeling_include.F90"
603#undef VOL7D_POLY_TYPE
604#undef VOL7D_POLY_TYPES
605#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
606#define VOL7D_POLY_TYPES c
607#include "modqc_peeling_include.F90"
608
609
610#undef VOL7D_POLY_SUBTYPE
611#undef VOL7D_POLY_SUBTYPES
612#define VOL7D_POLY_SUBTYPE INTEGER(kind=int_b)
613#define VOL7D_POLY_SUBTYPES b
614
615#undef VOL7D_POLY_TYPE
616#undef VOL7D_POLY_TYPES
617#define VOL7D_POLY_TYPE REAL
618#define VOL7D_POLY_TYPES r
619#include "modqc_peeling_include.F90"
620#undef VOL7D_POLY_TYPE
621#undef VOL7D_POLY_TYPES
622#define VOL7D_POLY_TYPE DOUBLE PRECISION
623#define VOL7D_POLY_TYPES d
624#include "modqc_peeling_include.F90"
625#undef VOL7D_POLY_TYPE
626#undef VOL7D_POLY_TYPES
627#define VOL7D_POLY_TYPE INTEGER
628#define VOL7D_POLY_TYPES i
629#include "modqc_peeling_include.F90"
630#undef VOL7D_POLY_TYPE
631#undef VOL7D_POLY_TYPES
632#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
633#define VOL7D_POLY_TYPES b
634#include "modqc_peeling_include.F90"
635#undef VOL7D_POLY_TYPE
636#undef VOL7D_POLY_TYPES
637#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
638#define VOL7D_POLY_TYPES c
639#include "modqc_peeling_include.F90"
640
641
642
643#undef VOL7D_POLY_SUBTYPE
644#undef VOL7D_POLY_SUBTYPES
645#define VOL7D_POLY_SUBTYPE CHARACTER(len=vol7d_cdatalen)
646#define VOL7D_POLY_SUBTYPES c
647
648#undef VOL7D_POLY_TYPE
649#undef VOL7D_POLY_TYPES
650#define VOL7D_POLY_TYPE REAL
651#define VOL7D_POLY_TYPES r
652#include "modqc_peeling_include.F90"
653#undef VOL7D_POLY_TYPE
654#undef VOL7D_POLY_TYPES
655#define VOL7D_POLY_TYPE DOUBLE PRECISION
656#define VOL7D_POLY_TYPES d
657#include "modqc_peeling_include.F90"
658#undef VOL7D_POLY_TYPE
659#undef VOL7D_POLY_TYPES
660#define VOL7D_POLY_TYPE INTEGER
661#define VOL7D_POLY_TYPES i
662#include "modqc_peeling_include.F90"
663#undef VOL7D_POLY_TYPE
664#undef VOL7D_POLY_TYPES
665#define VOL7D_POLY_TYPE INTEGER(kind=int_b)
666#define VOL7D_POLY_TYPES b
667#include "modqc_peeling_include.F90"
668#undef VOL7D_POLY_TYPE
669#undef VOL7D_POLY_TYPES
670#define VOL7D_POLY_TYPE CHARACTER(len=vol7d_cdatalen)
671#define VOL7D_POLY_TYPES c
672#include "modqc_peeling_include.F90"
673
674
675
676IF (.NOT.PRESENT(keep_attr) .AND. .NOT.PRESENT(delete_attr) .and. .not. optio_log(preserve)) THEN ! destroy all attributes
677 IF (ASSOCIATED(this%voldatiattrr)) DEALLOCATE(this%voldatiattrr)
678 IF (ASSOCIATED(this%voldatiattrd)) DEALLOCATE(this%voldatiattrd)
679 IF (ASSOCIATED(this%voldatiattri)) DEALLOCATE(this%voldatiattri)
680 IF (ASSOCIATED(this%voldatiattrb)) DEALLOCATE(this%voldatiattrb)
681 IF (ASSOCIATED(this%voldatiattrc)) DEALLOCATE(this%voldatiattrc)
682
683 CALL delete(this%datiattr)
684 CALL delete(this%dativarattr)
685END IF
686
687IF (PRESENT(keep_attr)) THEN ! set to missing non requested attributes and reform
688
689 if (optio_log(preserve)) call l4f_log(l4f_warn,"preserve parameter ignored: keep_attr passed")
690 CALL keep_var(this%datiattr%r)
691 CALL keep_var(this%datiattr%d)
692 CALL keep_var(this%datiattr%i)
693 CALL keep_var(this%datiattr%b)
694 CALL keep_var(this%datiattr%c)
695 CALL qc_reform(this,data_id, miss=.true., purgeana=purgeana)
696
697ELSE IF (PRESENT(delete_attr)) THEN ! set to missing requested attributes and reform
698
699 if (optio_log(preserve)) call l4f_log(l4f_warn,"preserve parameter ignored: delete_attr passed")
700 CALL delete_var(this%datiattr%r)
701 CALL delete_var(this%datiattr%d)
702 CALL delete_var(this%datiattr%i)
703 CALL delete_var(this%datiattr%b)
704 CALL delete_var(this%datiattr%c)
705 CALL qc_reform(this,data_id, miss=.true., purgeana=purgeana)
706
707ELSE IF (PRESENT(purgeana)) THEN
708
709 CALL qc_reform(this,data_id, purgeana=purgeana)
710
711ENDIF
712
713
714CONTAINS
716
718subroutine qc_reform(this,data_id,miss, purgeana)
719TYPE(vol7d),INTENT(INOUT) :: this
720integer,INTENT(inout),pointer,OPTIONAL :: data_id(:,:,:,:,:)
721logical,intent(in),optional :: miss
722logical,intent(in),optional :: purgeana
723
724integer,pointer :: data_idtmp(:,:,:,:,:)
725logical,allocatable :: llana(:)
726integer,allocatable :: anaind(:)
727integer :: i,j,nana
728
729if (optio_log(purgeana)) then
730 allocate(llana(size(this%ana)))
731 llana =.false.
732 do i =1,size(this%ana)
733 if (associated(this%voldatii)) llana(i)= llana(i) .or. any(c_e(this%voldatii(i,:,:,:,:,:)))
734 if (associated(this%voldatir)) llana(i)= llana(i) .or. any(c_e(this%voldatir(i,:,:,:,:,:)))
735 if (associated(this%voldatid)) llana(i)= llana(i) .or. any(c_e(this%voldatid(i,:,:,:,:,:)))
736 if (associated(this%voldatib)) llana(i)= llana(i) .or. any(c_e(this%voldatib(i,:,:,:,:,:)))
737 if (associated(this%voldatic)) llana(i)= llana(i) .or. any(c_e(this%voldatic(i,:,:,:,:,:)))
738
739#ifdef DEBUG
740 if (.not. llana(i)) call l4f_log(l4f_debug,"remove station"//t2c(i))
741#endif
742
743 end do
744
745 nana=count(llana)
747
748 allocate(anaind(nana))
749
750 j=0
751 do i=1,size(this%ana)
752 if (llana(i)) then
753 j=j+1
754 anaind(j)=i
755 end if
756 end do
757
758
759 if(present(data_id)) then
760 allocate(data_idtmp(nana,size(data_id,2),size(data_id,3),size(data_id,4),size(data_id,5)))
761 data_idtmp=data_id(anaind,:,:,:,:)
762 if (associated(data_id))deallocate(data_id)
763 data_id=>data_idtmp
764 end if
765
766 call vol7d_reform(this,miss=miss,lana=llana)
767
768 deallocate(llana,anaind)
769
770else
771
772 call vol7d_reform(this,miss=miss)
773
774end if
775
776end subroutine qc_reform
777
778
779SUBROUTINE keep_var(var)
780TYPE(vol7d_var),intent(inout),POINTER :: var(:)
781
782INTEGER :: i
783
784IF (ASSOCIATED(var)) THEN
785 if (size(var) == 0) then
786 var%btable = vol7d_var_miss%btable
787 else
788 DO i = 1, SIZE(var)
789 IF (all(var(i)%btable /= keep_attr(:))) THEN ! n.b. ALL((//)) = .TRUE.
790 var(i)%btable = vol7d_var_miss%btable
791 ENDIF
792 ENDDO
793 end if
794ENDIF
795
796END SUBROUTINE keep_var
797
798SUBROUTINE delete_var(var)
799TYPE(vol7d_var),intent(inout),POINTER :: var(:)
800
801INTEGER :: i
802
803IF (ASSOCIATED(var)) THEN
804 if (size(var) == 0) then
805 var%btable = vol7d_var_miss%btable
806 else
807 DO i = 1, SIZE(var)
808 IF (any(var(i)%btable == delete_attr(:))) THEN ! n.b. ANY((//)) = .FALSE.
809 var(i) = vol7d_var_miss
810 ENDIF
811 ENDDO
812 end if
813ENDIF
814
815END SUBROUTINE delete_var
816
817END SUBROUTINE vol7d_peeling
818
819
820end module modqc
Function to check whether a value is missing or not.
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
real data conversion
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.