JGedanken.sh

Go to the documentation of this file.

#!/bin/zsh
#
# \author mdejong
version=1.0
script=${0##*/}
 
# ------------------------------------------------------------------------------------------
#
#                        Auxiliary script for Gedanken experiment.
#
# ------------------------------------------------------------------------------------------
 
if [ -z $JPP_DIR ]; then
    echo "Variable JPP_DIR undefined."
    exit
fi
 
source $JPP_DIR/setenv.sh $JPP_DIR
 
zmodload zsh/mathfunc
 
set_variable:  DEBUG             SIRENE_DEBUG        2
set_variable:  WORKDIR           SIRENE_WORKDIR      ${TMPDIR:-/tmp}/
set_variable   NUMBER_OF_EVENTS  1000000
 
# detector
 
set_variable   R_M               25.0
set_variable   Z_M               50.0
set_variable   MODULE            $R_M 0.0 $Z_M               # position
 
# particle
 
set_variable   E_GEV             20                          # energy
set_variable   POS               0.0 0.0 0.0                 # position 
set_variable   DIR               0.0 0.0 1.0                 # direction
set_variable   TYPE              13                          # PDG
 
# graphics
 
set_variable:  FORMAT            GRAPHICS_FORMAT     gif
set_variable+  BATCH             GRAPHICS_BATCH      -B
 
if do_usage $*; then
    usage "$script [number of events]"
fi
 
case $# in 
    1) set_variable  NUMBER_OF_EVENTS  $1;;
esac
 
set_variable   DETECTOR          $WORKDIR/detector.detx
 
if [[ ! -f $DETECTOR ]]; then
 
    JModule \
        -D 1001                             \
        -M    1                             \
        -P "$MODULE"                        \
        -o $DETECTOR                        \
        -d $DEBUG
fi
 
if [[ ! -f $WORKDIR/gedanken.root ]]; then
    
    $JPP_DIR/examples/JSirene/JGedanken \
        -o $WORKDIR/gedanken.root           \
        -P "$POS"                           \
        -D "$DIR"                           \
        -E $E_GEV                           \
        -T $TYPE                            \
        -n $NUMBER_OF_EVENTS                \
        -a $DETECTOR \
        -d $DEBUG
fi
 
if [[ ! -f $WORKDIR/sirene.root ]]; then
 
    JSirene \
        -a $DETECTOR                        \
        -f $WORKDIR/gedanken.root           \
        -o $WORKDIR/sirene.root             \
        -F $JPP_DATA/I%p.dat                \
        -N 0                                \
        -d $DEBUG
fi
 
if [[ ! -f $WORKDIR/domino.root ]]; then
 
    $JPP_DIR/examples/JSirene/JDomino \
        -a $DETECTOR                        \
        -f $WORKDIR/sirene.root             \
        -o $WORKDIR/domino.root             \
        -d $DEBUG
fi
 
JPlot1D \
    -f $WORKDIR/"domino.root:h\[%\]"    \
    -y "1.0e-6 1.0e0" -Y                \
    -> "#Deltat [ns]"                  \
    -\^ "dP/dt npe/ns"                  \
    -T ""                               \
    -O "HIST]["                         \
    -o pdf.$FORMAT                      $BATCH
 
set_variable   PDF        $JPP_DATA/J%p.dat
set_variable   EPS        1.0E-6
set_variable   PI         $((acos(-1)))
 
# constrain angle between [epsilon, pi - epsilon]
#
function constrain()
{
    if   (( $1 > $PI - $EPS )); then
        echo $(($PI - $EPS))
    elif (( $1 < $EPS )); then
        echo $EPS
    else
        echo $1
    fi
}
 
for id in {1..31}; do
 
    echo "Generating graphics PMT $id"
 
    getPMT -a $DETECTOR -p $id  |  read ID X Y Z DX DY DZ T0 STATUS
 
    let "THETA = acos($DZ)"
 
    if   (( $DX > +$EPS )); then
        let "PHI = atan(fabs($DY/$DX))"
    elif (( $DX < -$EPS )); then
        let "PHI = $PI - atan(fabs($DY/$DX))"
    else
        let "PHI = 0.5*$PI"
    fi
 
    set_variable THETA  `constrain $THETA`
    set_variable PHI    `constrain $PHI`
 
    JPlotPDF \
        -f ${PDF/\%/1}                        \
        -f ${PDF/\%/2}                        \
        -f ${PDF/\%/3}                        \
        -f ${PDF/\%/4}                        \
        -f ${PDF/\%/5}                        \
        -f ${PDF/\%/6}                        \
        -D "$THETA $PHI"                      \
        -R $R_M -E $E_GEV                     \
        -H "860 -10 850"                      \
        -o $WORKDIR/f1\[${id}\].root
 
    set_variable  TITLE `printf "PMT[%d] (%5.3f,%5.3f)" $id $THETA $PHI`
 
    JPlot1D \
        -f $WORKDIR/"domino.root:h\[${id}\]"  \
        -f $WORKDIR/f1\[${id}\].root:h0       \
        -y "1.0e-6 1.0e0" -Y                  \
        -> "#Deltat [ns]"                    \
        -\^ "dP/dt npe/ns"                    \
        -T "$TITLE"                           \
        -O "HIST]["                           \
        -o $WORKDIR/pdf_${id}.$FORMAT         $BATCH
done
 
montage \
    -tile 8x4                                   \
    -geometry +0+0                              \
    $WORKDIR/pdf_{1..31}.$FORMAT                \
    PDF.$FORMAT                                 >& /dev/null