87 character (len=255),
parameter:: subcategorytem=
"QCtem"
89 integer,
parameter :: tem_nvar=1
90 CHARACTER(len=10) :: tem_btable(tem_nvar)=(/
"B12101"/)
92 real,
parameter :: tem_a(tem_nvar) = (/1.e5/)
94 real,
parameter :: tem_b(tem_nvar) = (/250./)
98 type (vol7d),
pointer :: v7d => null()
99 integer,
pointer :: data_id_in(:,:,:,:,:) => null()
100 integer,
pointer :: data_id_out(:,:,:,:,:) => null()
102 type (qcclitype) :: qccli
103 type (vol7d) :: clima
104 character(len=20):: operation
105 integer :: timeconfidence
111 module procedure qcteminit
116 module procedure qctemalloc
121 module procedure qctemdelete
134 subroutine qcteminit(qctem,v7d,var, timei, timef, coordmin, coordmax, data_id_in,extremepath,temporalpath,&
136 dsne,usere,passworde,&
137 dsntem,usertem,passwordtem,&
139 height2level,operation,timeconfidence,categoryappend)
141 type(qctemtype),
intent(in out) :: qctem
142 type (vol7d),
intent(in),
target:: v7d
143 character(len=*),
INTENT(in) :: var(:)
146 TYPE(geo_coord),
INTENT(inout),
optional :: coordmin,coordmax
148 TYPE(datetime),
INTENT(in),
optional :: timei, timef
149 integer,
intent(in),
optional,
target:: data_id_in(:,:,:,:,:)
150 character(len=*),
intent(in),
optional :: extremepath
151 character(len=*),
intent(in),
optional :: temporalpath
152 logical ,
intent(in),
optional :: height2level
153 character(len=*),
optional :: operation
154 integer,
intent(in),
optional :: timeconfidence
155 character(len=*),
INTENT(in),
OPTIONAL :: categoryappend
158 type (vol7d_dballe) :: v7d_dballetmp
159 character(len=*),
intent(in),
optional :: dsne
160 character(len=*),
intent(in),
optional :: usere
161 character(len=*),
intent(in),
optional :: passworde
162 character(len=*),
intent(in),
optional :: dsntem
163 character(len=*),
intent(in),
optional :: usertem
164 character(len=*),
intent(in),
optional :: passwordtem
165 character(len=512) :: ldsntem
166 character(len=512) :: lusertem
167 character(len=512) :: lpasswordtem
170 type (vol7d) :: v7dtmp
171 TYPE(vol7d_network):: network
173 character(len=512) :: filepathtem
174 character(len=512) :: a_name
175 character(len=9) ::netname(2)=(/
"qctemgndi",
"qctemsndi"/)
178 call l4f_launcher(a_name,a_name_append=trim(subcategorytem)//
"."//trim(categoryappend))
179 qctem%category=l4f_category_get(a_name)
181 nullify ( qctem%data_id_in )
182 nullify ( qctem%data_id_out )
187 qctem%operation=optio_c(operation,20)
188 filepathtem=optio_c(temporalpath,512)
191 if (qctem%operation /=
"gradient" .and. qctem%operation /=
"run")
then
192 call l4f_category_log(qctem%category,l4f_error,
"operation is wrong: "//qctem%operation)
197 if (
present(data_id_in))
then
198 qctem%data_id_in => data_id_in
201 qctem%timeconfidence = optio_i(timeconfidence)
204 call init(qctem%qccli,v7d,var, timei, timef, data_id_in,&
205 macropath=cmiss, climapath=cmiss, extremepath=extremepath, &
207 dsncli=cmiss,dsnextreme=dsne,user=usere,password=passworde,&
209 height2level=height2level,categoryappend=categoryappend)
216 if (qctem%operation ==
"run")
then
217 call init(network,netname(i))
220 call optio(dsntem,ldsntem)
221 call optio(usertem,lusertem)
222 call optio(passwordtem,lpasswordtem)
224 if (
c_e(filepathtem) .and. (
c_e(ldsntem).or.
c_e(lusertem).or.
c_e(lpasswordtem)))
then
225 call l4f_category_log(qctem%category,l4f_error,
"filepath defined together with dba options")
229 if (.not.
c_e(ldsntem))
then
233 if (.not.
c_e(filepathtem))
then
234 filepathtem=get_package_filepath(netname(i)//
'.v7d', filetype_data)
237 if (
c_e(filepathtem))
then
239 select case (trim(lowercase(suffixname(filepathtem))))
243 call import(v7dtmp,filename=filepathtem,unit=iuni)
248 call init(v7d_dballetmp,file=.true.,filename=filepathtem,categoryappend=trim(a_name)//
".clima")
249 call import(v7d_dballetmp,var=var, &
250 varkind=(/(
"r",i=1,
size(var))/),attr=(/
"*B33209"/),attrkind=(/
"b"/),network=network)
251 call copy(v7d_dballetmp%vol7d,v7dtmp)
252 call delete(v7d_dballetmp)
257 "file type not supported (use .v7d or .bufr suffix only): "//trim(filepathtem))
262 call l4f_category_log(qctem%category,l4f_warn,
"spatial clima volume not iniziatized: spatial QC will not be possible")
263 call init(qctem%clima)
264 call raise_fatal_error()
271 call init(v7d_dballetmp,dsn=ldsntem,user=lusertem,password=lpasswordtem,write=.false.,&
272 file=.false.,categoryappend=trim(a_name)//
".tem")
274 call import(v7d_dballetmp,var=var, &
275 varkind=(/("r
",i=1,size(var))/),attr=(/"*b33209
"/),attrkind=(/"b
"/),network=network)
277 call copy(v7d_dballetmp%vol7d,v7dtmp)
279 call delete(v7d_dballetmp)
283 call vol7d_merge(qctem%clima,v7dtmp)
290 end subroutine qcteminit
294 !>\brief Allocazioni di memoria
295 subroutine qctemalloc(qctem)
296 ! pseudo costruttore con distruttore automatico
298 type(qctemtype),intent(in out) :: qctem !< Oggetto per il controllo climatico
303 ! se ti sei dimenticato di deallocare ci penso io
304 call qctemdealloc(qctem)
306 if (associated(qctem%data_id_in))then
307 sh=shape(qctem%data_id_in)
308 allocate (qctem%data_id_out(sh(1),sh(2),sh(3),sh(4),sh(5)),stat=istatt)
310 call l4f_category_log(qctem%category,L4F_ERROR,"allocate error
")
311 call raise_error("allocate error
")
313 qctem%data_id_out=imiss
317 end subroutine qctemalloc
320 !>\brief Deallocazione della memoria
321 subroutine qctemdealloc(qctem)
324 type(qctemtype),intent(in out) :: qctem !< Oggetto per il controllo temporale
326 if (associated(qctem%data_id_out)) deallocate (qctem%data_id_out)
328 if (associated(qctem%data_id_out)) then
329 deallocate (qctem%data_id_out)
330 nullify (qctem%data_id_out)
333 end subroutine qctemdealloc
336 !>\brief Cancellazione
339 subroutine qctemdelete(qctem)
340 ! decostruttore a mezzo
341 type(qctemtype),intent(in out) :: qctem !< Oggetto per il controllo temporale
343 call qctemdealloc(qctem)
345 call delete(qctem%qccli)
348 call l4f_category_delete(qctem%category)
351 end subroutine qctemdelete
354 à
!>\brief Controllo di Qualit temporale.
355 èà
!!Questo il vero e proprio controllo di qualit temporale.
357 SUBROUTINE quacontem (qctem,battrinv,battrcli,battrout,&
358 anamask,timemask,levelmask,timerangemask,varmask,networkmask)
361 à
type(qctemtype),intent(in out) :: qctem !< Oggetto per il controllo di qualit
362 character (len=10) ,intent(in),optional :: battrinv !< attributo invalidated in input
363 character (len=10) ,intent(in),optional :: battrcli !< attributo con la confidenza climatologica in input
364 character (len=10) ,intent(in),optional :: battrout !< attributo con la confidenza temporale in output
365 logical ,intent(in),optional :: anamask(:) !< Filtro sulle anagrafiche
366 logical ,intent(in),optional :: timemask(:) !< Filtro sul tempo
367 logical ,intent(in),optional :: levelmask(:) !< Filtro sui livelli
368 logical ,intent(in),optional :: timerangemask(:) !< filtro sui timerange
369 logical ,intent(in),optional :: varmask(:) !< Filtro sulle variabili
370 logical ,intent(in),optional :: networkmask(:) !< Filtro sui network
372 !REAL(kind=fp_geo) :: lat,lon
375 integer :: indbattrinv,indbattrcli,indbattrout,grunit
376 logical :: anamaskl(size(qctem%v7d%ana)), timemaskl(size(qctem%v7d%time)), levelmaskl(size(qctem%v7d%level)), &
377 timerangemaskl(size(qctem%v7d%timerange)), varmaskl(size(qctem%v7d%dativar%r)), networkmaskl(size(qctem%v7d%network))
379 integer :: indana , indtime ,indlevel ,indtimerange ,inddativarr, indnetwork,indtimenear
380 integer :: indcana , indctime ,indclevel ,indctimerange ,indcdativarr, indcnetwork, indcnetworks, indcnetworkg
381 real :: datoqui,datoprima,datodopo,climaquii, climaquif
382 !integer, allocatable :: indcanav(:)
384 !TYPE(vol7d_ana) :: ana
385 TYPE(datetime) :: time,prima, ora, dopo
386 TYPE(vol7d_network):: network
387 type(timedelta) :: td
389 double precision :: gradprima,graddopo,grad
390 !call qctem_validate (qctem)
391 character(len=512) :: filename
395 !localize optional parameter
396 if (present(battrinv))then
397 indbattrinv = index_c(qctem%v7d%datiattr%b(:)%btable, battrinv)
399 indbattrinv = index_c(qctem%v7d%datiattr%b(:)%btable, qcattrvarsbtables(1))
402 if (present(battrcli))then
403 indbattrcli = index_c(qctem%v7d%datiattr%b(:)%btable, battrcli)
405 indbattrcli = index_c(qctem%v7d%datiattr%b(:)%btable, qcattrvarsbtables(2))
408 if (present(battrout))then
409 indbattrout = index_c(qctem%v7d%datiattr%b(:)%btable, battrout)
411 indbattrout = index_c(qctem%v7d%datiattr%b(:)%btable, qcattrvarsbtables(3))
415 ! some checks on input
416 .or..or.
!if (indbattrinv <=0 indbattrcli <= 0 indbattrout <= 0 ) then
417 if (indbattrout <= 0 ) then
419 call l4f_category_log(qctem%category,L4F_ERROR,"error finding attribute
index for output
")
424 if (qctem%operation == "gradient
") then
426 !check for gradient operation
427 .or.
if ( size(qctem%v7d%level) > 1 &
428 .or.
size(qctem%v7d%timerange) > 1 &
429 size(qctem%v7d%dativar%r) > 1 ) then
430 call l4f_category_log(qctem%category,L4F_ERROR,"gradient operation manage one level/timerange/var only
")
434 !check for data to check
435 if ( size(qctem%v7d%time) < 1 ) then
436 call l4f_category_log(qctem%category,L4F_INFO,"no data present for gradient operation
")
441 ! set other local variable from optional parameter
442 if(present(anamask)) then
447 if(present(timemask)) then
452 if(present(levelmask)) then
453 levelmaskl = levelmask
457 if(present(timerangemask)) then
458 timerangemaskl = timerangemask
460 timerangemaskl = .true.
462 if(present(varmask)) then
467 if(present(networkmask)) then
468 networkmaskl = networkmask
470 networkmaskl = .true.
473 ! do not touch data that will not in/validated by QC
474 qctem%v7d%voldatiattrb(:,:,:,:,:,:,indbattrout)=ibmiss
476 ! normalize data in space and time
477 call vol7d_normalize_data(qctem%qccli)
480 ! compute some index for temporal clima
481 !! compute the conventional generic datetime
482 !!cyclicdt = cyclicdatetime_new(chardate="/////////
") !TMMGGhhmm
483 time=cyclicdatetime_to_conventional(cyclicdatetime_new(chardate="/////////
")) !TMMGGhhmm
484 !!call init(time, year=1007, month=1, day=1, hour=01, minute=01)
487 if (qctem%operation == "run
") then
488 !!indcana = firsttrue(qctem%clima%ana == ana)
489 call init(network,"qctemsndi
")
490 indcnetworks = index(qctem%clima%network , network)
491 call init(network,"qctemgndi
")
492 indcnetworkg = index(qctem%clima%network , network)
493 indctime = index(qctem%clima%time , time)
496 do indana=1,size(qctem%v7d%ana)
497 .not.
if (anamaskl(indana)) cycle
498 call l4f_category_log(qctem%category,L4F_INFO,&
499 "check ana:
"//to_char(qctem%v7d%ana(indana)) )
501 do indnetwork=1,size(qctem%v7d%network)
502 do indlevel=1,size(qctem%v7d%level)
503 do indtimerange=1,size(qctem%v7d%timerange)
504 do inddativarr=1,size(qctem%v7d%dativar%r)
505 ind=index_c(tem_btable,qctem%v7d%dativar%r(inddativarr)%btable)
507 if (qctem%operation == "gradient
") then
508 ! open file to write gradient
510 !t2c(getilon(qctem%v7d%ana(indana)%coord))
511 !t2c(getilat(qctem%v7d%ana(indana)%coord))
513 filename=trim(to_char(qctem%v7d%level(indlevel)))//&
514 "_
"//trim(to_char(qctem%v7d%timerange(indtimerange)))//&
515 "_
"//trim(qctem%v7d%dativar%r(inddativarr)%btable)//&
518 call l4f_category_log(qctem%category,L4F_INFO,"try to open gradient file; filename below
")
519 call l4f_category_log(qctem%category,L4F_INFO,filename)
521 inquire(file=filename, exist=exist)
524 if (grunit /= -1) then
525 !open (unit=grunit, file=t2c(timei)//"_
"//t2c(timef)//".grad
",STATUS='UNKNOWN', form='FORMATTED')
526 open (grunit, file=filename ,STATUS='UNKNOWN', form='FORMATTED',position='APPEND')
528 ! say we have to write header in file
529 .not.
if ( exist) then
530 call l4f_category_log(qctem%category,L4F_INFO,"write header in gradient file
")
532 qctem%v7d%level(indlevel), &
533 qctem%v7d%timerange(indtimerange), &
534 qctem%v7d%dativar%r(inddativarr)
538 !!$ call l4f_log(L4F_INFO,"index:
"// t2c(indana)//t2c(indnetwork)//t2c(indlevel)//&
539 !!$ t2c(indtimerange)//t2c(inddativarr)//t2c(indtime))
542 do indtime=2,size(qctem%v7d%time)-1
544 .not..or..not..or.
if (timemaskl(indtime) levelmaskl(indlevel) &
545 .not..or..not..or..not.
timerangemaskl(indtimerange) varmaskl(inddativarr) networkmaskl(indnetwork)) cycle
548 datoqui = qctem%v7d%voldatir (indana ,indtime ,indlevel ,indtimerange ,inddativarr, indnetwork )
549 .not.
if ( c_e(datoqui)) cycle
550 ora = qctem%v7d%time (indtime)
553 if (indbattrinv > 0) then
554 if( invalidated(qctem%v7d%voldatiattrb&
555 (indana,indtime,indlevel,indtimerange,inddativarr,indnetwork,indbattrinv))) then
556 call l4f_category_log(qctem%category,L4F_WARN,&
557 "it
's better to do a reform on ana to v7d after peeling, before spatial QC")
563 if (indbattrcli > 0) then
564 .not.
if( vdge(qctem%v7d%voldatiattrb&
565 (indana,indtime,indlevel,indtimerange,inddativarr,indnetwork,indbattrcli))) then
566 call l4f_category_log(qctem%category,L4F_WARN,&
567 "It's better to do a reform on ana to v7d after peeling, before spatial qc
")
574 if (qctem%operation == "run
") then
576 indclevel = index(qctem%clima%level , qctem%v7d%level(indlevel))
577 indctimerange = index(qctem%clima%timerange , qctem%v7d%timerange(indtimerange))
579 ! attenzione attenzione TODO
580 è
! se leggo da bufr il default char e non reale
581 indcdativarr = index(qctem%clima%dativar%r, qctem%v7d%dativar%r(inddativarr))
584 call l4f_log(L4F_DEBUG,"qctem
index:
"// to_char(indctime)//to_char(indclevel)//&
585 to_char(indctimerange)//to_char(indcdativarr)//to_char(indcnetworks))
587 .or..or..or.
if ( indctime <= 0 indclevel <= 0 indctimerange <= 0 indcdativarr <= 0 &
588 .or.
indcnetworks <= 0 ) cycle
591 !!$ nintime=qctem%v7d%time(indtime)+timedelta_new(minute=30)
592 !!$ CALL getval(nintime, month=mese, hour=ora)
593 !!$ call init(time, year=1001, month=mese, day=1, hour=ora, minute=00)
596 !find the nearest data in time before
597 indtimenear=indtime-1
598 datoprima = qctem%v7d%voldatir (indana ,indtimenear ,indlevel ,indtimerange ,inddativarr, indnetwork )
599 prima = qctem%v7d%time (indtimenear)
602 if (indbattrinv > 0) then
603 if( invalidated(qctem%v7d%voldatiattrb&
604 (indana,indtimenear,indlevel,indtimerange,inddativarr,indnetwork,indbattrinv))) then
610 if (indbattrcli > 0) then
611 .not.
if( vdge(qctem%v7d%voldatiattrb&
612 (indana,indtimenear,indlevel,indtimerange,inddativarr,indnetwork,indbattrcli))) then
618 !find the nearest data in time after
619 indtimenear=indtime+1
620 datodopo = qctem%v7d%voldatir (indana ,indtimenear ,indlevel ,indtimerange ,inddativarr, indnetwork )
621 dopo = qctem%v7d%time (indtimenear)
624 if (indbattrinv > 0) then
625 if( invalidated(qctem%v7d%voldatiattrb&
626 (indana,indtimenear,indlevel,indtimerange,inddativarr,indnetwork,indbattrinv))) then
632 if (indbattrcli > 0) then
633 .not.
if( vdge(qctem%v7d%voldatiattrb&
634 (indana,indtimenear,indlevel,indtimerange,inddativarr,indnetwork,indbattrcli))) then
640 .NOT..and..NOT.
IF(C_E(datoprima) C_E(datodopo) ) cycle
646 !compute time gradient only inside timeconfidence
648 call getval(td,asec=asec)
649 .and..or.
if ((c_e(qctem%timeconfidence) asec <= qctem%timeconfidence) &
650 .not.
c_e(qctem%timeconfidence)) then
651 if (c_e(datoprima)) gradprima=(datoqui-datoprima) / dble(asec)
655 call getval(td,asec=asec)
656 .and..or.
if ((c_e(qctem%timeconfidence) asec <= qctem%timeconfidence) &
657 .not.
c_e(qctem%timeconfidence)) then
658 if (c_e(datodopo)) graddopo =(datodopo-datoqui ) / dble(asec)
663 call l4f_log(L4F_DEBUG,"qctem gradprima:
"// to_char(gradprima)//" graddopo:
"//to_char(graddopo))
665 ! we need some gradient
666 .NOT..and..NOT.
IF(C_E(gradprima) C_E(graddopo) ) cycle
669 ! for gap we set negative gradient
670 ! for spike positive gradinet
671 .NOT.
IF(C_E(gradprima) ) then
673 ! set gap for other one
674 grad= sign(abs(graddopo),-1.d0)
676 .NOT.
else IF(C_E(graddopo) ) then
678 ! set gap for other one
679 grad= sign(abs(gradprima),-1.d0)
683 if (abs(max(abs(gradprima),abs(graddopo))-min(abs(gradprima),abs(graddopo))) < &
684 .and.
max(abs(gradprima),abs(graddopo))/2. (sign(1.d0,gradprima)*sign(1.d0,graddopo)) < 0.) then
686 grad= min(abs(gradprima),abs(graddopo))
689 grad= sign(max(abs(gradprima),abs(graddopo)),-1.d0)
693 if (qctem%operation == "gradient
") then
697 È
!ATTENZIONE TODO : inddativarr UNA GRANDE SEMPLIFICAZIONE NON VERA SE TIPI DI DATO DIVERSI !!!!
698 if (qctem%operation == "run
") then
700 ! choice which network we have to use
703 call l4f_log(L4F_DEBUG,"qctem choice gradient type: spike
")
705 indcnetwork=indcnetworks
708 call l4f_log(L4F_DEBUG,"qctem choice gradient type: gradmax
")
710 indcnetwork=indcnetworkg
714 call l4f_log(L4F_DEBUG,"gradiente da confrontare con qctem clima:
"//t2c(grad))
716 do indcana=1,size(qctem%clima%ana)
718 climaquii=(qctem%clima%voldatir(indcana &
719 ,indctime,indclevel,indctimerange,indcdativarr,indcnetwork)&
720 -tem_b(ind))/tem_a(ind) ! denormalize
722 climaquif=(qctem%clima%voldatir(min(indcana+1,size(qctem%clima%ana)) &
723 ,indctime,indclevel,indctimerange,indcdativarr,indcnetwork)&
724 -tem_b(ind))/tem_a(ind) ! denormalize
727 call l4f_log(L4F_DEBUG,"qctem clima start:
"//t2c(climaquii))
728 call l4f_log(L4F_DEBUG,"qctem clima end:
"//t2c(climaquif))
730 .and.
if ( c_e(climaquii) c_e(climaquif )) then
732 .and..or.
if ( (grad >= climaquii grad < climaquif) &
733 .and..or.
(indcana == 1 grad < climaquii) &
734 .and.
(indcana == size(qctem%clima%ana) grad >= climaquif) ) then
737 call l4f_log(L4F_DEBUG,"qctem confidence:
"// t2c(qctem%clima%voldatiattrb&
738 (indcana,indctime,indclevel,indctimerange,indcdativarr,indcnetwork,1)))
741 qctem%v7d%voldatiattrb( indana, indtime, indlevel, indtimerange, inddativarr, indnetwork, indbattrout)=&
742 qctem%clima%voldatiattrb(indcana,indctime,indclevel,indctimerange,indcdativarr,indcnetwork,1 )
744 if ( associated ( qctem%data_id_in)) then
746 call l4f_log (L4F_DEBUG,"id:
"//t2c(&
747 qctem%data_id_in(indana,indtime,indlevel,indtimerange,indnetwork)))
749 qctem%data_id_out(indana,indtime,indlevel,indtimerange,indnetwork)=&
750 qctem%data_id_in(indana,indtime,indlevel,indtimerange,indnetwork)
762 if (qctem%operation == "gradient
") then
768 !!$print*,"risultato
"
769 !!$print *,qcspa%v7d%voldatiattrb(:,:,:,:,:,:,indbattrout)
770 !!$print*,"fine risultato
"
774 end subroutine quacontem
780 !> \example v7d_qctem.F90
781 !! Sample program for module qctem
Emit log message for a category with specific priority.
Check for problems return 0 if all check passed print diagnostics with log4f.
Utilities for CHARACTER variables.
Classi per la gestione delle coordinate temporali.
Utilities for managing files.
classe per la gestione del logging
Utilities and defines for quality control.
Controllo di qualità climatico.
Controllo di qualità temporale.
Classe per la gestione di un volume completo di dati osservati.
classe per import ed export di volumi da e in DB-All.e
Oggetto principale per il controllo di qualità