fragmentYields.C
Go to the documentation of this file.
1 #include "Riostream.h"
2 #include "TSystem.h"
3 #include "TInterpreter.h"
4 #include "TROOT.h"
5 #include "TApplication.h"
6 #include "TFile.h"
7 #include "TNtuple.h"
8 #include "TCanvas.h"
9 #include "TH1F.h"
10 #include "THStack.h"
11 #include "TCut.h"
12 #include "TString.h"
13 #include "TMath.h"
14 
15 /***************************
16  * Root script that produces the yield calculations for fragments
17  * below the angle of 10 degrees.
18  * This is the same numbers E.Haettner has calculated
19  * However, there is room for much improvement to be made here.
20  *
21  * @author Gillis Danielsen
22  * **************************/
23 
24 
26 
27 
28 // gROOT->SetStyle("clearRetro");
29  //this will be used as base for pulling the experimental data
30  TString dir = gSystem->UnixPathName(gInterpreter->GetCurrentMacroName());
31  dir.ReplaceAll("fragmentAngularDistribution.C","");
32  dir.ReplaceAll("/./","/");
33  ifstream in;
34 
35  TString pDepth, fragment, Znum, normToOneAtZeroAngle;
36  cout << "Enter phantom depth (eg. 27.9, see experimentalData directory for choices): ";
37  cin >> pDepth;
38 
39  TString simulationDataPath = "IAEA_" + pDepth + ".root";
40 
41  //Let's pull in the simulation-data
42  //TFile *MCData = TFile::Open("IAEA_" + pDepth + ".root");
43  TFile *MCData = TFile::Open(simulationDataPath);
44  TNtuple *fragments = (TNtuple*) MCData->Get("fragmentNtuple");
45 
46  //Block bellow pulls out the simulation's metadata from the metadata ntuple.
47  TNtuple *metadata = (TNtuple*) MCData->Get("metaData");
48  Float_t events, detectorDistance,waterThickness,beamEnergy,energyError,phantomCenterDistance;
49  metadata->SetBranchAddress("events",&events);
50  metadata->SetBranchAddress("waterThickness",&waterThickness);
51  metadata->SetBranchAddress("detectorDistance",&detectorDistance);
52  metadata->SetBranchAddress("beamEnergy",&beamEnergy);
53  metadata->SetBranchAddress("energyError",&energyError);
54  metadata->SetBranchAddress("phantomCenterDistance",&phantomCenterDistance);
55  metadata->GetEntry(0); //there is just one row to consider.
56 
57  //good to keep for ref. G4 might give weird units due to change.
58  metadata->Scan();
59 
60  //ALL UNITS ARE cm!
61  Double_t scatteringDistance = detectorDistance - phantomCenterDistance; //temporarily hard-coded, should be distance from target-center to detector
62 
63  Double_t r;
64  Double_t degrees = 10.0;
65  Double_t rMin;
66  Double_t rMax;
67  TString rMinString;
68  TString rMaxString;
69 
70  rMinString = "0.00";
71  rMaxString = Form("%f", scatteringDistance*TMath::ATan(degrees*TMath::DegToRad()));
72  std::cout << "Computes comparable figures to E.haettners yield calculations for zero to ten degree angles." << endl;
73  std::cout << "Firsta value is yield of Z=1 per C12 incident paritcle second is for Z=2 etc." << endl;
74  std::cout << "verboseness has been limited to make pasting easyer." << endl;
75  for(int i = 1; i <= 6; i++){
76  int numEntries = fragments->GetEntries("(Z == " + TString::Format("%i",i) + " && sqrt(posY^2 + posZ^2) < " + rMaxString + "&& sqrt(posY*posY + posZ*posZ) > " + rMinString + ")");
77  std::cout << numEntries / (events * .999) << endl; //< for calculation of FWHM
78  }
79 }
void fragmentYields()
OStream cout
Definition: OStream.cxx:6
ifstream in
Definition: comparison.C:7
TDirectory * dir
Definition: macro.C:5
TRandom3 r(0)
void Format(TGraph *gr, int lcol, int lsty, int lwid, int mcol, int msty, double msiz)
Definition: Style.cxx:154
void events(int which)
Definition: Cana.C:52