libsim Versione 7.2.4

◆ vol7d_normalize_data()

subroutine, public vol7d_normalize_data ( type(qcclitype), intent(inout) qccli)

Modulo 1: Calcolo dei parametri di normalizzazione dei dati I parametri di normalizzazione sono il 25°, il 50° e il 75° percentile (p25,p50,p75) oppure (p15.87,p50.,p84.13) Tali parametri verranno calcolati per ogni mese, per ogni ora, per ogni area.

Modulo 2: Normalizzazione dati Ciascun dato D verrà normalizzato come segue: DN = (D-p50)*2/(p75-p25) oppure DN = (D-p50)*2/(p16-p84) dove DN è il valore normalizzato. La scelta dei parametri di normalizzazione dipende dal mese, dall'ora, dall'area.

Parametri
[in,out]qcclivolume providing data to be computed, it is modified by the method

Definizione alla linea 785 del file modqccli.F90.

786 indcana=index(qccli%extreme%ana,ana)
787 call l4f_log(l4f_debug,"normalize data Index25:"//to_char(k)//to_char(indcana)//&
788 to_char(indctime)//to_char(indclevel)//&
789 to_char(indctimerange)//to_char(indcdativarr)//to_char(indcnetwork))
790 if (indcana > 0 )then
791 perc25=qccli%extreme%voldatir(indcana,indctime,indclevel,&
792 indctimerange,indcdativarr,indcnetwork)
793 end if
794
795 desc=50
796 write(ident,'("#",i2.2,2i3.3)')k,iarea,desc ! macro-area e descrittore
797 call init(ana,lat=0.d0,lon=0.d0,ident=ident)
798 indcana=index(qccli%extreme%ana,ana)
799 call l4f_log(l4f_debug,"normalize data Index50:"//to_char(k)//to_char(indcana)//&
800 to_char(indctime)//to_char(indclevel)//&
801 to_char(indctimerange)//to_char(indcdativarr)//to_char(indcnetwork))
802 if (indcana > 0 )then
803 perc50=qccli%extreme%voldatir(indcana,indctime,indclevel,&
804 indctimerange,indcdativarr,indcnetwork)
805 end if
806
807 desc=75
808 write(ident,'("#",i2.2,2i3.3)')k,iarea,desc ! macro-area e descrittore
809 call init(ana,lat=0.d0,lon=0.d0,ident=ident)
810 indcana=index(qccli%extreme%ana,ana)
811 call l4f_log(l4f_debug,"normalize data Index75:"//to_char(k)//to_char(indcana)//&
812 to_char(indctime)//to_char(indclevel)//&
813 to_char(indctimerange)//to_char(indcdativarr)//to_char(indcnetwork))
814 if (indcana > 0 )then
815 perc75=qccli%extreme%voldatir(indcana,indctime,indclevel,&
816 indctimerange,indcdativarr,indcnetwork)
817 end if
818
819 if ( c_e(perc25) .and. c_e(perc50) .and. c_e(perc75) .and. &
820 abs( perc75 - perc25 ) >= spacing( max(abs(perc75),abs(perc25))))then
821 ! normalize
822
823 call l4f_log(l4f_debug,"normalize data dato in : "//t2c(datoqui))
824 datoqui = (datoqui - perc50) / (perc75 - perc25) + &
825 base_value(qccli%v7d%dativar%r(inddativarr)%btable)
826 call l4f_log(l4f_debug,"normalize data dato out : "//t2c(datoqui))
827
828 else
829
830 datoqui=rmiss
831
832 endif
833
834 qccli%v7d%voldatir (indana ,indtime ,indlevel ,indtimerange ,&
835 inddativarr, indnetwork ) = datoqui
836
837 end do
838 end do
839 end do
840 end do
841 end do
842
843 read (qccli%v7d%ana(indana)%ident,'(a1,i2.2,2i3.3)') mycanc, clev, iarea, desc
844 if (mycanc == canc) then
845 clev=0
846 iarea=0
847 write (qccli%v7d%ana(indana)%ident,'(a1,i2.2,2i3.3)') canc, clev, iarea, desc
848 end if
849
850end do
851
852end SUBROUTINE vol7d_normalize_data
853
854
855real function base_value(btable)
856character (len=10) ,intent(in):: btable
857
858character (len=10) :: btables(1) =(/"B12101"/)
859real :: base_values(1) =(/273.15 /)
860integer :: ind
861
862ind = index_c(btables,btable)
863
864if (ind > 0) then
865 base_value = base_values(ind)
866else
867 call l4f_log(l4f_warn,"modqccli_base_value: variable "//btable//" do not have base value")
868 base_value = 0.
869end if
870
871return
872
873end function base_value
874
875
882SUBROUTINE quaconcli (qccli,battrinv,battrout,&
883 anamask,timemask,levelmask,timerangemask,varmask,networkmask)
884
885
886type(qcclitype),intent(in out) :: qccli
887character (len=10) ,intent(in),optional :: battrinv !< attributo invalidated in input/output
888character (len=10) ,intent(in),optional :: battrout !< attributo con la confidenza climatologica in output
889logical ,intent(in),optional :: anamask(:)
890logical ,intent(in),optional :: timemask(:)
891logical ,intent(in),optional :: levelmask(:)
892logical ,intent(in),optional :: timerangemask(:)
893logical ,intent(in),optional :: varmask(:)
894logical ,intent(in),optional :: networkmask(:)
895
896CHARACTER(len=vol7d_ana_lenident) :: ident
897REAL(kind=fp_geo) :: latc,lonc
898integer :: mese, ora
899 !local
900integer :: indbattrinv,indbattrout
901logical :: anamaskl(size(qccli%v7d%ana)), timemaskl(size(qccli%v7d%time)), levelmaskl(size(qccli%v7d%level)), &
902 timerangemaskl(size(qccli%v7d%timerange)), varmaskl(size(qccli%v7d%dativar%r)), networkmaskl(size(qccli%v7d%network))
903
904integer :: indana , indtime ,indlevel ,indtimerange ,inddativarr, indnetwork
905integer :: indcana, indctime,indclevel,indctimerange,indcdativarr,indcnetwork
906real :: datoqui,climaquii,climaquif, extremequii,extremequif,perc25,perc50,perc75
907integer :: iarea,desc,indn,indvar,k
908real :: height
909!integer, allocatable :: indcanav(:)
910
911
912TYPE(vol7d_network) :: network
913TYPE(vol7d_ana) :: ana
914TYPE(datetime) :: time, nintime
915TYPE(vol7d_level):: level
916type(vol7d_var) :: anavar
917integer :: iclv(size(qccli%v7d%ana))
918character(len=1) :: type
919
920
921!call qccli_validate (qccli)
922
923indbattrinv=0
924if (associated(qccli%v7d%datiattr%b))then
925 if (present(battrinv))then
926 indbattrinv = index_c(qccli%v7d%datiattr%b(:)%btable, battrinv)
927 else
928 indbattrinv = index_c(qccli%v7d%datiattr%b(:)%btable, qcattrvarsbtables(1))
929 end if
930end if
931
932if ( indbattrinv <= 0 ) then
933
934 call l4f_category_log(qccli%category,l4f_error,"error finding attribute index in/out "//qcattrvarsbtables(1))
935 call raise_error("error finding attribute index in/out "//qcattrvarsbtables(1))
936
937end if
938
939
940if (present(battrout))then
941 indbattrout = index_c(qccli%v7d%datiattr%b(:)%btable, battrout)
942else
943 indbattrout = index_c(qccli%v7d%datiattr%b(:)%btable, qcattrvarsbtables(2))
944end if
945
946if ( indbattrout <= 0 ) then
947
948 call l4f_category_log(qccli%category,l4f_error,"error finding attribute index in/out "//qcattrvarsbtables(2))
949 call raise_error("error finding attribute index in/out "//qcattrvarsbtables(2))
950
951end if
952
953if(present(anamask)) then
954 anamaskl = anamask
955else
956 anamaskl = .true.
957endif
958if(present(timemask)) then
959 timemaskl = timemask
960else
961 timemaskl = .true.
962endif
963if(present(levelmask)) then
964 levelmaskl = levelmask
965else
966 levelmaskl = .true.
967endif
968if(present(timerangemask)) then
969 timerangemaskl = timerangemask
970else
971 timerangemaskl = .true.
972endif
973if(present(varmask)) then
974 varmaskl = varmask
975else
976 varmaskl = .true.
977endif
978if(present(networkmask)) then
979 networkmaskl = networkmask
980else
981 networkmaskl = .true.
982endif
983
984qccli%v7d%voldatiattrb(:,:,:,:,:,:,indbattrout)=ibmiss
985call init(anavar,"B07030" )
986type=cmiss
987indvar = index(qccli%v7d%anavar, anavar, type=type)
988
989do indana=1,size(qccli%v7d%ana)
990
991! if (qccli%height2level) then
992! iarea= supermacroa(qccli%in_macroa(indana))
993! else
994 iarea= qccli%in_macroa(indana)
995
996 if (.not. c_e(iarea)) cycle
997
998! end if
999! write(ident,'("BOX",2i3.3)')iarea,desc ! macro-area e descrittore
1000 !lat=0.0d0
1001 !lon=0.0d0
1002 !write(ident,'("BOX-",2i2.2)')iarea,lperc ! macro-area e percentile
1003 !call init(ana,lat=lat,lon=lon,ident=ident)
1004
1005 !allocate (indcanav(count(match(qccli%clima%ana(:)%ident,ident))))
1006 !indcanav=match(qccli%clima%ana(:)%ident,ident))))
1007
1008 latc=0.d0
1009 lonc=0.d0
1010
1011
1012 ! use conventional level starting from station height
1013 if (qccli%height2level) then
1014
1015 height=rmiss
1016
1017 ! here we take the height fron any network (the first network win)
1018 do indn=1,size(qccli%v7d%network)
1019
1020 if( indvar > 0 .and. indn > 0 ) then
1021 select case (type)
1022 case("d")
1023 height=realdat(qccli%v7d%volanad(indana,indvar,indn),qccli%v7d%anavar%d(indvar))
1024 case("r")
1025 height=realdat(qccli%v7d%volanar(indana,indvar,indn),qccli%v7d%anavar%r(indvar))
1026 case ("i")
1027 height=realdat(qccli%v7d%volanai(indana,indvar,indn),qccli%v7d%anavar%i(indvar))
1028 case("b")
1029 height=realdat(qccli%v7d%volanab(indana,indvar,indn),qccli%v7d%anavar%b(indvar))
1030 case("c")
1031 height=realdat(qccli%v7d%volanac(indana,indvar,indn),qccli%v7d%anavar%c(indvar))
1032 case default
1033 height=rmiss
1034 end select
1035 end if
1036
1037 if (c_e(height)) exit
1038 end do
1039
Index method.

Generated with Doxygen.