#! /bin/tcsh -f
#
#	Automatically set up the aperture widths on 8/12.3.1
#
#
#
# characterize and re-define slits
set beamline = `beamline.com`

#goto horiz

echo "opening slits all the way..."
# open everything WIDE
blmotor.com Aperture Upper -1000 > /dev/null
blmotor.com Aperture Lower -1000  > /dev/null
blmotor.com Aperture Inboard -1000  > /dev/null
blmotor.com Aperture Outboard -1000  > /dev/null

sleep 1

# record these for posterity
set pos = `blmotor.com Aperture Upper | awk '/is now at/{print $NF}'` 
change.com found Aperture Upper home at $pos
set pos = `blmotor.com Aperture Lower | awk '/is now at/{print $NF}'` 
change.com found Aperture Lower home at $pos
set pos = `blmotor.com Aperture Inboard | awk '/is now at/{print $NF}'` 
change.com found Aperture Inboard home at $pos
set pos = `blmotor.com Aperture Outboard | awk '/is now at/{print $NF}'` 
change.com found Aperture Outboard home at $pos

# home them
if(0 && "$beamline" != "1231") then
    echo "setpos Aperture Upper 0\r" | sock_exchange.tcl beamline 10001 1
    echo "setpos Aperture Lower 0\r" | sock_exchange.tcl beamline 10001 1
    echo "setpos Aperture Line 9 0\r" | sock_exchange.tcl beamline 10001 1
    echo "setpos Aperture Line 7 0\r" | sock_exchange.tcl beamline 10001 1
endif

###set Aperture = Upper


####################################################
#first set up apertures - Upper and Lower (Vertical)
####################################################

foreach Aperture ( Upper Lower )

set aperture = `echo $Aperture | awk -F "[_]" '{print $1,$2}'`
echo "scanning Aperture $aperture for outer edge of beam"
set corner = ""

foreach step ( 1 0.1 )

wait_for_beam.com

set pos = `blmotor.com Aperture $aperture | awk '/is now at/{print $NF}'`
if($corner != "") set pos = `echo $corner | awk '{print $1-2}'` 
sleep 2
set Izero = `blmotor.com Izero | awk '{print $NF}'`
set maxIzero = `blmotor.com Izero | awk '{print $NF}'`
if(-e Ap_${Aperture}.log) mv Ap_${Aperture}.log Ap_${Aperture}.log.old
echo -n "" >! Ap_${Aperture}.log
set sloping = 0
while ( `echo $pos | awk '{print ( $1 < 100 )}'` )
set pos = `echo $pos $step | awk '{print $1+$2}'`

blmotor.com Aperture $aperture $pos  > /dev/null

sleep 1
set lastIzero = "$Izero"
set Izero = `blmotor.com Izero | awk '{print $NF}'`

set prct  = `echo $Izero $maxIzero | awk '{printf "%d", 100*$1/$2}'`
set mean = `awk '{sum+=$2;++n} END{if(n) print sum/n}' Ap_${Aperture}.log`
set sigma = `head -10 Ap_${Aperture}.log | awk -v mean=$maxIzero '{sum+=($2-mean)^2;++n} END{if(n) print sqrt(sum/n)}'`
set sigdev = `awk -v mean="$maxIzero" -v sigma="$sigma" '{if(sigma+0) print int(($2-mean)/sigma)}' Ap_${Aperture}.log | tail -1`
if("$sigdev" == "") set sigdev = 0
if($sigdev < -2) then
    if ("$corner" == "") set corner = $pos
    @ sloping = ( $sloping + 1 )
endif
if( $prct < 20 ) break
echo "$pos $Izero $prct    $sigdev $sloping" | tee -a Ap_${Aperture}.log
if($sloping > 10 && $prct < 90) break
if($sloping > 10 && $prct < 95 && $step == 0.1) break
if( $prct < 70 ) break
end

if("$corner" == "") set corner = `awk 'NR==1{max=$2} $2<0.9*max{print $1;exit}' Ap_${Aperture}.log`
set corner = `echo "$corner $step" | awk '{print $1-3*$2}'`
if($step == 1) set slope = `echo $Izero $lastIzero $step | awk '{print ($1-$2)/$3}'`

set aperture = `echo $Aperture | awk -F "[_]" '{print $1,$2}'`
rm -f fit.log >& /dev/null
gnuplot -persist << EOF
knife(x) = x*((1+sgn(x))/2)
slope = $slope
corner = $corner
base = $maxIzero
scan(x) = base+slope*knife(x-corner)

