|
◆ qc_compute_percentile()
subroutine, public modqccli::qc_compute_percentile |
( |
type(qcclitype), intent(inout) |
this, |
|
|
real, dimension(:), intent(in) |
perc_vals, |
|
|
type(cyclicdatetime), intent(in) |
cyclicdt, |
|
|
real, optional |
presentperc, |
|
|
integer, optional |
presentnumb |
|
) |
| |
- Parametri
-
[in,out] | this | volume providing data to be computed and output data. Input data are not modified by the method, apart from performing a vol7d_alloc_vol on it |
[in] | perc_vals | percentile values to use in compute, between 0. and 100. |
[in] | cyclicdt | cyclic date and time |
| presentperc | percentual of data present for compute (default=not used) |
| presentnumb | number of data present for compute (default=not used) $logical, optional :: height2level !< use conventional level starting from station height |
Definizione alla linea 1618 del file modqccli.F90.
1620 do inddativarr=1, size(this%v7d%dativar%r)
1621 do indtimerange=1, size(this%v7d%timerange)
1622 do indlevel=1, size(this%v7d%level)
1628 mask = spread(spread((this%v7d%time == cyclicdt ),1, size(this%v7d%ana)),3, size(this%v7d%network))
1631 CALL l4f_category_log(this%category, l4f_debug, 'count has value '//t2c(count &
1632 (mask .and. c_e(this%v7d%voldatir(:,:, indlevel, indtimerange, inddativarr,:)))))
1636 do j=1, size(mask,1)
1637 if (areav(j) /= area(i)) mask(j,:,:) =.false.
1640 if (this%height2level) then
1641 allocate (maskplus( size(this%v7d%ana), size(this%v7d%time), size(this%v7d%network)))
1644 CALL l4f_category_log(this%category, l4f_debug, 'k has value '//t2c(k))
1647 do iana=1, size(mask,1)
1648 if (iclv(iana) /= k) maskplus(iana,:,:) =.false.
1649 if (iclv(iana) == k) maskplus(iana,:,:) = mask(iana,:,:)
1652 call sub_perc(maskplus)
1655 deallocate(maskplus)
1667 deallocate (perc,mask,area)
1671 subroutine sub_perc(mymask)
1673 logical :: mymask(:,:,:)
1677 ( c_e(lpresentperc) .and. ((float(count &
1678 (mymask .and. c_e(this%v7d%voldatir(:,:, indlevel, indtimerange, inddativarr,:)))&
1680 float(count(mymask))) < lpresentperc)) &
1684 ( c_e(lpresentnumb) .and. (count &
1685 (mymask .and. c_e(this%v7d%voldatir(:,:, indlevel, indtimerange, inddativarr,:))) < lpresentnumb)&
1690 perc= stat_percentile(&
1691 pack(this%v7d%voldatir(:,:, indlevel, indtimerange, inddativarr,:), &
1696 do j=1, size(perc_vals)
1697 indana=(k-1)* size(perc_vals)*narea + (j-1)*narea + i
1698 this%extreme%voldatir(indana, indtime, indlevel, indtimerange, inddativarr, indnetwork)=&
1703 end subroutine sub_perc
1706 end SUBROUTINE qc_compute_percentile
1709 SUBROUTINE qc_compute_normalizeddensityindex(this, perc_vals,cyclicdt,presentperc, presentnumb,data_normalized)
1711 TYPE(qcclitype), INTENT(inout) :: this
1715 real, intent(in) :: perc_vals(:)
1716 TYPE(cyclicdatetime), INTENT(in) :: cyclicdt
1717 real, optional :: presentperc
1718 integer, optional :: presentnumb
1719 logical, optional, intent(in) :: data_normalized
1721 integer :: indana,indtime,indvar,indnetwork,indlevel ,indtimerange ,inddativarr, indattr
1722 integer :: i,j,k,iana,narea
1724 TYPE(vol7d_var) :: var
1725 character(len=vol7d_ana_lenident) :: ident
1726 character(len=1) :: type
1727 integer :: areav(size(this%v7d%ana)),iclv(size(this%v7d%ana))
1728 logical, allocatable :: mask(:,:,:),maskplus(:,:,:), maskarea(:)
1729 integer, allocatable :: area(:)
1730 REAL, DIMENSION(:), allocatable :: ndi,limbins
1731 real :: lpresentperc
1732 integer :: lpresentnumb
1735 lnorm = optio_log(data_normalized)
1736 lpresentperc=optio_r(presentperc)
1737 lpresentnumb=optio_i(presentnumb)
1739 allocate (ndi( size(perc_vals)-1),limbins( size(perc_vals)))
1740 call delete(this%clima)
1741 CALL init(this%clima, time_definition=this%v7d%time_definition)
1742 call init(var, btable= "B01192")
1744 if (.NOT.(lnorm)) then
1747 indvar = index(this%v7d%anavar, var, type=type)
1753 areav=integerdat(this%v7d%volanad(:,indvar,indnetwork),this%v7d%anavar%d(indvar))
1755 areav=integerdat(this%v7d%volanar(:,indvar,indnetwork),this%v7d%anavar%r(indvar))
1757 areav=integerdat(this%v7d%volanai(:,indvar,indnetwork),this%v7d%anavar%i(indvar))
1759 areav=integerdat(this%v7d%volanab(:,indvar,indnetwork),this%v7d%anavar%b(indvar))
1761 areav=integerdat(this%v7d%volanac(:,indvar,indnetwork),this%v7d%anavar%c(indvar))
1767 allocate(maskarea( size(this%v7d%ana)))
1768 maskarea(:)= areav(:) /= imiss
1769 narea=count_distinct(areav,maskarea)
1770 allocate(area(narea))
1771 area=pack_distinct(areav,narea,maskarea)
1772 deallocate(maskarea)
1774 if (this%height2level) then
1775 call vol7d_alloc(this%clima,nana=narea*( size(perc_vals)-1)*cli_nlevel)
1777 call vol7d_alloc(this%clima,nana=narea*( size(perc_vals)-1))
1784 do j=1, size(perc_vals)-1
1785 if (this%height2level) then
1787 write(ident, '("#",i2.2,2i3.3)')k,area(i),nint(perc_vals(j))
1788 call init(this%clima%ana((k-1)*( size(perc_vals)-1)*narea + (j-1)*narea + i),ident=ident,lat=0d0,lon=0d0)
1792 write(ident, '("#",i2.2,2i3.3)')k,area(i),nint(perc_vals(j))
1793 call init(this%clima%ana((j-1)*narea+i),ident=ident,lat=0d0,lon=0d0)
1802 if (this%height2level) then
1804 call init(var, btable= "B07030")
1807 indvar = index(this%v7d%anavar, var, type=type)
1809 do k=1, size(this%v7d%ana)
1813 do indnetwork=1, size(this%v7d%network)
1815 if( indvar > 0 .and. indnetwork > 0 ) then
1818 height=realdat(this%v7d%volanad(k,indvar,indnetwork),this%v7d%anavar%d(indvar))
1820 height=realdat(this%v7d%volanar(k,indvar,indnetwork),this%v7d%anavar%r(indvar))
1822 height=realdat(this%v7d%volanai(k,indvar,indnetwork),this%v7d%anavar%i(indvar))
1824 height=realdat(this%v7d%volanab(k,indvar,indnetwork),this%v7d%anavar%b(indvar))
1826 height=realdat(this%v7d%volanac(k,indvar,indnetwork),this%v7d%anavar%c(indvar))
1832 if (c_e(height)) exit
1835 if (c_e(height)) then
1836 iclv(k)=firsttrue(cli_level1 <= height .and. height <= cli_level2 )
1842 CALL l4f_log(l4f_debug, 'qc_compute_NormalizedDensityIndex height has value '//t2c(height, "Missing"))
1843 CALL l4f_log(l4f_debug, 'for k having number '//t2c(k)//&
1844 ' iclv has value '//t2c(iclv(k)))
1855 call vol7d_alloc(this%clima,nana= size(perc_vals)-1)
1856 do j=1, size(perc_vals)-1
1857 write(ident, '("#",i2.2,2i3.3)')0,0,nint(perc_vals(j))
1858 call init(this%clima%ana(j),ident=ident,lat=0d0,lon=0d0)
1868 call vol7d_alloc(this%clima,nlevel= size(this%v7d%level), ntimerange= size(this%v7d%timerange), &
1869 ndativarr= size(this%v7d%dativar%r), nnetwork=1,ntime=1,ndativarattrr= size(this%v7d%dativar%r),ndatiattrr=1)
1871 this%clima%level=this%v7d%level
1872 this%clima%timerange=this%v7d%timerange
1873 this%clima%dativar%r=this%v7d%dativar%r
1874 this%clima%dativarattr%r=this%clima%dativar%r
1875 call init(this%clima%datiattr%r(1), btable= "*B33209")
1876 this%clima%time(1)=cyclicdatetime_to_conventional(cyclicdt)
1878 call l4f_category_log(this%category,l4f_info, "vol7d_compute_ndi conventional datetime "//to_char(this%clima%time(1)))
1879 call init(this%clima%network(1),name= "qcclima-ndi")
1881 call vol7d_alloc_vol(this%clima,inivol=.true.)
1883 allocate (mask( size(this%v7d%ana), size(this%v7d%time), size(this%v7d%network)))
1888 do inddativarr=1, size(this%v7d%dativar%r)
1889 do indtimerange=1, size(this%v7d%timerange)
1890 do indlevel=1, size(this%v7d%level)
1891 if (.NOT.(lnorm)) then
|