Attached are instructions and scripts for carrying out a parameter space scan with ShelfMC, the simulation package for the ARIANNA detector.
Because some of the plotted outputs looked like colored stripes and did not offer any insight into how effective volume changed with some variables, I made some changes to the simulation and plotting scripts so that different maximum, minimum, and increment values can be chosen for each variable. Now rather than having fixed, hard-coded values for all variables, the parameter space scan and plotting is more flexible for use with variable inputs.
| Quote: |
|
I am trying to write a script that will plot a 2d histogram of effective volume versus two of ShelfMC's parameters.
The script prompts the user for which two parameters (out of five that we vary in our parameter space scan) to plot along the x- and y-axes, as well as what values to hold the other 3 parameters constant at. It then collects the necessary root files from simulation results, generates a plotting script, and runs the plotting script to produce a plot in pdf form.
After many struggles I have the script written to the point where it functions, but the plots don't look right. Some plots look like they could be actual data (like Veff_A_I), and others just look flat-out wrong (like Veff_R_S).
I have yet to pin down the cause of this, but hopefully will be able to sometime in the near future.
|
|
#!/bin/bash
#04/18/2017
# script will collect root files and plot as dictated by user
#Hannah Hasan
shelfmc=$HOME/ShelfMC/git_shelfmc ##set to your ShelfMC directory - should have the parameter space scan directory and a directory called "outputs". outputs may be empty
AttenMin=500
AttenMax=1000
AttenInc=100
RadiusMin=3
RadiusMax=31
RadiusInc=7
IceMin=500
IceMax=2900
IceInc=400
FirnMin=60
FirnMax=140
FirnInc=20
StDepthMin=0
StDepthMax=200
StDepthInc=50
line1='{'
line2=' TFile *f = new TFile("'"Veff_${xvar}_${yvar}.root"'");' #set later
line3=' TTree *SimParam = (TTree*)f->Get("sim_parameters");'
line4=' TCanvas *c = new TCanvas("c", "V_{eff}", 2000, 1000);'
line5=' TPad* p = new TPad("p","p",0.0,0.0,1.0,1.0);'
line6=' p->Draw();'
line7=' p->SetLeftMargin(0.15);'
line8=' p->SetRightMargin(0.15);'
line9=' p->cd();'
line10=' TH2F *h = new TH2F("''h''", "''V_{eff} vs '"$var1"' vs '"$var2"'"'",$xbins,$xmin,$xmax"', '"$ybins,$ymin,$ymax);" set later
line11=' double x;'
line12=' double y;'
line13=' double z;'
line14=' double veff=0;'
line15=' double energy;'
line16=" double avgs[$SIZE];" #set later; SIZE varies depending on which variables we are plotting
line17=' SimParam->SetBranchAddress("'"$xbranch"'", &x );' #set later
line18=' SimParam->SetBranchAddress("'"$ybranch"'", &y );' #set later
line19=' SimParam->SetBranchAddress("Veff", &z);'
line20=' int v=0;'
line21=" for(int k=1; k<=$SIZE; k++){" #set later for SIZE
line22=' veff = 0;'
line23=' for(v; v<(k*20.0); v++){' #change 20 to number of processors ShelfMC runs were split across
line24=' SimParam->GetEntry(v);'
line25=' veff += z;'
line26=' }'
line27=' avgs[(k-1)] = (veff/20.0);' #change 20 to number of processors
line28=' }'
line29=" for(int i=0; i<$SIZE; i++){" #set later for SIZE
line30=' SimParam->GetEntry((i*20));' #change 20 to number of processors
line31=' h->Fill(x,y,avgs[i]);'
line32=' }'
line33=' gStyle->SetPalette(1,0);'
line34=' h->GetXaxis()->SetTitle("'"$var1"' (m)");' #set later
line35=' h->GetYaxis()->SetTitle("'"$var2"' (m)");' #set later
line36=' h->GetZaxis()->SetTitle("V_{eff} (km^{3} str)");'
line37=' h->GetXaxis()->SetTitleOffset(1.2);'
line38=' h->GetYaxis()->SetTitleOffset(1.0);'
line39=' h->GetZaxis()->SetTitleOffset(1.3);'
line40=' gStyle->SetOptStat(0);'
line41=' h->Draw("COLZ");'
line42=' double const1;'
line43=' double const2;'
line44=' double const3;'
line45=' SimParam->SetBranchAddress("'"$const1br"'", &const1);' #set later
line46=' SimParam->SetBranchAddress("'"$const2br"'", &const2);' #set later
line47=' SimParam->SetBranchAddress("'"$const3br"'", &const3);' #set later
line48=' SimParam->SetBranchAddress("exponent", &energy);'
line49=' SimParam->GetEntry(0);'
line50=' char const_1[30];'
line51=' char const_2[30];'
line52=' char const_3[30];'
line53=' char enrgy[30];'
line54=' sprintf(enrgy, "Energy (exponent): %d", energy);'
line55=' sprintf(const_1, "'"$const1br"' %d m", const1);'
line56=' sprintf(const_2, "'"$const2br"' %d m", const2);'
line57=' sprintf(const_3, "'"$const3br"' %d m", const3);' #set later
line58=' TPaveText *t = new TPaveText(0,0,0.07,0.3,"brNDC");'
line59=' t->SetLabel("Simulation Parameters:");'
line60=' t->AddText(enrgy);'
line61=' t->AddText(const_1);'
line62=' t->AddText(const_2);'
line63=' t->AddText(const_3);'
line64=' t->Draw();'
line65=' c->Modified();'
line66=' c->Update();'
line67=' c->SaveAs("'"Veff_${xvar}_${yvar}.pdf"'");' #set later
line68='}'
echo
echo
echo 'This script will allow you to choose two variables to plot effective volume against.'
echo
echo 'Enter the name of the run you wish to plot data for.'
echo
read -p 'Name of run: ' runName
echo
echo 'Thank you!'
echo
echo 'Please make your selections from the following list:'
echo '(A) - Attenuation Length'
echo '(R) - Antenna Radius'
echo '(I) - Ice Thickness'
echo '(F) - Firn Thickness'
echo '(S) - Station Depth'
echo
echo '** please use capital letters in your selection **'
echo
read -p 'Select the variable you wish to plot along the X-axis: ' xvar
read -p 'Now select the variable you wish to plot along the Y-axis: ' yvar
echo
cd $shelfmc/outputs
mkdir ${runName}_${xvar}_${yvar}
cd $shelfmc/$runName
attenprompt='Select what value at which to hold constant Attenuation Length'
atteninc=' > Choose only values between 500-1000, at increments of 100: '
radiusprompt='Select what value at which to hold constant Antenna Radius'
radiusinc=' > Choose only values between 3-31, at increments of 7: '
iceprompt='Select what value at which to hold constant Ice Thickness'
iceinc=' > Choose only values between 500-2900, at increments of 400: '
firnprompt='Select what value at which to hold constant Firn Thickness'
firninc=' > Choose only values between 60-140, at increments of 20: '
stprompt='Select what value at which to hold constant Station Depth'
stinc=' > Choose only values between 0-200, at increments of 50: '
if ([ "$xvar" = "A" ] && [ "$yvar" = "R" ]) || ([ "$yvar" = "A" ] && [ "$xvar" = "R" ]); then
echo $iceprompt
echo $iceinc
read const1
const1br='icethick'
echo $firnprompt
echo $firninc
read const2
const2br='firndepth'
echo $stprompt
echo $stinc
read const3
const3br='station_depth'
echo
echo 'Adding files...'
echo
for ((L=$AttenMin;L<=$AttenMax;L+=$AttenInc))
do
cd Atten_Up$L
for ((AR=$RadiusMin;AR<=$RadiusMax;AR+=$RadiusInc))
do
cd AntennaRadius$AR/IceThick$const1/FirnThick$const2/StationDepth$const3
cp *.root $shelfmc/outputs/${runName}_${xvar}_${yvar}/Result_${L}_${AR}.root
cd ../../../../
done
cd ..
done
elif ([ "$xvar" = "A" ] && [ "$yvar" = "I" ]) || ([ "$yvar" = "A" ] && [ "$xvar" = "I" ]); then
echo $radiusprompt
echo $radiusinc
read const1
const1br='st4_r'
echo $firnprompt
echo $firninc
read const2
const2br='firndepth'
echo $stprompt
echo $stinc
read const3
const3br='station_depth'
echo
for ((L=$AttenMin;L<=$AttenMax;L+=$AttenInc))
do
cd Atten_Up$L/AntennaRadius$const1
for ((Ice=$IceMin;Ice<=$IceMax;Ice+=$IceInc))
do
cd IceThick$Ice/FirnThick$const2/StationDepth$const3
cp *.root $shelfmc/outputs/${runName}_${xvar}_${yvar}/Result_${L}_${Ice}.root
cd ../../../
done
cd ../../
done
elif ([ "$xvar" = "A" ] && [ "$yvar" = "F" ]) || ([ "$yvar" = "A" ] && [ "$xvar" = "F" ]); then
echo $radiusprompt
echo $radiusinc
read const1
const1br='st4_r'
echo $iceprompt
echo $iceinc
read const2
const2br='icethick'
echo $stprompt
echo $stinc
read const3
const3br='station_depth'
echo
for ((L=$AttenMin;L<=$AttenMax;L+=$AttenInc))
do
cd Atten_Up$L/AntennaRadius$const1/IceThick$const2
for ((FT=$FirnMin;FT<=$FirnMax;FT+=$FirnInc))
do
cd FirnThick$FT/StationDepth$const3
cp *.root $shelfmc/outputs/${runName}_${xvar}_${yvar}/Result_${L}_${FT}.root
cd ../..
done
cd ../../..
done
elif ([ "$xvar" = "A" ] && [ "$yvar" = "S" ]) || ([ "$yvar" = "A" ] && [ "$xvar" = "S" ]); then
echo $radiusprompt
echo $radiusinc
read const1
const1br='st4_r'
echo $iceprompt
echo $iceinc
read const2
const2br='icethick'
echo $firnprompt
echo $firninc
read const3
const3br='firndepth'
echo
for ((L=$AttenMin;L<=$AttenMax;L+=$AttenInc))
do
cd Atten_Up$L/AntennaRadius$const1/IceThick$const2/FirnThick$const3
for ((SD=$StDepthMin;SD<=$StDepthMax;SD+=$StDepthInc))
do
cd StationDepth$SD
cp *.root $shelfmc/outputs/${runName}_${xvar}_${yvar}/Result_${L}_${SD}.root
cd ../
done
cd ../../../..
done
elif ([ "$xvar" = "R" ] && [ "$yvar" = "I" ]) || ([ "$yvar" = "R" ] && [ "$xvar" = "I" ]); then
echo $attenprompt
echo $atteninc
read const1
const1br='atten_up'
echo $firnprompt
echo $firninc
read const2
const2br='firndepth'
echo $stprompt
echo $stinc
read const3
const3br='station_depth'
echo
for ((AR=$RadiusMin;AR<=$RadiusMax;AR+=$RadiusInc))
do
cd Atten_Up$const1/AntennaRadius$AR
for ((Ice=$IceMin;Ice<=$IceMax;Ice+=$IceInc))
do
cd IceThick$Ice/FirnThick$const2/StationDepth$const3
cp *.root $shelfmc/outputs/${runName}_${xvar}_${yvar}/Result_${AR}_${Ice}.root
cd ../../..
done
cd ../..
done
elif ([ "$xvar" = "R" ] && [ "$yvar" = "F" ]) || ([ "$yvar" = "R" ] && [ "$xvar" = "F" ]); then
echo $attenprompt
echo $atteninc
read const1
const1br='atten_up'
echo $iceprompt
echo $iceinc
read const2
const2br='icethick'
echo $stprompt
echo $stinc
read const3
const3br='station_depth'
echo
for ((AR=$RadiusMin;AR<=$RadiusMax;AR+=$RadiusInc))
do
cd Atten_Up$const1/AntennaRadius$AR/IceThick$const2
for ((FT=$FirnMin;FT<=$FirnMax;FT+=$FirnInc))
do
cd FirnThick$FT/StationDepth$const3
cp *.root $shelfmc/outputs/${runName}_${xvar}_${yvar}/Result_${AR}_${FT}.root
cd ../..
done
cd ../../..
done
elif ([ "$xvar" = "R" ] && [ "$yvar" = "S" ]) || ([ "$yvar" = "R" ] && [ "$xvar" = "S" ]); then
... 215 more lines ...
|