JOMGsim.sh

Go to the documentation of this file.

#!/usr/bin/env zsh
#
#
# \author mdejong
#
version=1.0
script=${0##*/}
 
# ------------------------------------------------------------------------------------------
#
#                         Utility script to process OMGsim output.
#
# ------------------------------------------------------------------------------------------
 
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       OMGSIM_WORKDIR         2
set_variable:  WORKDIR     OMGSIM_WORKDIR         ./
set_variable:  FACTOR      OMGSIM_FACTOR          0.5
set_variable:  RECYCLING   OMGSIM_RECYCLING       1
 
if do_usage $*; then
    usage "$script (<detector file>|URL) <file name pattern> <OMGsim data card>"
fi
 
if (( $# == 3 )); then
    if [[ $argv[1] == "http"* ]]; then
        set_variable  URL                 $argv[1]
        set_variable  DETECTOR            singleDOM_OMGsim.detx
    else
        set_variable  URL
        set_variable  DETECTOR            $argv[1]
    fi
    set_variable  FILE_NAME_PATTERN   $argv[2]
    set_variable  OMGSIM              $argv[3]
else
    fatal "Wrong number of arguments."
fi
 
function getv() { JPrintProperties -f $OMGSIM -k $1 -@ "cc=#" }
 
typeset -A BUFFER
 
let "BUFFER[K40Glass]   = `getv glass.K40  ` * `getv mass.glass`"
let "BUFFER[Th232Glass] = `getv glass.Th232` * `getv mass.glass`"
let "BUFFER[U238Glass]  = `getv glass.U238 ` * `getv mass.glass`"
let "BUFFER[K40Water]   = `getv water.K40  ` * (4.0/3.0) * acos(-1.0) * (`getv radius.K40 `**3 - `getv radius.glass`**3)"
let "BUFFER[U238Water]  = `getv water.U238 ` * (4.0/3.0) * acos(-1.0) * (`getv radius.U238`**3 - `getv radius.glass`**3)"
 
echo Radio activities:
 
for KEY VALUE in ${(kv)BUFFER}; do
    printf "%-12s %12.3e [Bq]\n" $KEY $VALUE 
done
 
if [[ "$URL" != "" ]]; then
 
    echo Downloading data from $URL:
 
    for KEY in ${(k)BUFFER}; do
        if [[ ! -f ${FILE_NAME_PATTERN/\%/${KEY}} ]]; then
            wget $URL/${FILE_NAME_PATTERN/\%/${KEY}}
        fi
    done
 
    echo Downloading detector file from $URL:
    if [[ ! -f $DETECTOR ]]; then
        wget $URL/$DETECTOR
    fi
fi
 
#ideal PMT simulation
cat>$WORKDIR/PMT_parameters.A.txt<<EOF
%.TTS_ns           =  1.5
%.gain             =  1.0
%.gainSpread       =  0.0
%.PunderAmplified  =  0.0
%.slewing          =  0
EOF
 
#realistic PMT simulation (R12199)
cat>$WORKDIR/PMT_parameters.B.txt<<EOF
%.TTS_ns           = -1.0
%.gain             =  1.0
%.gainSpread       =  0.4
%.PunderAmplified  =  0.05
%.slewing          =  1
EOF
 
OPTIONS=("BELL_SHAPE = 2.0;")
 
for TYPE in A B; do
 
    rm -f $TYPE.log
 
    for KEY VALUE in ${(kv)BUFFER}; do
 
        set_variable  INPUT_FILE   ${FILE_NAME_PATTERN/\%/$KEY}
        set_variable  OUTPUT_FILE  \[${KEY}\].$TYPE.root
        
        echo "Processing type $A $INPUT_FILE -> $KEY.$TYPE.log"
 
        nohup \
                $JPP_DIR/examples/JCalibrate/JOMGsim  \
                -a $DETECTOR                          \
                -f $INPUT_FILE                        \
                -o $WORKDIR/$OUTPUT_FILE              \
                -R $VALUE                             \
                -P $WORKDIR/PMT_parameters.$TYPE.txt  \
                -P "QE = $FACTOR"                     \
            -N $RECYCLING                         \
                -d $DEBUG --!                         >&  $KEY.$TYPE.log &
    done
 
    wait
 
    KEYS=(${(k)BUFFER})
 
    JOpera1D \
        -f$WORKDIR/\[${^KEYS}\].$TYPE.root:'^M$'  \
        -u Add                                    \
        -O =                                      \
        -o $WORKDIR/M.$TYPE.root                  >& /dev/null
 
    JOpera2D \
        -f$WORKDIR/\[${^KEYS}\].$TYPE.root:1.2R   \
        -u Add                                    \
        -O =                                      \
        -o $WORKDIR/\[Total\].$TYPE.root          >& /dev/null
 
    H2=1.2R
 
    JPrintRange2D -f $WORKDIR/\[Total\].$TYPE.root:$H2                             |  read X1 Y1 Z1 X2 Y2 Z2
    JPrintResult  -f $WORKDIR/\[Total\].$TYPE.root:$H2  -F GetNbinsX -F GetNbinsY  |  read NX NY
    JPrintResult  -f $WORKDIR/\[Total\].$TYPE.root:$H2  -F GetSumOfWeights         |  read W
 
    echo "Processing $KEYS[*] -> $TYPE.log"
 
    printf "Total rate %7.3f [Hz]\n"  $(($W * ($Y2 - $Y1) / ($NY - 1)))            >>  $TYPE.log
 
    echo "KM3NeT.hh"                                                               >>  $TYPE.log
    $JPP_DIR/examples/JGizmo/JPrint1D \
        -f $WORKDIR/M.$TYPE.root:M  |  awk 'NR==3,NR==6 { printf ", %6.2f", $4 }'  >>  $TYPE.log
    echo                                                                           >>  $TYPE.log
 
    JFitK40 \
        -a $DETECTOR                              \
        -f $WORKDIR/\[Total\].$TYPE.root          \
        -@${^OPTIONS}                             \
        -o $WORKDIR/fit.$TYPE.root                \
        -w                                        \
        -M                                                                         >>  $TYPE.log
done