*file member=anxx_search_lowocc_cut library=snoman language=fortran77
*file date=08:Nov:2004
<<< subroutine anxx_search_lowocc_cut
* ANX: ANXX Banks Processor: Search for low occ cuts for each period
* Contact: M. Huang, Texas
* Parameters:-
* ==========
* Revision History:-
* ================
* 5.01 M. Huang First version.
implicit none
include 'zunit.inc'
include 'anxx_count_com.inc'
* Arguments Declarations:-
* ======================
* Local Variable Declarations:-
* ============================
integer num_valley, num_bin_for_peak1, num_bin_one_peak
integer ipmt, ihd, ibin, ibin_zero, iperiod, tot_num,
+ bin_peak1, bin1, bin2, min_bin, num_period,
+ bin_min_count_0
parameter (num_valley = 20, num_bin_for_peak1 = 7,
+ num_bin_one_peak = num_bin_for_peak1-2)
parameter (bin1 = 1, bin2 = 60)
real low_occ_stat, bin_size_low_occ,
+ min_rate_valley, count_peak1, min_count, min_count_0
integer num_allow_diff, tot_num_diff
real count_allow_diff, count_diff
parameter (num_allow_diff = 20, count_allow_diff = 10.0)
logical warn_short
* --- general initialization --------------------------------------------
warn_short = .false.
bin_size_low_occ = (frac2 - frac1)/float(num_bin_low_occ)
ibin_zero = int((0.0 - frac1)/bin_size_low_occ) + 1
* --- compute num_period; if too short, don't deal with low occ -------
num_period = int(elapse_time/time_interval)
if (num_period .le. 0) then
if (.not.warn_short) then
write(iqprnt,90004)
write(iqlog, 90004)
warn_short = .true.
endif
else
if (num_period .gt. nperiod) num_period = nperiod
do iperiod = 1, nperiod+1
low_occ_rate_cut(iperiod) = -9999.0
do ibin = 1, max_num_bin
low_occ_rate_count(ibin,iperiod) = 0.0
end do
end do
* --- count low occ hits and fill into different period --------------
bin_size_low_occ = (frac2 - frac1)/float(num_bin_low_occ)
do ipmt = 0, npmt
if (occ_count(ipmt) .gt. 0.0) then
low_occ_stat = occ_count(ipmt)/float(num_period)
ibin = (low_occ_stat-frac1)/bin_size_low_occ + 1
if (ibin.gt.num_bin_low_occ) ibin = num_bin_low_occ
low_occ_rate_count(ibin,num_period+1) =
+ low_occ_rate_count(ibin,num_period+1) + 1.0
ihd = 3000 + nperiod + 1
call hf1 (ihd,low_occ_stat,1.0)
do iperiod = 1, num_period
low_occ_stat = low_occ_count(ipmt,iperiod)
ibin = (low_occ_stat-frac1)/bin_size_low_occ + 1
if (ibin.gt.num_bin_low_occ) ibin = num_bin_low_occ
low_occ_rate_count(ibin,iperiod) =
+ low_occ_rate_count(ibin,iperiod) + 1.0
ihd = 3000 + iperiod
call hf1 (ihd,low_occ_stat,1.0)
end do
endif
end do ! end of 1st ipmt loop
* --- look for criterion cuts for each period -------------------------
do iperiod = 1, num_period+1
* * ... look for where is the min_count & min_bin .......... *
min_count = 1.0e30
do ibin = bin1, bin2
if ( low_occ_rate_count(ibin,iperiod) .le.
+ min_count) then
min_count = low_occ_rate_count(ibin,iperiod)
min_bin = ibin
endif
end do
* * ... check if there is any bins after min_bin ....... *
tot_num_diff = 0
do ibin = min_bin+1, num_bin_low_occ
count_diff = abs(min_count -
+ low_occ_rate_count(ibin,iperiod) )
if (count_diff .ge. count_allow_diff)
+ tot_num_diff = tot_num_diff + 1
end do
* * ... for cases that have valley not near zero ....... *
if (tot_num_diff .ge. num_allow_diff) then
min_rate_valley = frac1 +
+ ( float(min_bin) - 0.5)*bin_size_low_occ
else
tot_num = 0
bin_peak1 = ibin_zero
count_peak1 = low_occ_rate_count(ibin_zero,iperiod)
min_count_0 = -999.0
bin_min_count_0 = -999
do ibin = ibin_zero+1, ibin_zero+num_bin_for_peak1
if ( low_occ_rate_count(ibin,iperiod) .gt.
+ count_peak1 ) then
bin_peak1 = ibin
endif
if ( low_occ_rate_count(ibin, iperiod) .ge.
+ low_occ_rate_count(ibin-1,iperiod) )
+ tot_num = tot_num + 1
if ( (low_occ_rate_count(ibin,iperiod).lt.
+ low_occ_rate_count(ibin-1,iperiod) ) .and.
+ (low_occ_rate_count(ibin,iperiod).lt.
+ low_occ_rate_count(ibin+1,iperiod) ) ) then
min_count_0 = low_occ_rate_count(ibin,iperiod)
bin_min_count_0 = ibin
endif
end do ! end of ibin loop
if ( (tot_num.ge.num_bin_one_peak) .and.
+ ((bin_peak1-ibin_zero).ge.num_bin_one_peak) .and.
+ (min_count_0 .lt. 0.0) .and.
+ (bin_min_count_0 .lt. ibin_zero ) ) then
min_rate_valley = -99.0
else if ( (min_count_0 .ge. 0.0) .and.
+ (bin_min_count_0 .gt. ibin_zero) ) then
min_rate_valley = frac1 +
+ ( float(bin_min_count_0) - 0.5)*bin_size_low_occ
endif
endif ! end of if (tot_num_diff >= num_allow_diff)
* *... calculate the final criterion cut for each period ... *
if (iperiod .le. num_period) then
low_occ_rate_cut(iperiod) = min_rate_valley/2.0
else
low_occ_rate_cut(iperiod) = min_rate_valley
endif
end do ! end of iperiod loop
endif ! end of if (num_period <= 0 )
* ---------------------------------------------------------------------
90004 format('Neutrino run is too short for low occupancy check')
return
end
*endfile member=anxx_search_lowocc_cut