set title "Izero vs Aperture $aperture"
set xlabel "Aperture $aperture"
set ylabel "Izero (nA)"
plot 'Ap_${Aperture}.log', scan(x)
fit scan(x) 'Ap_${Aperture}.log' via base
plot 'Ap_${Aperture}.log', scan(x)
fit scan(x) 'Ap_${Aperture}.log' via slope
plot 'Ap_${Aperture}.log', scan(x)
fit scan(x) 'Ap_${Aperture}.log' via base,corner
plot 'Ap_${Aperture}.log', scan(x)
fit scan(x) 'Ap_${Aperture}.log' via base,slope,corner
fit scan(x) 'Ap_${Aperture}.log' via base,slope,corner
plot 'Ap_${Aperture}.log', scan(x)
EOF

# move the slit over to the "corner"
echo "corner = $corner"
set corner = `awk '$1=="corner" && $4=="+/-"{print $3}' fit.log | tail -1`
set slope = `awk '$1=="slope" && $4=="+/-"{print $3}' fit.log | tail -1`
echo "corner = $corner"
blmotor.com Aperture $aperture $corner
end

change.com found beam edge at Aperture $aperture $corner

end

# make this the "outer limits"
set pos = `blmotor.com Vertical Aperture Position | awk '/is now at/{print $NF}'`
change.com found beam at Vertical Aperture Position $pos
if("$beamline" != "1231") then
    echo "setpos Vertical Aperture Position 0\r" | sock_exchange.tcl beamline 10001 1
endif



###################################################################
# now find the point where the aperture paddles cross
set direction = Vertical

set corner = ""

set step = 0.1
foreach step ( 0.1 )

wait_for_beam.com

set pos = `blmotor.com $direction Aperture Size | awk '{print $NF+1}'`
if($corner != "") then
    set pos = `echo $corner | awk '{print $1+2}'`
    blmotor.com $direction Aperture Size $pos
    blmotor.com $direction Aperture Size $pos
endif
sleep 2
set maxIzero = `blmotor.com Izero | awk '{print $NF}'`
set startpos = $pos
set Izero = `blmotor.com Izero | awk '{print $NF}'`
if(-e ${direction}_Aperture_Size.log) mv ${direction}_Aperture_Size.log ${direction}_Aperture_Size.log.old
echo -n "" >! ${direction}_Aperture_Size.log
set flat = 0
while ( $flat < 15 )
set pos = `echo $pos $step | awk '{print $1-$2}'`

blmotor.com $direction Aperture Size $pos > /dev/null

sleep 1
set lastIzero = "$Izero"
set Izero = `blmotor.com Izero | awk '{print $NF}'`

set mean = `tail -10 ${direction}_Aperture_Size.log | awk '{sum+=$2;++n} END{if(n) print sum/n}'`
set sigma = `tail -10 ${direction}_Aperture_Size.log | awk -v mean=$mean '{sum+=($2-mean)^2;++n} END{if(n) print sqrt(sum/n)}'`
set sigdev = `echo $maxIzero $sigma | awk '{if($2+0) print int($1/$2)}'`
if("$sigdev" == "") set sigdev = 0
set belowmax = `echo $maxIzero $mean | awk '{print ($2 != "" && $2/$1 < 0.8)}'`
if($sigdev > 30 && $belowmax) then
    if ("$corner" == "") set corner = $pos
    @ flat = ( $flat + 1 )
endif
echo "$pos $Izero    $belowmax $sigma $sigdev $flat" | tee -a ${direction}_Aperture_Size.log
end

# filter out upper edge
awk '$3>0' ${direction}_Aperture_Size.log >! crap
mv crap ${direction}_Aperture_Size.log

set slope = `awk '{x=$1;y=$2;getline;print ($2-y)/($1-x);exit}' ${direction}_Aperture_Size.log`
rm -f fit.log >& /dev/null
gnuplot -persist << EOF
knife(x) = x*((1+sgn(x))/2)
slope = $slope
corner = $corner
base = 0.0001
scan(x) = base+slope*knife(x-corner)

set title "Izero vs $direction Aperture"
set xlabel "$direction Aperture Size"
set ylabel "Izero (nA)"
fit scan(x) '${direction}_Aperture_Size.log' via base
plot '${direction}_Aperture_Size.log', scan(x)
fit scan(x) '${direction}_Aperture_Size.log' via base,corner
#fit scan(x) '${direction}_Aperture_Size.log' via base,slope,corner
fit scan(x) '${direction}_Aperture_Size.log' via base,slope,corner
plot '${direction}_Aperture_Size.log', scan(x)
EOF

