Well that was only a part of it. Similar if statements in my script posed no problem. I guess my question is whether there are any other cases will ignore an if statement simply because filling a histogram takes priority, or some other similar thing. If this is not the case, here is my entire code:
void plots(int Nev=1000, double binSize=0.5, int rMin=0, char *fName="./comparison_tests/test_internal_gammas/test_internal_gammas.root", double ringsize=2, const bool radYes=0, const bool ravYes=0, const bool sYes=1, const bool centYes=1) {
//load library definition
gSystem->Load("./tools/ConfigHandler/lib/libEventData.so");
gSystem->Load("./tools/ConfigHandler/lib/libConfigHandler.so");
//open file
TFile* fdata = new TFile(fName);
//~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
// CONFIGURATION
//~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
TTree* config = (TTree*)(fdata->Get("ConfigOptions"));
//find the vessel size, livetime, mass, etc.
char acx[30]="", acy[30]="", acz[30]="";
//char acxx[30]="", acyy[30]="", aczz[30]="";
char lar_top_d[30]="", lar_bot_d[30]="", lar_h[30]="";
char ppd[30]="", pgen[30]="";
config->SetBranchAddress("active_center_z",acz);
config->SetBranchAddress("active_center_x",acx);
config->SetBranchAddress("active_center_y",acy);
config->SetBranchAddress("lar_top_d",lar_top_d);
config->SetBranchAddress("lar_bot_d",lar_bot_d);
config->SetBranchAddress("lar_h",lar_h);
config->GetEntry(0);
const Double_t active_center_z = atof(acz);
const Double_t active_center_x = atof(acx);
const Double_t active_center_y = atof(acy);
const double rtop = atof(lar_top_d)/2.;
const double rbot = atof(lar_bot_d)/2.;
const double h = atof(lar_h);
const int overBin= (int)(1/binSize+0.5); //1 over binSize
const int nBins=rtop*overBin+1; //number of bins used to split the radius of the LAr chamber
const double printBins=(rtop-rMin)*overBin; //number of bins printed
const int lessRings = nBins-(1+ringsize*overBin); //determines number of bins to cut off to simulate upper ring
const double maxR=rtop; //outer edge of detector, or outer radius to be studied for average energy deposits in a radius
//~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
// Prepare Histogram
//~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
if (radYes==1) {
//radgen plots the average energy deposit made by a particle in each bin (dr)
TH1D *radgen = new TH1D("radgen", "Average Energy Deposit in dr (for all events)", printBins, rMin, rtop);
radgen->GetXaxis()->SetNdivisions(rtop/2);
radgen->GetXaxis()->CenterTitle();
radgen->GetXaxis()->SetTitle("Radius in cm");
radgen->GetYaxis()->CenterTitle();
radgen->GetYaxis()->SetTitle("Average Energy Deposit in keV");
}
if (ravYes==1) {
//ravgen plots the average accumulation of edeps looping through the radius (for all events)
TH1D *ravgen = new TH1D("ravgen", "Average Total Energy Deposit at Radius R", printBins,rMin,rtop);
ravgen->GetXaxis()->SetNdivisions(rtop/2);
ravgen->GetXaxis()->CenterTitle();
ravgen->GetXaxis()->SetTitle("Radius in cm");
ravgen->GetYaxis()->CenterTitle();
ravgen->GetYaxis()->SetTitle("Accumulation of Edeps in keV");
}
if (sYes==1) {
//s2s1 plots the average ratio of s2/s1 with increasing ring sizes
TH1D *s2s1 = new TH1D("s2s1", "S2/S1 vs. Rrid Width", printBins, rMin, rtop);
s2s1->GetXaxis()->SetNdivisions(rtop/2);
s2s1->GetXaxis()->CenterTitle();
s2s1->GetXaxis()->SetTitle("Grid Width in cm");
s2s1->GetYaxis()->CenterTitle();
s2s1->GetYaxis()->SetTitle("Average S2/S21");
}
if (centYes==1) {
//centroid plots the s2/s1 ratio as a function of R_centroid on an event by event basis
TH2D *centroid = new TH2D("centroid", "s2/s1 vs. R_centroid per Event", 100, 0, rtop, 100, 0, 1.1);
centroid->GetXaxis()->CenterTitle();
centroid->GetXaxis()->SetTitle("r centroid in cm");
centroid->GetYaxis()->CenterTitle();
centroid->GetYaxis()->SetTitle("s1/s2");
}
//get data tree from file
TTree *Events = (TTree*)(fdata->Get("Events"));
EventData* evt = 0;
Events->SetBranchAddress("particle_data",&evt);
//Define
double S1=0;
double S2[nBins];
int Neff=0; //number of events used to average s2s1
double s2=0; //these are used for histogram of s1/s2 vs r_centroid
double ratio=0; //this is the ratio of s2/s1 used in centroid
double En[7][nBins]; //stores edep information as well as particle counts in radial bins
for(int q=0; q<7; q++) //fill En with zeros
{
if (q==0) {
for (int W=0; W<nBins; W++) {
En[q][W]=(W+1)*binSize;
}
}
if (q>0){
for (int w=0; w<nBins; w++) {
En[q][w]=0;
}
}
}
//Ensure no repeated events by only processing the number of events in the tree
if (Nev>Events->GetEntries()){
std::cout<<"Too Many Events Requested"<<endl<<"Changing Number Of Events To "<<Events->GetEntries()<<" "<<endl;
Nev=Events->GetEntries();
}
cout << "This tree contains " << Events->GetEntries() << " events; only print " << Nev << endl;
//loop through NeV events in the data tree
for(int event = 0; event < Nev; event++)
{
//print every once in a while, just so I can see that the code is running
if((event%5000)==0)cout << "At event " << event << "." << endl;
//get the next event
Events->GetEntry(event);
//define variables pertaining to this specific event and set them equal to zero
double r;
S1=0;
for (int k=0; k<nBins; k++) S2[k]=0;
double EdepActive=0, EdepTot=0;
double xCentroid=0, yCentroid=0, zCentroid=0;
double rCentroid=0;
ratio=0;
for(int j=0; j<evt->n_particles; j++) //loop through the events and get info from each
{
r=sqrt((evt->spawn_x.at(j))*(evt->spawn_x.at(j)) + (evt->spawn_y.at(j))*(evt->spawn_y.at(j)));
EdepActive+=evt->edep_active.at(j);
xCentroid+=evt->spawn_x.at(j)*evt->edep_active.at(j);
yCentroid+=evt->spawn_y.at(j)*evt->edep_active.at(j);
zCentroid+=evt->spawn_z.at(j)*evt->edep_active.at(j);
for (int k=0; k<nBins; k++) //loops through R in binSize incriments and places edep data in En and S2 arrays
{
if ((r>=rMin) && (r>=rMin+(k*binSize)) && (r<rMin+(k+1)*binSize)) {
En[1][k]+=evt->edep_active.at(j); //store energy deposited at this radius
En[2][k]+=1; //store number of events that have occured in dr thus far
}
if (r<rMin+(k+1)*binSize) {
En[4][k]+=evt->edep_active.at(j); //accumulates edeptotals as dr increases
En[5][k]+=1; //keeps track of the number of events that have occured by this radius
S2[k]+=evt->edep_active.at(j); //records S2 signal if ringsize was to this bin
}
}
}
if(EdepActive>0) { //calculate rCentroid
xCentroid /= EdepActive;
yCentroid /= EdepActive;
zCentroid /= EdepActive;
rCentroid=sqrt(xCentroid*xCentroid+yCentroid*yCentroid);
}
S1=S2[nBins-1]; //S1 is total edep in LAr chamber
if ((centYes==1) && (S1!=0) && (EdepActive>0)) {
ratio=S2[lessRings]/S1;
centroid->Fill(rCentroid, ratio);
}
if ((S1>0) && (sYes==1)) {
Neff+=1;
for (int z=0; z<nBins; z++) {
S2[z]=S2[z]/S1;
s2s1->Fill((rtop-(rMin+(z+.05)*binSize)), S2[z]);
}
}
}
if (sYes==1) s2s1->Scale(1./Neff);
if ((radYes==1) || (ravYes==1)) {
for (int J=0; J<=(nBins-1); J++) { //loop over J and produce the En[3][J] which is the average energy
//deposited in that bin for all events
if (radYes==1) En[3][J]=En[1][J]/Nev;
if (ravYes==1) En[6][J]=En[4][J]/Nev;
}
for (int n=0; n<=((printBins+rMin*overBin)-1); n++) {
if (radYes==1) radgen->Fill(n*binSize, En[3][n]);
if (ravYes==1) ravgen->Fill(n*binSize, En[6][n]);
}
}
TFile* out = new TFile("splotter.root", "recreate");
if (radYes==1) radgen->Write();
if (ravYes==1) ravgen->Write();
if (centYes==1)centroid->Write();
if (sYes==1) s2s1->Write();
if (radYes==1) {
if (gROOT->FindObject("c1")) delete gROOT->FindObject("c1");
TCanvas* c1=new TCanvas("c1");
radgen->SetTitle("Radial Energy Distribution");
radgen->SetLineColor(1);
radgen->Draw();
}
if (ravYes==1) {
if (gROOT->FindObject("c2")) delete gROOT->FindObject("c2");
TCanvas* c2=new TCanvas("c2");
ravgen->SetTitle("Radial Energy Distribution");
ravgen->SetLineColor(1);
ravgen->Draw();
}
if (sYes==1) {
if (gROOT->FindObject("c3")) delete gROOT->FindObject("c3");
TCanvas* c3=new TCanvas("c3");
s2s1->SetLineColor(1);
s2s1->GetXaxis()->SetRangeUser(0,5);
s2s1->Draw();
}
if (centYes==1) {
if (gROOT->FindObject("c4")) delete gROOT->FindObject("c4");
TCanvas* c4=new TCanvas("c4");
centroid->Draw();
}
}