libsim Versione 7.2.4
modqctem.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
19! derived from a work of:
20!!$CC**********************************************************************CC
21!!$CC**********************************************************************CC
22!!$CC CC
23!!$CC SERVIZIO METEOROLOGICO REGIONE EMILA ROMAGNA CC
24!!$CC E.R.S.A. CC
25!!$CC CC
26!!$CC CC
27!!$CC PAOLO PATRUNO PIER PAOLO ALBERONI CC
28!!$CC CC
29!!$CC BOLOGNA 1992 CC
30!!$CC**********************************************************************CC
31!!$CC**********************************************************************CC
32
33#include "config.h"
34
66
67
68module modqctem
69
75use modqc
76use modqccli
77
78!use array_utilities
79!use io_units
80#ifdef HAVE_DBALLE
82#endif
83
84implicit none
85
86
87character (len=255),parameter:: subcategorytem="QCtem"
88
89integer, parameter :: tem_nvar=1
90CHARACTER(len=10) :: tem_btable(tem_nvar)=(/"B12101"/)
92real, parameter :: tem_a(tem_nvar) = (/1.e5/)
94real, parameter :: tem_b(tem_nvar) = (/250./)
95
97type :: qctemtype
98 type (vol7d),pointer :: v7d => null()
99 integer,pointer :: data_id_in(:,:,:,:,:) => null()
100 integer,pointer :: data_id_out(:,:,:,:,:) => null()
101 integer :: category
102 type (qcclitype) :: qccli
103 type (vol7d) :: clima
104 character(len=20):: operation
105 integer :: timeconfidence
106end type qctemtype
107
108
110interface init
111 module procedure qcteminit
112end interface
113
115interface alloc
116 module procedure qctemalloc
117end interface
118
120interface delete
121 module procedure qctemdelete
122end interface
123
124PRIVATE
125PUBLIC tem_nvar, tem_btable, tem_a, tem_b, qctemtype, init, alloc, delete, &
126 quacontem
127
128
129contains
130
133
134subroutine qcteminit(qctem,v7d,var, timei, timef, coordmin, coordmax, data_id_in,extremepath,temporalpath,&
135#ifdef HAVE_DBALLE
136 dsne,usere,passworde,&
137 dsntem,usertem,passwordtem,&
138#endif
139 height2level,operation,timeconfidence,categoryappend)
140
141type(qctemtype),intent(in out) :: qctem
142type (vol7d),intent(in),target:: v7d
143character(len=*),INTENT(in) :: var(:)
146TYPE(geo_coord),INTENT(inout),optional :: coordmin,coordmax
148TYPE(datetime),INTENT(in),optional :: timei, timef
149integer,intent(in),optional,target:: data_id_in(:,:,:,:,:)
150character(len=*),intent(in),optional :: extremepath
151character(len=*),intent(in),optional :: temporalpath
152logical ,intent(in),optional :: height2level
153character(len=*), optional :: operation
154integer,intent(in),optional :: timeconfidence
155character(len=*),INTENT(in),OPTIONAL :: categoryappend
156
157#ifdef HAVE_DBALLE
158type (vol7d_dballe) :: v7d_dballetmp
159character(len=*),intent(in),optional :: dsne
160character(len=*),intent(in),optional :: usere
161character(len=*),intent(in),optional :: passworde
162character(len=*),intent(in),optional :: dsntem
163character(len=*),intent(in),optional :: usertem
164character(len=*),intent(in),optional :: passwordtem
165character(len=512) :: ldsntem
166character(len=512) :: lusertem
167character(len=512) :: lpasswordtem
168#endif
169
170type (vol7d) :: v7dtmp
171TYPE(vol7d_network):: network
172integer :: iuni,i
173character(len=512) :: filepathtem
174character(len=512) :: a_name
175character(len=9) ::netname(2)=(/"qctemgndi","qctemsndi"/)
176
177
178call l4f_launcher(a_name,a_name_append=trim(subcategorytem)//"."//trim(categoryappend))
179qctem%category=l4f_category_get(a_name)
180
181nullify ( qctem%data_id_in )
182nullify ( qctem%data_id_out )
183
184! riporto il volume dati nel mio oggetto
185qctem%v7d => v7d
186
187qctem%operation=optio_c(operation,20)
188filepathtem=optio_c(temporalpath,512)
189
190!check options
191if (qctem%operation /= "gradient" .and. qctem%operation /= "run") then
192 call l4f_category_log(qctem%category,l4f_error,"operation is wrong: "//qctem%operation)
193 call raise_error()
194end if
195
196!inglobe id
197if (present(data_id_in))then
198 qctem%data_id_in => data_id_in
199end if
200
201qctem%timeconfidence = optio_i(timeconfidence)
202
203! load extreme
204call init(qctem%qccli,v7d,var, timei, timef, data_id_in,&
205 macropath=cmiss, climapath=cmiss, extremepath=extremepath, &
206#ifdef HAVE_DBALLE
207 dsncli=cmiss,dsnextreme=dsne,user=usere,password=passworde,&
208#endif
209 height2level=height2level,categoryappend=categoryappend)
210
211
212! now load temporal clima
213
214do i=1,size(netname)
215
216 if (qctem%operation == "run") then
217 call init(network,netname(i))
218
219#ifdef HAVE_DBALLE
220 call optio(dsntem,ldsntem)
221 call optio(usertem,lusertem)
222 call optio(passwordtem,lpasswordtem)
223
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")
226 call raise_error()
227 end if
228
229 if (.not. c_e(ldsntem)) then
230
231#endif
232
233 if (.not. c_e(filepathtem)) then
234 filepathtem=get_package_filepath(netname(i)//'.v7d', filetype_data)
235 end if
236
237 if (c_e(filepathtem))then
238
239 select case (trim(lowercase(suffixname(filepathtem))))
240
241 case("v7d")
242 iuni=getunit()
243 call import(v7dtmp,filename=filepathtem,unit=iuni)
244 close (unit=iuni)
245
246#ifdef HAVE_DBALLE
247 case("bufr")
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)
253#endif
254
255 case default
256 call l4f_category_log(qctem%category,l4f_error,&
257 "file type not supported (use .v7d or .bufr suffix only): "//trim(filepathtem))
258 call raise_error()
259 end select
260
261 else
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()
265 end if
266
267#ifdef HAVE_DBALLE
268 else
269
270 call l4f_category_log(qctem%category,l4f_debug,"init v7d_dballetem")
271 call init(v7d_dballetmp,dsn=ldsntem,user=lusertem,password=lpasswordtem,write=.false.,&
272 file=.false.,categoryappend=trim(a_name)//".tem")
273 call l4f_category_log(qctem%category,l4f_debug,"import v7d_dballetmp")
274 call import(v7d_dballetmp,var=var, &
275 varkind=(/("r",i=1,size(var))/),attr=(/"*B33209"/),attrkind=(/"b"/),network=network)
276
277 call copy(v7d_dballetmp%vol7d,v7dtmp)
279 call delete(v7d_dballetmp)
281 end if
282#endif
283 call vol7d_merge(qctem%clima,v7dtmp)
284
285 end if
287end do
289return
290end subroutine qcteminit
295subroutine qctemalloc(qctem)
296 ! pseudo costruttore con distruttore automatico
297
298type(qctemtype),intent(in out) :: qctem
299
300integer :: istatt
301integer :: sh(5)
302
303! se ti sei dimenticato di deallocare ci penso io
304call qctemdealloc(qctem)
305
306if (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)
309 if (istatt /= 0)then
310 call l4f_category_log(qctem%category,l4f_error,"allocate error")
311 call raise_error("allocate error")
312 else
313 qctem%data_id_out=imiss
314 end if
315end if
316
317end subroutine qctemalloc
318
319
321subroutine qctemdealloc(qctem)
322 ! pseudo distruttore
323
324type(qctemtype),intent(in out) :: qctem
325
326if (associated(qctem%data_id_out)) deallocate (qctem%data_id_out)
327
328if (associated(qctem%data_id_out)) then
329 deallocate (qctem%data_id_out)
330 nullify (qctem%data_id_out)
331end if
332
333end subroutine qctemdealloc
334
335
337
338
339subroutine qctemdelete(qctem)
340 ! decostruttore a mezzo
341type(qctemtype),intent(in out) :: qctem
342
343call qctemdealloc(qctem)
344
345call delete(qctem%qccli)
346
347!delete logger
348call l4f_category_delete(qctem%category)
349
350return
351end subroutine qctemdelete
352
353
356
357SUBROUTINE quacontem (qctem,battrinv,battrcli,battrout,&
358 anamask,timemask,levelmask,timerangemask,varmask,networkmask)
359
360
361type(qctemtype),intent(in out) :: qctem
362character (len=10) ,intent(in),optional :: battrinv
363character (len=10) ,intent(in),optional :: battrcli
364character (len=10) ,intent(in),optional :: battrout
365logical ,intent(in),optional :: anamask(:)
366logical ,intent(in),optional :: timemask(:)
367logical ,intent(in),optional :: levelmask(:)
368logical ,intent(in),optional :: timerangemask(:)
369logical ,intent(in),optional :: varmask(:)
370logical ,intent(in),optional :: networkmask(:)
371
372 !REAL(kind=fp_geo) :: lat,lon
373integer :: asec
374 !local
375integer :: indbattrinv,indbattrcli,indbattrout,grunit
376logical :: 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))
378
379integer :: indana , indtime ,indlevel ,indtimerange ,inddativarr, indnetwork,indtimenear
380integer :: indcana , indctime ,indclevel ,indctimerange ,indcdativarr, indcnetwork, indcnetworks, indcnetworkg
381real :: datoqui,datoprima,datodopo,climaquii, climaquif
382 !integer, allocatable :: indcanav(:)
383
384 !TYPE(vol7d_ana) :: ana
385TYPE(datetime) :: time,prima, ora, dopo
386TYPE(vol7d_network):: network
387type(timedelta) :: td
388
389double precision :: gradprima,graddopo,grad
390 !call qctem_validate (qctem)
391character(len=512) :: filename
392logical :: exist
393integer :: ind
394
395!localize optional parameter
396if (present(battrinv))then
397 indbattrinv = index_c(qctem%v7d%datiattr%b(:)%btable, battrinv)
398else
399 indbattrinv = index_c(qctem%v7d%datiattr%b(:)%btable, qcattrvarsbtables(1))
400end if
401
402if (present(battrcli))then
403 indbattrcli = index_c(qctem%v7d%datiattr%b(:)%btable, battrcli)
404else
405 indbattrcli = index_c(qctem%v7d%datiattr%b(:)%btable, qcattrvarsbtables(2))
406end if
407
408if (present(battrout))then
409 indbattrout = index_c(qctem%v7d%datiattr%b(:)%btable, battrout)
410else
411 indbattrout = index_c(qctem%v7d%datiattr%b(:)%btable, qcattrvarsbtables(3))
412end if
413
414
415! some checks on input
416!if (indbattrinv <=0 .or. indbattrcli <= 0 .or. indbattrout <= 0 ) then
417if (indbattrout <= 0 ) then
418
419 call l4f_category_log(qctem%category,l4f_error,"error finding attribute index for output")
420 call raise_error()
421
422end if
423
424if (qctem%operation == "gradient") then
425
426 !check for gradient operation
427 if ( size(qctem%v7d%level) > 1 .or.&
428 size(qctem%v7d%timerange) > 1 .or.&
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")
431 call raise_error()
432 end if
433
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")
437 return
438 end if
439end if
440
441 ! set other local variable from optional parameter
442if(present(anamask)) then
443 anamaskl = anamask
444else
445 anamaskl = .true.
446endif
447if(present(timemask)) then
448 timemaskl = timemask
449else
450 timemaskl = .true.
451endif
452if(present(levelmask)) then
453 levelmaskl = levelmask
454else
455 levelmaskl = .true.
456endif
457if(present(timerangemask)) then
458 timerangemaskl = timerangemask
459else
460 timerangemaskl = .true.
461endif
462if(present(varmask)) then
463 varmaskl = varmask
464else
465 varmaskl = .true.
466endif
467if(present(networkmask)) then
468 networkmaskl = networkmask
469else
470 networkmaskl = .true.
471endif
472
473! do not touch data that will not in/validated by QC
474qctem%v7d%voldatiattrb(:,:,:,:,:,:,indbattrout)=ibmiss
475
476! normalize data in space and time
477call vol7d_normalize_data(qctem%qccli)
478
479
480! compute some index for temporal clima
481!! compute the conventional generic datetime
482!!cyclicdt = cyclicdatetime_new(chardate="/////////") !TMMGGhhmm
483time=cyclicdatetime_to_conventional(cyclicdatetime_new(chardate="/////////")) !TMMGGhhmm
484!!call init(time, year=1007, month=1, day=1, hour=01, minute=01)
485
486
487if (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)
494end if
495
496do indana=1,size(qctem%v7d%ana)
497 if (.not.anamaskl(indana)) cycle
498 call l4f_category_log(qctem%category,l4f_info,&
499 "Check ana:"//to_char(qctem%v7d%ana(indana)) )
500
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)
506
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))
512
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)//&
516 ".grad"
517
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)
520
521 inquire(file=filename, exist=exist)
522
523 grunit=getunit()
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')
527 end if
528 ! say we have to write header in file
529 if (.not. exist) then
530 call l4f_category_log(qctem%category,l4f_info,"write header in gradient file")
531 write (grunit,*) &
532 qctem%v7d%level(indlevel), &
533 qctem%v7d%timerange(indtimerange), &
534 qctem%v7d%dativar%r(inddativarr)
535 end if
536 end if
537
538!!$ call l4f_log(L4F_INFO,"Index:"// t2c(indana)//t2c(indnetwork)//t2c(indlevel)//&
539!!$ t2c(indtimerange)//t2c(inddativarr)//t2c(indtime))
540
541
542 do indtime=2,size(qctem%v7d%time)-1
543
544 if (.not.timemaskl(indtime).or. .not. levelmaskl(indlevel).or. &
545 .not. timerangemaskl(indtimerange) .or. .not. varmaskl(inddativarr) .or. .not. networkmaskl(indnetwork)) cycle
546
547
548 datoqui = qctem%v7d%voldatir (indana ,indtime ,indlevel ,indtimerange ,inddativarr, indnetwork )
549 if (.not. c_e(datoqui)) cycle
550 ora = qctem%v7d%time (indtime)
551
552 ! invalidated
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")
558 cycle
559 end if
560 end if
561
562 ! gross error check
563 if (indbattrcli > 0) then
564 if( .not. 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")
568 cycle
569 end if
570 end if
571
572
573
574 if (qctem%operation == "run") then
575
576 indclevel = index(qctem%clima%level , qctem%v7d%level(indlevel))
577 indctimerange = index(qctem%clima%timerange , qctem%v7d%timerange(indtimerange))
578
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))
582
583#ifdef DEBUG
584 call l4f_log(l4f_debug,"QCtem Index:"// to_char(indctime)//to_char(indclevel)//&
585 to_char(indctimerange)//to_char(indcdativarr)//to_char(indcnetworks))
586#endif
587 if ( indctime <= 0 .or. indclevel <= 0 .or. indctimerange <= 0 .or. indcdativarr <= 0 &
588 .or. indcnetworks <= 0 ) cycle
589 end if
590
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)
594!!$
595
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)
600
601 ! invalidated
602 if (indbattrinv > 0) then
603 if( invalidated(qctem%v7d%voldatiattrb&
604 (indana,indtimenear,indlevel,indtimerange,inddativarr,indnetwork,indbattrinv))) then
605 datoprima=rmiss
606 end if
607 end if
608
609 ! gross error check
610 if (indbattrcli > 0) then
611 if( .not. vdge(qctem%v7d%voldatiattrb&
612 (indana,indtimenear,indlevel,indtimerange,inddativarr,indnetwork,indbattrcli))) then
613 datoprima=rmiss
614 end if
615 end if
616
617
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)
622
623 ! invalidated
624 if (indbattrinv > 0) then
625 if( invalidated(qctem%v7d%voldatiattrb&
626 (indana,indtimenear,indlevel,indtimerange,inddativarr,indnetwork,indbattrinv))) then
627 datodopo=rmiss
628 end if
629 end if
630
631 ! gross error check
632 if (indbattrcli > 0) then
633 if( .not. vdge(qctem%v7d%voldatiattrb&
634 (indana,indtimenear,indlevel,indtimerange,inddativarr,indnetwork,indbattrcli))) then
635 datodopo=rmiss
636 end if
637 end if
638
639
640 IF(.NOT.c_e(datoprima) .and. .NOT.c_e(datodopo) ) cycle
641
642 gradprima=dmiss
643 graddopo=dmiss
644 grad=dmiss
645
646 !compute time gradient only inside timeconfidence
647 td=ora-prima
648 call getval(td,asec=asec)
649 if ((c_e(qctem%timeconfidence) .and. asec <= qctem%timeconfidence) .or. &
650 .not. c_e(qctem%timeconfidence)) then
651 if (c_e(datoprima)) gradprima=(datoqui-datoprima) / dble(asec)
652 end if
653
654 td=dopo-ora
655 call getval(td,asec=asec)
656 if ((c_e(qctem%timeconfidence) .and. asec <= qctem%timeconfidence) .or. &
657 .not. c_e(qctem%timeconfidence)) then
658 if (c_e(datodopo)) graddopo =(datodopo-datoqui ) / dble(asec)
659 end if
660
661
662#ifdef DEBUG
663 call l4f_log(l4f_debug,"QCtem gradprima:"// to_char(gradprima)//" graddopo:"//to_char(graddopo))
664#endif
665 ! we need some gradient
666 IF(.NOT.c_e(gradprima) .and. .NOT.c_e(graddopo) ) cycle
667
668
669 ! for gap we set negative gradient
670 ! for spike positive gradinet
671 IF(.NOT.c_e(gradprima) ) then
672
673 ! set gap for other one
674 grad= sign(abs(graddopo),-1.d0)
675
676 else IF(.NOT.c_e(graddopo) ) then
677
678 ! set gap for other one
679 grad= sign(abs(gradprima),-1.d0)
680
681 else
682
683 if (abs(max(abs(gradprima),abs(graddopo))-min(abs(gradprima),abs(graddopo))) < &
684 max(abs(gradprima),abs(graddopo))/2. .and. (sign(1.d0,gradprima)*sign(1.d0,graddopo)) < 0.) then
685 ! spike
686 grad= min(abs(gradprima),abs(graddopo))
687 else
688 ! gap
689 grad= sign(max(abs(gradprima),abs(graddopo)),-1.d0)
690 end if
691 end IF
692
693 if (qctem%operation == "gradient") then
694 write(grunit,*)grad
695 end if
696
697 !ATTENZIONE TODO : inddativarr È UNA GRANDE SEMPLIFICAZIONE NON VERA SE TIPI DI DATO DIVERSI !!!!
698 if (qctem%operation == "run") then
699
700 ! choice which network we have to use
701 if (grad >= 0) then
702#ifdef DEBUG
703 call l4f_log(l4f_debug,"QCtem choice gradient type: spike")
704#endif
705 indcnetwork=indcnetworks
706 else
707#ifdef DEBUG
708 call l4f_log(l4f_debug,"QCtem choice gradient type: gradmax")
709#endif
710 indcnetwork=indcnetworkg
711 end if
712
713 grad=abs(grad)
714 call l4f_log(l4f_debug,"gradiente da confrontare con QCtem clima:"//t2c(grad))
715
716 do indcana=1,size(qctem%clima%ana)
717
718 climaquii=(qctem%clima%voldatir(indcana &
719 ,indctime,indclevel,indctimerange,indcdativarr,indcnetwork)&
720 -tem_b(ind))/tem_a(ind) ! denormalize
721
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
725
726#ifdef DEBUG
727 call l4f_log(l4f_debug,"QCtem clima start:"//t2c(climaquii))
728 call l4f_log(l4f_debug,"QCtem clima end:"//t2c(climaquif))
729#endif
730 if ( c_e(climaquii) .and. c_e(climaquif )) then
731
732 if ( (grad >= climaquii .and. grad < climaquif) .or. &
733 (indcana == 1 .and. grad < climaquii) .or. &
734 (indcana == size(qctem%clima%ana) .and. grad >= climaquif) ) then
735
736#ifdef DEBUG
737 call l4f_log(l4f_debug,"QCtem confidence:"// t2c(qctem%clima%voldatiattrb&
738 (indcana,indctime,indclevel,indctimerange,indcdativarr,indcnetwork,1)))
739#endif
740
741 qctem%v7d%voldatiattrb( indana, indtime, indlevel, indtimerange, inddativarr, indnetwork, indbattrout)=&
742 qctem%clima%voldatiattrb(indcana,indctime,indclevel,indctimerange,indcdativarr,indcnetwork,1 )
743
744 if ( associated ( qctem%data_id_in)) then
745#ifdef DEBUG
746 call l4f_log (l4f_debug,"id: "//t2c(&
747 qctem%data_id_in(indana,indtime,indlevel,indtimerange,indnetwork)))
748#endif
749 qctem%data_id_out(indana,indtime,indlevel,indtimerange,indnetwork)=&
750 qctem%data_id_in(indana,indtime,indlevel,indtimerange,indnetwork)
751 end if
752 end if
753 end if
754 end do
755 end if
756 end do
757 end do
758 end do
759 end do
760 end do
761
762 if (qctem%operation == "gradient") then
763 close (unit=grunit)
764 end if
765
766end do
767
768!!$print*,"risultato"
769!!$print *,qcspa%v7d%voldatiattrb(:,:,:,:,:,:,indbattrout)
770!!$print*,"fine risultato"
771
772return
773
774end subroutine quacontem
775
776
777end module modqctem
778
779
Set of functions that return a trimmed CHARACTER representation of the input variable.
Set of functions that return a CHARACTER representation of the input variable.
Operatore di valore assoluto di un intervallo.
Restituiscono il valore dell'oggetto nella forma desiderata.
Index method.
Emit log message for a category with specific priority.
Test di dato invalidato.
Definition modqc.F90:411
Check data validity based on gross error check.
Definition modqc.F90:406
Allocazione di memoria.
Definition modqctem.F90:303
Cancellazione.
Definition modqctem.F90:308
Inizializzazione.
Definition modqctem.F90:298
Lettura da file.
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.
Definition modqc.F90:357
Controllo di qualità climatico.
Definition modqccli.F90:295
Controllo di qualità temporale.
Definition modqctem.F90:256
Classe per la gestione di un volume completo di dati osservati.
classe per import ed export di volumi da e in DB-All.e
Class for expressing an absolute time value.
Class for expressing a relative time interval.
Oggetto principale per il controllo di qualità
Definition modqccli.F90:322
Oggetto principale per il controllo di qualità
Definition modqctem.F90:285
Definisce un oggetto contenente i volumi anagrafica e dati e tutti i descrittori delle loro dimension...

Generated with Doxygen.