set test = `awk '$1=="corner" && $4=="+/-"{print (sqrt($3*$3)>$5)}' fit.log | tail -1`
if(! $test) then
    echo "AAAAK! "
    set in = ( $< )
endif

# move slits back so they are just barely closed
set corner = `awk '$1=="corner" && $4=="+/-"{print $3}' fit.log | tail -1`
set slope = `awk '$1=="slope" && $4=="+/-"{print $3}' fit.log | tail -1`
blmotor.com $direction Aperture Size $corner
blmotor.com $direction Aperture Size $corner


# define this as the "closed" position
set pos = `blmotor.com Vertical Aperture Position | awk '/is now at/{print $NF}'`
change.com found beam at Vertical Aperture Position $pos
echo "setpos $direction Aperture Position 0\r" | sock_exchange.tcl beamline 10001 1
# define this as the "closed" position
change.com paddles crossed at $direction Aperture Size $corner
echo "setpos $direction Aperture Size 0\r" | sock_exchange.tcl beamline 10001 1


# open it back up
blmotor.com $direction Convergance Angle 1 > /dev/null


# check that divergence actually makes sense























horiz:
#########################################################
#now for the other set of apertures 11 and 1 (horizontal)
#########################################################

set Aperture = Outboard
foreach Aperture ( Outboard Inboard )

set aperture = `echo $Aperture | awk -F "[_]" '{print $1,$2}'`
echo "scanning Aperture $aperture for outer edge of beam"
set corner = ""

set step = 1
foreach step ( 1 0.1 )

wait_for_beam.com

set pos = `blmotor.com Aperture $aperture | awk '/is now at/{print $NF}'`
if($corner != "") then
    set pos = `echo $corner | awk '{print $1-4}'` 
endif
sleep 2
set Izero = `blmotor.com Izero | awk '{print $NF}'`
set maxIzero = `blmotor.com Izero | awk '{print $NF}'`
if(-e Ap_${Aperture}.log) mv Ap_${Aperture}.log Ap_${Aperture}.log.old
echo -n "" >! Ap_${Aperture}.log
set sloping = 0
while ( `echo $pos | awk '{print ( $1 < 100 )}'` )
set pos = `echo $pos $step | awk '{print $1+$2}'`

blmotor.com Aperture $aperture $pos  > /dev/null

sleep 1
set lastIzero = "$Izero"
set Izero = `blmotor.com Izero | awk '{print $NF}'`

set slope = `echo $Izero $lastIzero $maxIzero $step | awk '{printf "%d", 1000*sqrt((($2-$1)/$3/$4)^2)}'`
set myslope = `echo $Izero $lastIzero $maxIzero $step | awk '{printf "%d", sqrt((($2-$1)/$4)^2)}'`
set prct  = `echo $Izero $maxIzero | awk '{printf "%d", 100*$1/$2}'`
if($slope > 15) then
    if ("$corner" == "") set corner = $pos
    @ sloping = ( $sloping + 1 )
endif
if($sloping > 10 && $prct < 90) break
if($sloping > 10 && $prct < 95 && $step == 0.1) break
if( $prct < 70 ) break
echo "$pos $Izero $prct $slope $sloping $myslope" | tee -a Ap_${Aperture}.log
end

if("$corner" == "") set corner = `awk 'NR==1{max=$2} $2<0.9*max{print $1;exit}' Ap_${Aperture}.log`
set slope = -0.1
set myslope = -0.1
set aperture = `echo $Aperture | awk -F "[_]" '{print $1,$2}'`
rm -f fit.log >& /dev/null
gnuplot -persist << EOF
knife(x) = x*((1+sgn(x))/2)
slope=$myslope
corner = $corner
base = $maxIzero
scan(x) = base+slope*knife(x-corner)

set title "Izero vs Aperture $aperture"
set xlabel "Aperture $aperture"
set ylabel "Izero (nA)"
fit scan(x) 'Ap_${Aperture}.log' via base
plot 'Ap_${Aperture}.log', scan(x)
fit scan(x) 'Ap_${Aperture}.log' via slope
plot 'Ap_${Aperture}.log', scan(x)
fit scan(x) 'Ap_${Aperture}.log' via base,corner
plot 'Ap_${Aperture}.log', scan(x)
fit scan(x) 'Ap_${Aperture}.log' via base,slope,corner
fit scan(x) 'Ap_${Aperture}.log' via base,slope,corner
plot 'Ap_${Aperture}.log', scan(x)
EOF

