Skip to content
Snippets Groups Projects
Commit 1040c4af authored by Chad Lantz's avatar Chad Lantz
Browse files

Finally made things work as intended

parent 44159e9f
No related branches found
No related tags found
No related merge requests found
#include "TFile.h"
#include "TH1F.h"
#include "TGraph.h"
#include "TF1.h"
#include "TRandom.h"
#include <TTree.h>
......@@ -38,14 +39,14 @@ void cutZDCoutput(string fileName, string rod, double xOffset, double zOffset){
inTree->SetBranchAddress("Py",&Pyf);
inTree->SetBranchAddress("Pz",&Pzf);
//inTree->SetBranchAddress("EventNo",&EventNf);
//intree->SetBranchAddress("NCherenkovs",&NCherenkovs);
// intree->SetBranchAddress("NCherenkovs",&NCherenkovs);
//Output root file
TFile* outFile = new TFile( Form("output/Out_%s",fileName.c_str() ), "RECREATE");
TFile* outFile = new TFile( Form("output/%s_%s.root",fileName.substr(0,fileName.find_last_of(".")).c_str(), rod.c_str() ), "RECREATE");
TTree* outTree = new TTree( "tree", "tree" );
outTree->Branch("X0",&x_0);
outTree->Branch("Z0",&z_0);
// outTree->Branch("X0",&x_0);
// outTree->Branch("Z0",&z_0);
outTree->Branch("X",&x);
outTree->Branch("Z",&z);
outTree->Branch("Px",&px);
......@@ -53,11 +54,12 @@ void cutZDCoutput(string fileName, string rod, double xOffset, double zOffset){
outTree->Branch("Pz",&pz);
outTree->Branch("Energy",&Energy);
outTree->Branch("EventNo",&EventNo);
// outTree->SetBranchAddress("NCherenkovs",&NCherenkovs);
// output txt file
ofstream outputFile( Form("output/Out_%s.txt",fileName.substr(0,fileName.find_last_of("," ) ).c_str() ) );
ofstream outputFile( Form("output/%s_%s.txt",fileName.substr(0,fileName.find_last_of("." ) ).c_str(), rod.c_str() ) );
//// Get the transmission data
//// Get the tr.ansmission data
////
vector< double > transData;
string str;
......@@ -76,11 +78,14 @@ void cutZDCoutput(string fileName, string rod, double xOffset, double zOffset){
// Load that data into a histogram so we can use the GetBin function
TH1F* hTrans = new TH1F( "TransData", "TransData", nBins, 197.2888, 1027.23596);
TH1F* hTrans2 = new TH1F( "TransData2", "TransData2", nBins, 197.2888, 1027.23596);
float trans, rand;
for(int bin = 0; bin < nBins; bin++){
hTrans->SetBinContent( bin, transData[bin] );
}
// Apply a modified baseline to the fused quartz spectrum
if(rod == "fusedQuartz"){
TF1* bLine = new TF1("baseline","0.0883*log(0.013087*(x-156.8) ) + 0.802",197,1028);
hTrans->Multiply(bLine);
......@@ -115,13 +120,14 @@ void cutZDCoutput(string fileName, string rod, double xOffset, double zOffset){
// get the data from the appropriate bin.
wl = 1240./energy;
if(wl < 197.2888){
trans = hTrans->GetBinContent(1);
trans = hTrans->GetBinContent(10);
}else{
trans = hTrans->GetBinContent( hTrans->GetBin( wl ) );
trans = hTrans->GetBinContent( hTrans->FindBin( wl ) );
}
// Give the photon a random chance to make it weighted
// by transmission %
if( rnd.Rndm() < trans ){
// if( rnd.Rndm() < trans ){
if( rnd.Rndm() < 1 ){
// If the photon is transmitted, add it to the vector
X = Xf->at(k) + xOffset;
......@@ -214,17 +220,15 @@ void PMTcuts(string fileName, int PMTmodel = 6091){
}
}
// Load that data into a histogram so we can use the GetBin function
TH1F* hQE = new TH1F( "TransData", "TransData", wavelength.size(), wavelength.front(), wavelength.back());
for(int bin = 0; bin < wavelength.size(); bin++){
hQE->SetBinContent(bin,efficiency[bin]);
}
TVectorD TvWL(wavelength.size(),&wavelength[0]);
TVectorD TvQE(efficiency.size(),&efficiency[0]);
TGraph* gQE = new TGraph(TvWL,TvQE);
//// begin loop over events
////
TRandom rnd;
int nCut, nPhotons, totalCut=0, totalPhotons=0;
double X,Z,wl,trans;
double X,Z,wl,trans,rand;
int nEvents = inTree->GetEntries();
for (int ev = 0; ev < nEvents ; ev++){
inTree->GetEntry(ev);
......@@ -243,9 +247,10 @@ void PMTcuts(string fileName, int PMTmodel = 6091){
nCut++;
} else{
// Get transmission % from data
trans = hQE->GetBinContent( hQE->GetBin( wl ) );
trans = gQE->Eval( wl );
// Give the photon a random chance to make it past the
// Quantum Efficiency check
rand = 100*rnd.Rndm();
if( 100*rnd.Rndm() > trans ){
nCut++;
}
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment