icemc
 All Classes Namespaces Files Functions Variables Friends Macros Pages
secondaries.hh
Go to the documentation of this file.
1 
2 //class Secondaries:
4 #include <algorithm>
5 #include <numeric>
6 #include <iostream>
7 #include <string>
8 #include <sstream>
9 #include <fstream>
10 #include <iomanip>
11 #include <vector>
12 #include "Primaries.h"
13 class TH1F;
14 class Vector;
15 class Settings;
16 class TRandom3;
17 class Primaries;
18 
19 using std::vector;
20 using std::ifstream;
21 using std::string;
22 
24 class Secondaries {
25 
26 private:
27  TRandom3 Rand3;
28  void Picky(double *y_cumulative,int NPROB_MAX,double rnd,double& y);
29  // secondaries
30  static const int NPROB_MAX=100; // this sets the maximum length of the arrays
31  int NPROB; // this counts how many of those array elements are actually used
32  double plepton; // energy of charged lepton after interaction(s)
33  double em_secondaries_max; // em energy due to secondary interactions
34  double had_secondaries_max; // had energy due to secondary interactions
35  // first index=energy from 10^18 to 10^21 in half-decades
36  // second index=y (100 bins)
37  double dsdy_muon_brems[7][NPROB_MAX]; // probability distribution vs. y for brem
38  double dsdy_muon_epair[7][NPROB_MAX]; // prob. dist. vs. y for pair production
39  double dsdy_muon_pn[7][NPROB_MAX]; // photonuclear interactions
40 
41  double y_muon_brems[7][NPROB_MAX]; // inelasticity corresponding to each bin,brem
42  double y_muon_epair[7][NPROB_MAX]; // same for pair production
43  double y_muon_pn[7][NPROB_MAX]; // same for photonuclear
44 
45  vector<double> vy_muon_brems[7]; // vectors containing inelasticities distributed such
46  //that they follow the dsdy curves above.
47  vector<double> vy_muon_epair[7];
48  vector<double> vy_muon_pn[7];
49 
50  double y_cumulative_muon_brems[7][NPROB_MAX];
51  double y_cumulative_muon_epair[7][NPROB_MAX];
52  double y_cumulative_muon_pn[7][NPROB_MAX];
53 
54  double int_muon_brems[7]; // integral of prob. dist.=# of expected interactions.
55  // for each energy
56  double int_muon_epair[7]; // same for pair prod.
57  double int_muon_pn[7]; // same for photnuclear
58 
59  double max_muon_brems; // maximum value of prob. dist., brems
60  double max_muon_epair; // max value of prob. dist., pair prod.
61  double max_muon_pn; // same for photonucl.
62 
63  double min_muon_brems; // minimum value of prob. dist., brems
64  double min_muon_epair; // min value of prob. dist., pair prod.
65  double min_muon_pn; // same for photonucl.
66 
67  double dsdy_tauon_brems[7][NPROB_MAX]; // same as above, but for taus
68  double dsdy_tauon_epair[7][NPROB_MAX];
69  double dsdy_tauon_pn[7][NPROB_MAX];
70  double dsdy_tauon_hadrdecay[7][NPROB_MAX]; // hadronic decay of taus
71  double dsdy_tauon_edecay[7][NPROB_MAX]; // tau decay to electrons
72  double dsdy_tauon_mudecay[7][NPROB_MAX]; // tau decay to muons.
73 
74  double y_tauon_brems[7][NPROB_MAX];
75  double y_tauon_epair[7][NPROB_MAX];
76  double y_tauon_pn[7][NPROB_MAX];
77  double y_tauon_hadrdecay[7][NPROB_MAX];
78  double y_tauon_edecay[7][NPROB_MAX];
79  double y_tauon_mudecay[7][NPROB_MAX];
80 
81  double y_cumulative_tauon_brems[7][NPROB_MAX];
82  double y_cumulative_tauon_epair[7][NPROB_MAX];
83  double y_cumulative_tauon_pn[7][NPROB_MAX];
84  double y_cumulative_tauon_hadrdecay[7][NPROB_MAX];
85  double y_cumulative_tauon_edecay[7][NPROB_MAX];
86  double y_cumulative_tauon_mudecay[7][NPROB_MAX];
87 
88  vector<double> vy_tauon_brems[7]; // vectors containing inelasticities distributed such
89  //that they follow the dsdy curves above.
90  vector<double> vy_tauon_epair[7];
91  vector<double> vy_tauon_pn[7];
92  vector<double> vy_tauon_hadrdecay[7];
93  vector<double> vy_tauon_edecay[7];
94  vector<double> vy_tauon_mudecay[7];
95 
96  double int_tauon_brems[7];
97  double int_tauon_epair[7];
98  double int_tauon_pn[7];
99  double int_tauon_hadrdecay[7];
100  double int_tauon_edecay[7];
101  double int_tauon_mudecay[7];
102 
103  double max_tauon_brems;
104  double max_tauon_epair;
105  double max_tauon_pn;
106  double max_tauon_hadrdecay;
107  double max_tauon_edecay;
108  double max_tauon_mudecay;
109 
110  double min_tauon_brems; // minimum value of prob. dist., brems
111  double min_tauon_epair; // min value of prob. dist., pair prod.
112  double min_tauon_pn; // same for photonucl.
113  double min_tauon_hadrdecay; // minimum value of prob. dist., brems
114  double min_tauon_edecay; // min value of prob. dist., pair prod.
115  double min_tauon_mudecay; // same for photonucl.
116 
117 
118 // double weight_doublebang=0;
119 // double dtime_doublebang=0;
120 
121 
122  static const int N_TAUOLA=10001;
123  double tauola[N_TAUOLA][6];//tau decay energy fractions, [0]=tau neutrino, [1]=mu neutrino, [2]=Electron neutrino, [3]=hadrons, [4]=electromagnetic, [5]=muons
124  ifstream tauolainfile;
125 
126 
127  void readData(string,string,double (*)[NPROB_MAX], double (*)[NPROB_MAX]);
128  void ReadSecondaries(); // read in prob. dist. for secondary interactions
129 
130  //For Tau Weight & Survival probability equations.
131  //n.b. not in SI units.
132  //from Tau neutrino propagaiton and tau energy loss 2005 Dutta, Huang, & Reno.
133  //Equation 16 & used in Equation 30.
134  double B0,B1,E0;//parameterization using a logarithmic dependence on energy
135  //for B, the tau elecromagnetic energy loss parameter.
136  double mT;//mass of the Tau in Gev
137  double cT;//Tau Decay length in cm
138  double p;//Density of Standard Rock. g/cm^3
139 
140  //Used in Connolly Calc 2011.(d_dzPsurvNu())
141  //p, the Density of Standard Rock. g/cm^3
142  double A;
143  double Mn;// nucleon/ proton mass in grams,also equal to 0.938 GeV.
144 
145 
146 public:
147  Secondaries();
148 
149 
150  int SECONDARIES;
151  int TAUDECAY; // is tau decay counted as a secondary interaction
152  double TAUFRAC; //fraction of tau neutrino-cc current events where the primare interaction point is the first bang
153  int count_nfb;
154  int secondary_e_noncons;
155 
156  void GetSecondaries(Settings *settings1,string,double,double&,double&,int&,TH1F*);
157 
158  void InitTauola();
159  void GetTauDecay(string nuflavor,string current,string& taudecay, double& emfrac_db, double& hadfrac_db);
160 
161  void GetEMFracDB(double& emfrac_db, double& hadfrac_db);
162  double GetDBViewAngle(const Vector &refr, const Vector &nnu);
163  //void GetFirstBang(const Position &r_in, const Vector &nnu, Position &posnu, double len_int_kgm2, double d1, double &nuentrancelength);
164  double NFBWeight(double ptau, double taulength);
165 
166  int GetEMFrac(Settings *settings1,
167  string nuflavor,
168  string current,
169  string taudecay,
170  double y,
171  TH1F *hy,
172  double pnu,
173  int inu,
174 
175  double& emfrac,
176  double& hadfrac,
177  int& n_interactions, int taumodes1,double ptauf);
178 
179  bool secondbang;
180  static const bool interestedintaus=false;
181 
182  //string flavors[3]={"nue","numu","nutau"}; // the gps path of the anita-lite flight
183 string flavors[3]; // the gps path of the anita-lite flight
184 
185 }; //class Secondaries
186 
187 
Generated on Sat Nov 15 2014 22:18:21 for icemc by   doxygen 1.8.3.1