# move the slit over to the "corner"
echo "corner = $corner"
set corner = `awk '$1=="corner" && $4=="+/-"{print $3}' fit.log | tail -1`
echo "corner = $corner"
blmotor.com Aperture $aperture $corner
end

change.com found beam edge at Aperture $aperture $corner

end

# make these the "outer limits"
# define this as the "closed" position
set pos = `blmotor.com Horizontal Aperture Position | awk '/is now at/{print $NF}'`
change.com found beam at Horizontal Aperture Position $pos
if("$beamline" != "1231") then
    echo "setpos Horizontal Aperture Position 0\r" | sock_exchange.tcl beamline 10001 1
endif

# now find the point where the aperture paddles cross
set direction = Horizontal

set corner = ""

set step = 1
foreach step ( 1 0.1 )

wait_for_beam.com

set pos = `blmotor.com $direction Aperture Size | awk '{print $NF}'`
if($corner != "") set pos = `echo $corner | awk '{print $1+2}'` 
sleep 2
if($corner == "") set maxIzero = `blmotor.com Izero | awk '{print $NF}'`
set startpos = $pos
set Izero = `blmotor.com Izero | awk '{print $NF}'`
if(-e ${direction}_Aperture_Size.log) mv ${direction}_Aperture_Size.log ${direction}_Aperture_Size.log.old
echo -n "" >! ${direction}_Aperture_Size.log
set flat = 0
while ( $flat < 15 )
set pos = `echo $pos $step | awk '{print $1-$2}'`

blmotor.com $direction Aperture Size $pos > /dev/null

sleep 1
set lastIzero = "$Izero"
set Izero = `blmotor.com Izero | awk '{print $NF}'`

set mean = `tail -10 ${direction}_Aperture_Size.log | awk '{sum+=$2;++n} END{if(n) print int(sum/n)}'`
set sigma = `tail -10 ${direction}_Aperture_Size.log | awk -v mean=$mean '{sum+=($2-mean)^2;++n} END{if(n) print sqrt(sum/n)}'`
set sigdev = `echo $maxIzero $sigma | awk '{if($2+0) print int($1/$2)}'`
if("$sigdev" == "") set sigdev = 0
set belowmax = `echo $maxIzero $mean | awk '{print ($2 != "" && $2/$1 < 0.8)}'`
if($sigdev > 30 && $mean < 1) then
    if ("$corner" == "") set corner = $pos
    @ flat = ( $flat + 1 )
endif
echo "$pos $Izero    $belowmax $sigma $sigdev $flat" | tee -a ${direction}_Aperture_Size.log
end

# filter out upper edge
awk '$3>0' ${direction}_Aperture_Size.log >! crap
mv crap ${direction}_Aperture_Size.log

set slope = 0.5
rm -f fit.log >& /dev/null
gnuplot -persist << EOF
knife(x) = x*((1+sgn(x))/2)
slope = $slope
corner = $corner
base = 0.0001
scan(x) = base+slope*knife(x-corner)

set title "Izero vs $direction Aperture"
set xlabel "$direction Aperture Size"
set ylabel "Izero (nA)"
fit slope*x+intercept '${direction}_Aperture_Size.log' via slope,intercept
slope=slope*3
fit scan(x) '${direction}_Aperture_Size.log' via corner
plot '${direction}_Aperture_Size.log', scan(x)
fit scan(x) '${direction}_Aperture_Size.log' via base,corner
#fit scan(x) '${direction}_Aperture_Size.log' via base,slope,corner
fit scan(x) '${direction}_Aperture_Size.log' via base,slope,corner
plot '${direction}_Aperture_Size.log', scan(x)
EOF

# move slits back so they are just barely closed
set corner = `awk '$1=="corner" && $4=="+/-"{print $3}' fit.log | tail -1`
blmotor.com $direction Aperture Size $corner
if($status) then
    sleep 1
    blmotor.com $direction Aperture Size $corner
endif
if($status) then
    sleep 1
    blmotor.com $direction Aperture Size $corner
endif

end


