*To estimate the cosmic ray background
*
* EHT
*
h/file 1  /disk/upmu/2007nov_officialdataset/ntuple/2007nov_allmuons_nearhorizontal_costh0-0.1.hbk

opt liny
set ylab .8
set xlab .8
set vsiz .3
set tsiz .3
set asiz .3
set xval .4
application comis
	real function func(x)
	common /pawpar/par(3)
	func = par(1) + par(2)*exp(-x*par(3))
	return
	end
exit
vec/cre par(3) R 10 10 10
opt fit
set fit 1111

nt/cut $1 fit_pid=1&&fit_mom>1.6
nt/cut $2 fit_pid=2&&fit_len>700
nt/cut $3 fit_dir(3)>-.08&&fit_dir(3)<0
nt/cut $4 fit_dir(3)>-.08
nt/cut $5 sh_delta>.5&&sh_chi1p>50
nt/cut $6 mod(atan2(fit_dir(2),fit_dir(1))*(180/3.14)+360,360.0)>60
nt/cut $7 mod(atan2(fit_dir(2),fit_dir(1))*(180/3.14)+360,360.0)<240
nt/cut $8 mod(atan2(fit_dir(2),fit_dir(1))*(180/3.14)+360,360.0)<310

h/del 50,51
1d 50 'zenith stop' 10 -.08 0
1d 51 'zenith thru' 10 -.08 0
1d 60 'zenith stop' 10 -.08 0
1d 61 'zenith thru' 10 -.08 0
nt/proj 50 1.fit_dir(3) $3&&$1&&$6&&$7
nt/proj 51 1.fit_dir(3) $3&&$2&&$6&&$8
nt/proj 60 1.fit_dir(3) $3&&$1&&!($6&&$7)
nt/proj 61 1.fit_dir(3) $3&&$2&&!($6&&$8)

zone 2 2
set mscf 1
set mtyp 1
h/fit 50 func 0 3 par
h/pl 50 PE
set pmci 3
h/pl 60 PSE
h/fit 51 func 0 3 par
set pmci 1
h/pl 50 PE
set pmci 3
h/pl 61 PSE
set pmci 1

set mtyp 1
opt nstat
graph/hplot/null 0 360 -1 .08
nt/pl 1.(-fit_dir(3))%mod(atan2(fit_dir(2),fit_dir(1))*(180/3.14)+360,360.0) _
	$4&&$1 ! ! ! S
graph/hplot/atit 'stop [g]' '-z'
set plci 2
set dmod 2
line 0 0 360 0
set plci 4
line 60 -1 60 0.08
line 310 -1 310 0.08

graph/hplot/null 0 360 -1 .08
nt/pl 1.(-fit_dir(3))%mod(atan2(fit_dir(2),fit_dir(1))*(180/3.14)+360,360.0) _
	$4&&$2 ! ! ! S
graph/hplot/atit 'thru [g]' '-z'
set plci 2
line 0 0 360 0
set plci 4
line 60 -1 60 0.08
line 240 -1 240 0.08

h/del 70,71,80,81
1d 70 'stopping phi' 36 0 360
1d 80 'stopping phi' 36 0 360
1d 71 'thru phi' 36 0 360
1d 81 'thru phi' 36 0 360
nt/proj 70 1.mod(atan2(fit_dir(2),fit_dir(1))*(180/3.14)+360,360.0) _
	$4&&$1
nt/proj 71 1.mod(atan2(fit_dir(2),fit_dir(1))*(180/3.14)+360,360.0) _
	$4&&$2
h/file 2 /disk/upmu/2007nov_officialdataset/ntuple/2007nov_upmu_final_eventsample.hbk
cd //lun2
nt/proj 80 1.mod(atan2(fit_dir(2),fit_dir(1))*(180/3.14)+360,360.0) _
        $1&&fit_dir(3)>0
nt/proj 81 1.mod(atan2(fit_dir(2),fit_dir(1))*(180/3.14)+360,360.0) _
        $2&&fit_dir(3)>0

opt logy
set dmod 1
h/pl 70
set dmod 2
h/pl 80 S
set dmod 1
h/pl 71
set dmod 2
h/pl 80 S

cdir  //lun1
set plci 1
set dmod 1