# make these the "zero" position
set pos = `blmotor.com $direction Aperture Position | awk '/is now at/{print $NF}'`
change.com found beam at $direction Aperture Position $pos
echo "setpos $direction Aperture Position 0\r" | sock_exchange.tcl beamline 10001 1
# define this as the "closed" position
change.com paddles crossed at $direction Aperture Size $corner
echo "setpos $direction Aperture Size 0\r" | sock_exchange.tcl beamline 10001 1


# open it back up
blmotor.com $direction Convergance Angle 4 > /dev/null

tune_M2.com
tune_chi2.com

# record these for posterity?


# now check the divergence vs flux curves
diode.com in
shutter.com open
divergence.com 3.4 0.7 > /dev/null
autoscale.com
echo -n "" >! Vdiv.log
divergence.com 3.4 0.0001 > /dev/null
sleep 1
foreach div ( `awk 'BEGIN{for(x=-0.1;x<=0.7;x+=0.02) print x}'` )
blmotor.com Vertical Convergance Angle $div  > /dev/null
sleep 1
set Iend = `blmotor.com Iend | awk '{print $NF}'`
set Izero = `blmotor.com Izero | awk '{print $NF}'`
echo "$div $Iend $Izero" | tee -a Vdiv.log
end

echo -n "" >! Hdiv.log
divergence.com 0.001 0.7  > /dev/null
sleep 1
foreach div ( `awk 'BEGIN{for(x=-0.5;x<=4;x+=0.1) print x}'` )
blmotor.com Horizontal Convergance Angle $div  > /dev/null
sleep 1
set Iend = `blmotor.com Iend | awk '{print $NF}'`
set Izero = `blmotor.com Izero | awk '{print $NF}'`
echo "$div $Iend $Izero" | tee -a Hdiv.log
end
divergence.com 2 0.35
diode.com out
shutter.com close


cat << EOF >! variables.txt
a=1
b=1
z=1
EOF

gnuplot -persist << EOF
set grid
set xrange [*:*]
set autoscale y
set autoscale x
set xtics auto
set title "flux vs divergence $beamline"
set xlabel "Horizontal Divergence"
set ylabel "Izero (nA)"
c(x)=column(x)
a=1
b=1
fit a*x 'Hdiv.log' using 1:(c(1)<2.5 && c(1)>0 ? c(3) : 1/0) via a
fit b 'Hdiv.log' using 1:(c(1)>3.5 ? c(3) : 1/0) via b
fit z 'Hdiv.log' using 1:(c(1)<0 ? c(3) : 1/0) via z
plot 'Hdiv.log' using 1:3 ti "horizontal", a*x ti "", b ti "", z ti ""
print "horiz crossover: ", b/a
print "horiz zero crossover: ", z/a
update "variables.txt"
EOF
set horz_corner_Izero = `awk '/^a/{a=$3} /^b/{b=$3} /^z/{z=$3} END{print b/a}' variables.txt`
echo "$horz_corner_Izero" >! /data/calibrations/max_Hdiv_Izero.txt
change.com horizontal divergence max flux corner at $horz_corner_Izero for Izero

gnuplot -persist << EOF
set grid
set xrange [*:*]
set autoscale y
set autoscale x
set xtics auto
set title "flux vs divergence $beamline"
set xlabel "Vertical Divergence"
set ylabel "Izero (nA)"
c(x)=column(x)
a=1
b=1
fit a*x 'Vdiv.log' using 1:(c(1)<0.22 && c(1)>0 ? c(3) : 1/0) via a
fit b 'Vdiv.log' using 1:(c(1)>0.5 ? c(3) : 1/0) via b
fit z 'Vdiv.log' using 1:(c(1)<0 ? c(3) : 1/0) via z
plot 'Vdiv.log' using 1:3 ti "vertical", a*x ti "", b ti "", z ti ""
print "vert crossover: ", b/a
print "vert zero crossover: ", z/a
update "variables.txt"
EOF
set vert_corner_Izero = `awk '/^a/{a=$3} /^b/{b=$3} /^z/{z=$3} END{print b/a}' variables.txt`
echo "$vert_corner_Izero" >! /data/calibrations/max_Vdiv_Izero.txt
change.com vertical divergence max flux corner at $vert_corner_Izero for Izero


gnuplot -persist << EOF
set grid
set xrange [*:*]
set autoscale y
set autoscale x
set xtics auto
set title "flux vs divergence $beamline"
set xlabel "Horizontal Divergence"
set ylabel "Iend (nA)"
c(x)=column(x)
a=1
b=1
fit a*x 'Hdiv.log' using 1:(c(1)<2.5 && c(1)>0 ? c(2) : 1/0) via a
fit b 'Hdiv.log' using 1:(c(1)>3.5 ? c(2) : 1/0) via b
fit z 'Hdiv.log' using 1:(c(1)<0 ? c(2) : 1/0) via z
plot 'Hdiv.log' using 1:2 ti "horizontal", a*x ti "", b ti "", z ti ""
print "horiz crossover: ", b/a
print "horiz zero crossover: ", z/a
update "variables.txt"
EOF
set horz_corner_Iend = `awk '/^a/{a=$3} /^b/{b=$3} /^z/{z=$3} END{print b/a}' variables.txt`
echo "$horz_corner_Iend" >! /data/calibrations/max_Hdiv_Iend.txt
change.com horizontal divergence max flux corner at $horz_corner_Iend for Iend


gnuplot -persist << EOF
set grid
set xrange [*:*]
set autoscale y
set autoscale x
set xtics auto
set title "flux vs divergence $beamline"
set xlabel "Vertical Divergence"
set ylabel "Iend (nA)"
c(x)=column(x)
a=1
b=1
set xlabel "Vertical Divergence"
fit a*x 'Vdiv.log' using 1:(c(1)<0.22 && c(1)>0 ? c(2) : 1/0) via a
fit b 'Vdiv.log' using 1:(c(1)>0.5 ? c(2) : 1/0) via b
fit z 'Vdiv.log' using 1:(c(1)<0 ? c(2) : 1/0) via z
plot 'Vdiv.log' using 1:2 ti "vertical", a*x ti "", b ti "", z ti ""
print "vert crossover: ", b/a
print "vert zero crossover: ", z/a
update "variables.txt"
EOF
set vert_corner_Iend = `awk '/^a/{a=$3} /^b/{b=$3} /^z/{z=$3} END{print b/a}' variables.txt`
echo "$vert_corner_Iend" >! /data/calibrations/max_Vdiv_Iend.txt
change.com vertical divergence max flux corner at $vert_corner_Iend for Iend



exit

gnuplot 
set grid
set xrange [*:*]
set autoscale y
set autoscale x
set xtics auto
set title "flux vs divergence $beamline"
set xlabel "Horizontal Divergence"
set ylabel "Izero (nA)"
a=1
b=1
fit a*x 'Hdiv.log' using 1:($1<2.5 && $1>0 ? $3 : 1/0) via a
fit b 'Hdiv.log' using 1:($1>3.5 ? $3 : 1/0) via b
fit z 'Hdiv.log' using 1:($1<0 ? $3 : 1/0) via z
plot 'Hdiv.log' using 1:3 ti "horizontal", a*x ti "", b ti "", z ti ""
print "horiz crossover: ", b/a
print "horiz zero crossover: ", z/a


set xlabel "Vertical Divergence"
fit a*x 'Vdiv.log' using 1:($1<0.22 && $1>0 ? $3 : 1/0) via a
fit b 'Vdiv.log' using 1:($1>0.5 ? $3 : 1/0) via b
fit z 'Vdiv.log' using 1:($1<0 ? $3 : 1/0) via z
plot 'Vdiv.log' using 1:3 ti "vertical", a*x ti "", b ti "", z ti ""
print "vert crossover: ", b/a
print "vert zero crossover: ", z/a


set ylabel "Iend (nA)"
fit a*x 'Hdiv.log' using 1:($1<2.5 && $1>0 ? $2 : 1/0) via a
fit b 'Hdiv.log' using 1:($1>3.5 ? $2 : 1/0) via b
fit z 'Hdiv.log' using 1:($1<0 ? $2 : 1/0) via z
plot 'Hdiv.log' using 1:2 ti "horizontal", a*x ti "", b ti "", z ti ""
print "horiz crossover: ", b/a
print "horiz zero crossover: ", z/a


set xlabel "Vertical Divergence"
fit a*x 'Vdiv.log' using 1:($1<0.22 && $1>0 ? $2 : 1/0) via a
fit b 'Vdiv.log' using 1:($1>0.5 ? $2 : 1/0) via b
fit z 'Vdiv.log' using 1:($1<0 ? $2 : 1/0) via z
plot 'Vdiv.log' using 1:2 ti "vertical", a*x ti "", b ti "", z ti ""
print "vert crossover: ", b/a
print "vert zero crossover: ", z/a



knife(x) = x*((1+sgn(x))/2)
slope = -1
corner = 1
base
scan(x) = base+slope*knife(x-corner)