GRSISort "v4.0.0.5"
An extension of the ROOT analysis Framework
Loading...
Searching...
No Matches
MakeGriffinHistograms.cxx
Go to the documentation of this file.
1#include "TRuntimeObjects.h"
2
3#include "TRunInfo.h"
4
5#include "TGriffin.h"
6#include "TGriffinBgo.h"
7#include "TSceptar.h"
8#include "TZeroDegree.h"
9#include "TDescant.h"
10
11const Double_t ps = 1E-3;
12const Double_t ns = 1.;
13const Double_t us = 1E3;
14const Double_t ms = 1E6;
15const Double_t s = 1E9;
16
17const Double_t keV = 1.;
18const Double_t MeV = 1E3;
19
20const Double_t ggTLow = -200 * ns;
21const Double_t ggTHigh = 200 * ns;
22const Double_t gbgoTLow = -250 * ns;
23const Double_t gbgoTHigh = 250 * ns;
24const Double_t gbTLow = -150 * ns;
25const Double_t gbTHigh = 300 * ns;
26const Double_t gzTLow = -50 * ns;
27const Double_t gzTHigh = 50 * ns;
28
29const Double_t bELow = 50 * keV;
30const Double_t zELow = 200 * keV;
31const ULong64_t kCycleLength = 34.5 * s;
32
34{
35 return (((hit2->GetTime() * ns - hit1->GetTime() * ns) >= ggTLow) && ((hit2->GetTime() * ns - hit1->GetTime() * ns) <= ggTHigh));
36}
37
38bool PromptCoincidence(TBgoHit* bgo_hit, TGriffinHit* grif_hit)
39{
40 return (((grif_hit->GetTime() * ns - bgo_hit->GetTime() * ns) >= gbgoTLow) && ((grif_hit->GetTime() * ns - bgo_hit->GetTime() * ns) <= gbgoTHigh));
41}
42
44{
45 return (((grif_hit->GetTime() * ns - sc_hit->GetTime() * ns) >= gbTLow) && ((grif_hit->GetTime() * ns - sc_hit->GetTime() * ns) <= gbTHigh) && (sc_hit->GetEnergy() > bELow));
46}
47
49{
50 return (((grif_hit->GetTime() * ns - sc_hit->GetTime() * ns) >= gzTLow) && ((grif_hit->GetTime() * ns - sc_hit->GetTime() * ns) <= gzTHigh) && (sc_hit->GetEnergy() > zELow) && (sc_hit->GetDetector() == 1));
51}
52
54{
55 auto grif = obj.GetDetector<TGriffin>();
56 auto grif_bgo = obj.GetDetector<TGriffinBgo>().get();
57 auto scep = obj.GetDetector<TSceptar>();
58 auto zds = obj.GetDetector<TZeroDegree>();
59 auto descant = obj.GetDetector<TDescant>();
60
61 static int kRunNumber = 0;
62 if(kRunNumber == 0) {
63 gROOT->cd();
64 kRunNumber = TRunInfo::RunNumber();
65 }
66
67 static ULong64_t less_background = 0.5 * s;
68
69 if(kRunNumber >= 12012) {
70 less_background = 0.0 * s;
71 } else {
72 less_background = 0.5 * s;
73 }
74
75 // static TPPG* ppg = TPPG::Get();
76 // static const ULong64_t kCycleLength = ppg->OdbCycleLength();//in us
77
78 if(grif != nullptr) {
79 for(auto g1 = 0; g1 < grif->GetMultiplicity(); ++g1) {
80 obj.FillHistogram("GRIFFIN", "gE_Cyc", 5000, 0, 5000, grif->GetGriffinHit(g1)->GetEnergy(), kCycleLength * 20 / s, 0, kCycleLength / s, (((ULong64_t)(grif->GetHit(g1)->GetTimeStamp() * ns)) % (kCycleLength - less_background) + less_background) / s);
81 if(scep != nullptr) {
82 bool beta_found = false;
83 for(auto b1 = 0; b1 < scep->GetMultiplicity(); ++b1) {
84 obj.FillHistogram("TimeDiffs", "gbT", 2000, -2000, 2000, grif->GetHit(g1)->GetTime() * ns - scep->GetHit(b1)->GetTime() * ns, 20, 1, 21, scep->GetHit(b1)->GetDetector());
85 if(PromptCoincidence(scep->GetSceptarHit(b1), grif->GetGriffinHit(g1)) && !beta_found) {
86 obj.FillHistogram("GRIFFIN", "gEB", 8000, 0, 8000, grif->GetHit(g1)->GetEnergy());
87 beta_found = true;
88 }
89 }
90 }
91 if(zds != nullptr) {
92 bool beta_found = false;
93 for(auto z1 = 0; z1 < zds->GetMultiplicity(); ++z1) {
94 obj.FillHistogram("TimeDiffs", "gzT", 2000, -2000, 2000, grif->GetHit(g1)->GetTime() * ns - zds->GetHit(z1)->GetTime() * ns, 2, 1, 3, zds->GetHit(z1)->GetDetector());
95 if(PromptCoincidence(zds->GetZeroDegreeHit(z1), grif->GetGriffinHit(g1)) && !beta_found) {
96 obj.FillHistogram("GRIFFIN", "gEZ", 8000, 0, 8000, grif->GetHit(g1)->GetEnergy());
97 beta_found = true;
98 }
99 }
100 }
101 for(auto g2 = g1 + 1; g2 < grif->GetMultiplicity(); ++g2) {
102 obj.FillHistogram("TimeDiffs", "ggT", 2000, -2000, 2000, grif->GetHit(g2)->GetTime() * ns - grif->GetHit(g1)->GetTime() * ns);
103 if(PromptCoincidence(grif->GetGriffinHit(g1), grif->GetGriffinHit(g2))) {
104 obj.FillHistogram("GRIFFIN", "ggE", 5000, 0, 5000, grif->GetHit(g1)->GetEnergy(), 5000, 0, 5000, grif->GetHit(g2)->GetEnergy());
105 obj.FillHistogram("GRIFFIN", "ggE", 5000, 0, 5000, grif->GetHit(g2)->GetEnergy(), 5000, 0, 5000, grif->GetHit(g1)->GetEnergy());
106 obj.FillHistogram("HitPatterns", "gg_HP", 65, 0, 65, grif->GetGriffinHit(g1)->GetChannel()->GetNumber(), 65, 0, 65, grif->GetGriffinHit(g2)->GetChannel()->GetNumber());
107 }
108 }
109 if(grif_bgo) {
110 for(auto bgo1 = 0; bgo1 < grif_bgo->GetMultiplicity(); ++bgo1) {
111 obj.FillHistogram("TimeDiffs", "gbgoT", 2000, -2000, 2000, grif_bgo->GetHit(bgo1)->GetTime() * ns - grif->GetHit(g1)->GetTime() * ns, 64 * 5, 1, 64 * 5 + 1, grif_bgo->GetHit(bgo1)->GetArrayNumber());
112 obj.FillHistogram("HitPatterns", "ggbgo_HP", 65, 0, 65, grif->GetGriffinHit(g1)->GetArrayNumber(), 64 * 5, 1, 64 * 5 + 1, grif_bgo->GetBgoHit(bgo1)->GetArrayNumber());
113
114 obj.FillHistogram("HitPatterns", "gbgoDet", 17, 0, 17, grif->GetGriffinHit(g1)->GetDetector(), 17, 0, 17, grif_bgo->GetHit(bgo1)->GetDetector());
115 if((grif_bgo->GetHit(bgo1)->GetEnergy() > 50) && PromptCoincidence(grif_bgo->GetBgoHit(bgo1), grif->GetGriffinHit(g1))) {
116 obj.FillHistogram("HitPatterns", "gbgoDet_ETCut", 17, 0, 17, grif->GetGriffinHit(g1)->GetDetector(), 17, 0, 17, grif_bgo->GetHit(bgo1)->GetDetector());
117 }
118 if(grif->GetHit(g1)->GetDetector() == grif_bgo->GetHit(bgo1)->GetDetector()) {
119 obj.FillHistogram("GRIFFIN", "bgoEgE", 2000, 0, 2000, grif->GetHit(g1)->GetEnergy(), 2000, 0, 2000, grif_bgo->GetHit(bgo1)->GetEnergy());
120 obj.FillHistogram("GRIFFIN", "bgoTgE", 2000, 0, 2000, grif->GetHit(g1)->GetEnergy(), 2000, -1000, 1000, grif_bgo->GetHit(bgo1)->GetTime() * ns - grif->GetHit(g1)->GetTime() * ns);
121 if(std::fabs(grif_bgo->GetHit(bgo1)->GetTime() * ns - grif->GetHit(g1)->GetTime() * ns) < 300 * ns) {
122 obj.FillHistogram("GRIFFIN", "bgoEgE_TCut", 2000, 0, 2000, grif->GetHit(g1)->GetEnergy(), 2000, 0, 2000, grif_bgo->GetHit(bgo1)->GetEnergy());
123 }
124 if(grif_bgo->GetHit(bgo1)->GetEnergy() * keV > 50 * keV) {
125 obj.FillHistogram("GRIFFIN", "bgoTgE_ECut", 2000, 0, 2000, grif->GetHit(g1)->GetEnergy(), 2000, -1000, 1000, grif_bgo->GetHit(bgo1)->GetTime() * ns - grif->GetHit(g1)->GetTime() * ns);
126 }
127 }
128 }
129 }
130 }
131 for(auto a1 = 0; a1 < grif->GetAddbackMultiplicity(); ++a1) {
132 obj.FillHistogram("GRIFFIN", "aE", 5000, 0, 5000, grif->GetAddbackHit(a1)->GetEnergy());
133 obj.FillHistogram("GRIFFIN", "aE_Cyc", 5000, 0, 5000, grif->GetAddbackHit(a1)->GetEnergy(), kCycleLength * 20 / s, 0, kCycleLength / s, (((ULong64_t)(grif->GetHit(a1)->GetTimeStamp() * ns)) % (kCycleLength - less_background) + less_background) / s);
134 if(zds != nullptr) {
135 bool beta_found = false;
136 for(auto z1 = 0; z1 < zds->GetMultiplicity(); ++z1) {
137 if(PromptCoincidence(zds->GetZeroDegreeHit(z1), grif->GetAddbackHit(a1)) && !beta_found) {
138 obj.FillHistogram("GRIFFIN", "aEZ_Cyc", 8000, 0, 8000, grif->GetAddbackHit(a1)->GetEnergy(), kCycleLength * 20 / s, 0, kCycleLength / s, (((ULong64_t)(grif->GetHit(a1)->GetTimeStamp() * ns)) % (kCycleLength - less_background) + less_background) / s);
139 beta_found = true;
140 }
141 }
142 }
143 }
144 for(auto as1 = 0; as1 < grif->GetSuppressedAddbackMultiplicity(grif_bgo); ++as1) {
145 obj.FillHistogram("GRIFFIN", "aES", 8000, 0, 8000, grif->GetSuppressedAddbackHit(as1)->GetEnergy());
146 bool zds_or_scep_found = false;
147 if(scep != nullptr) {
148 bool beta_found = false;
149 for(auto b1 = 0; b1 < scep->GetMultiplicity(); ++b1) {
150 if(PromptCoincidence(scep->GetSceptarHit(b1), grif->GetSuppressedAddbackHit(as1)) && !beta_found) {
151 obj.FillHistogram("GRIFFIN", "aESB_Cyc", 8000, 0, 8000, grif->GetSuppressedAddbackHit(as1)->GetEnergy(), kCycleLength * 20 / s, 0, kCycleLength / s, (((ULong64_t)(grif->GetSuppressedAddbackHit(as1)->GetTimeStamp() * ns)) % (kCycleLength - less_background) + less_background) / s);
152 beta_found = true;
153
154 if(!zds_or_scep_found) {
155 obj.FillHistogram("GRIFFIN", "aESbeta_Cyc", 8000, 0, 8000, grif->GetSuppressedAddbackHit(as1)->GetEnergy(), kCycleLength * 20 / s, 0, kCycleLength / s, (((ULong64_t)(grif->GetSuppressedAddbackHit(as1)->GetTimeStamp() * ns)) % (kCycleLength - less_background) + less_background) / s);
156 zds_or_scep_found = true;
157 }
158 }
159 }
160 }
161 if(zds != nullptr) {
162 bool beta_found = false;
163 for(auto z1 = 0; z1 < zds->GetMultiplicity(); ++z1) {
164 if(PromptCoincidence(zds->GetZeroDegreeHit(z1), grif->GetSuppressedAddbackHit(as1)) && !beta_found) {
165 obj.FillHistogram("GRIFFIN", "aESZ_Cyc", 8000, 0, 8000, grif->GetSuppressedAddbackHit(as1)->GetEnergy(), kCycleLength * 20 / s, 0, kCycleLength / s, (((ULong64_t)(grif->GetSuppressedAddbackHit(as1)->GetTimeStamp() * ns)) % (kCycleLength - less_background) + less_background) / s);
166 beta_found = true;
167 if(!zds_or_scep_found) {
168 obj.FillHistogram("GRIFFIN", "aESbeta_Cyc", 8000, 0, 8000, grif->GetSuppressedAddbackHit(as1)->GetEnergy(), kCycleLength * 20 / s, 0, kCycleLength / s, (((ULong64_t)(grif->GetSuppressedAddbackHit(as1)->GetTimeStamp() * ns)) % (kCycleLength - less_background) + less_background) / s);
169 zds_or_scep_found = true;
170 }
171 }
172 }
173 }
174 for(auto as2 = as1 + 1; as2 < grif->GetSuppressedAddbackMultiplicity(grif_bgo); ++as2) {
175 obj.FillHistogram("GRIFFIN", "aaES", 8000, 0, 8000, grif->GetSuppressedAddbackHit(as1)->GetEnergy(), 8000, 0, 8000, grif->GetSuppressedAddbackHit(as2)->GetEnergy());
176 obj.FillHistogram("GRIFFIN", "aaES", 8000, 0, 8000, grif->GetSuppressedAddbackHit(as2)->GetEnergy(), 8000, 0, 8000, grif->GetSuppressedAddbackHit(as1)->GetEnergy());
177
178 bool combined_beta_found = false;
179 if(scep != nullptr) {
180 bool beta_found = false;
181 for(auto b1 = 0; b1 < scep->GetMultiplicity(); ++b1) {
182 if((PromptCoincidence(scep->GetSceptarHit(b1), grif->GetSuppressedAddbackHit(as1)) || PromptCoincidence(scep->GetSceptarHit(b1), grif->GetSuppressedAddbackHit(as2))) && !beta_found) {
183 obj.FillHistogram("GRIFFIN", "aaESb", 8000, 0, 8000, grif->GetSuppressedAddbackHit(as2)->GetEnergy(), 8000, 0, 8000, grif->GetSuppressedAddbackHit(as1)->GetEnergy());
184 obj.FillHistogram("GRIFFIN", "aaESb", 8000, 0, 8000, grif->GetSuppressedAddbackHit(as1)->GetEnergy(), 8000, 0, 8000, grif->GetSuppressedAddbackHit(as2)->GetEnergy());
185 beta_found = true;
186
187 if(!combined_beta_found) {
188 obj.FillHistogram("GRIFFIN", "aaESbeta", 8000, 0, 8000, grif->GetSuppressedAddbackHit(as2)->GetEnergy(), 8000, 0, 8000, grif->GetSuppressedAddbackHit(as1)->GetEnergy());
189 obj.FillHistogram("GRIFFIN", "aaESbeta", 8000, 0, 8000, grif->GetSuppressedAddbackHit(as1)->GetEnergy(), 8000, 0, 8000, grif->GetSuppressedAddbackHit(as2)->GetEnergy());
190 combined_beta_found = true;
191 }
192 }
193 }
194 }
195 if(zds != nullptr) {
196 bool beta_found = false;
197 for(auto z1 = 0; z1 < zds->GetMultiplicity(); ++z1) {
198 if((PromptCoincidence(zds->GetZeroDegreeHit(z1), grif->GetSuppressedAddbackHit(as1)) || PromptCoincidence(zds->GetZeroDegreeHit(z1), grif->GetSuppressedAddbackHit(as2))) && !beta_found) {
199 obj.FillHistogram("GRIFFIN", "aaESz", 8000, 0, 8000, grif->GetSuppressedAddbackHit(as2)->GetEnergy(), 8000, 0, 8000, grif->GetSuppressedAddbackHit(as1)->GetEnergy());
200 obj.FillHistogram("GRIFFIN", "aaESz", 8000, 0, 8000, grif->GetSuppressedAddbackHit(as1)->GetEnergy(), 8000, 0, 8000, grif->GetSuppressedAddbackHit(as2)->GetEnergy());
201 beta_found = true;
202 if(!combined_beta_found) {
203 obj.FillHistogram("GRIFFIN", "aaESbeta", 8000, 0, 8000, grif->GetSuppressedAddbackHit(as2)->GetEnergy(), 8000, 0, 8000, grif->GetSuppressedAddbackHit(as1)->GetEnergy());
204 obj.FillHistogram("GRIFFIN", "aaESbeta", 8000, 0, 8000, grif->GetSuppressedAddbackHit(as1)->GetEnergy(), 8000, 0, 8000, grif->GetSuppressedAddbackHit(as2)->GetEnergy());
205 combined_beta_found = true;
206 }
207 }
208 }
209 }
210 }
211 }
212
213 for(auto sg1 = 0; sg1 < grif->GetSuppressedMultiplicity(grif_bgo); ++sg1) {
214 obj.FillHistogram("HitPatterns", "sMultgMult", 100, 0, 100, grif->GetMultiplicity(), 100, 0, 100, grif->GetSuppressedMultiplicity(grif_bgo));
215 obj.FillHistogram("GRIFFIN", "gES", 8000, 0, 8000, grif->GetSuppressedHit(sg1)->GetEnergy(), 65, 0, 65, grif->GetSuppressedHit(sg1)->GetArrayNumber());
216 for(auto sg2 = sg1 + 1; sg2 < grif->GetSuppressedMultiplicity(grif_bgo); ++sg2) {
217 obj.FillHistogram("GRIFFIN", "ggES", 8000, 0, 8000, grif->GetSuppressedHit(sg1)->GetEnergy(), 8000, 0, 8000, grif->GetSuppressedHit(sg2)->GetEnergy());
218 obj.FillHistogram("GRIFFIN", "ggES", 8000, 0, 8000, grif->GetSuppressedHit(sg2)->GetEnergy(), 8000, 0, 8000, grif->GetSuppressedHit(sg1)->GetEnergy());
219 }
220
221 if(scep != nullptr) {
222 bool beta_found = false;
223 for(auto b1 = 0; b1 < scep->GetMultiplicity(); ++b1) {
224 if(PromptCoincidence(scep->GetSceptarHit(b1), grif->GetSuppressedHit(sg1)) && !beta_found) {
225 obj.FillHistogram("GRIFFIN", "gEBS", 8000, 0, 8000, grif->GetSuppressedHit(sg1)->GetEnergy());
226 beta_found = true;
227 }
228 }
229 }
230 if(grif_bgo) {
231 for(auto bgo1 = 0; bgo1 < grif_bgo->GetMultiplicity(); ++bgo1) {
232 obj.FillHistogram("HitPatterns", "ggbgoS_HP", 65, 0, 65, grif->GetSuppressedHit(sg1)->GetArrayNumber(), 64 * 5, 1, 64 * 5 + 1, grif_bgo->GetBgoHit(bgo1)->GetArrayNumber());
233 }
234 }
235 }
236 }
237 if(grif_bgo != nullptr) {
238 for(auto bgo1 = 0; bgo1 < grif_bgo->GetMultiplicity(); ++bgo1) {
239 obj.FillHistogram("BGO", "bgoE", 3000, 0, 3000, grif_bgo->GetHit(bgo1)->GetEnergy(), 64 * 5, 1, 64 * 5 + 1, grif_bgo->GetBgoHit(bgo1)->GetArrayNumber());
240 obj.FillHistogram("HitPatterns", "bgoHP", 20 * 16 + 1, 0, 20 * 16 + 1, grif_bgo->GetHit(bgo1)->GetArrayNumber());
241 }
242 }
243 if(zds != nullptr) {
244 for(auto z1 = 0; z1 < zds->GetMultiplicity(); ++z1) {
245 obj.FillHistogram("ZDS", "ZE", 2500, 0, 65000, zds->GetHit(z1)->GetEnergy(), 2, 1, 3, zds->GetHit(z1)->GetDetector());
246 obj.FillHistogram("ZDS", "z_Cyc", 1000, 0, 10000, zds->GetHit(z1)->GetEnergy(), kCycleLength * 20 / s, 0, kCycleLength / s, (((ULong64_t)(zds->GetHit(z1)->GetTimeStamp() * ns)) % (kCycleLength - less_background) + less_background) / s);
247 }
248 }
249 if(scep != nullptr) {
250 for(auto b1 = 0; b1 < scep->GetMultiplicity(); ++b1) {
251 obj.FillHistogram("SCEPTAR", "bE", 10000, 0, 6e7, scep->GetHit(b1)->GetEnergy(), 20, 1, 21, scep->GetHit(b1)->GetDetector());
252 obj.FillHistogram("HitPatterns", "bHP", 21, 0, 21, scep->GetHit(b1)->GetDetector());
253 obj.FillHistogram("SCEPTAR", "b_Cyc", kCycleLength * 20 / s, 0, kCycleLength / s, (((ULong64_t)(scep->GetHit(b1)->GetTimeStamp() * ns)) % (kCycleLength - less_background) + less_background) / s);
254 }
255 }
256 if(descant != nullptr) {
257 for(auto d1 = 0; d1 < descant->GetMultiplicity(); ++d1) {
258 obj.FillHistogram("HitPatterns", "dHP", 10, 0, 10, descant->GetHit(d1)->GetDetector());
259 }
260 }
261 if(grif != nullptr) {
262 for(auto g1 = 0; g1 < grif->GetLowGainMultiplicity(); ++g1) {
263 obj.FillHistogram("HitPatterns", "gHP_lg", 65, 0, 65, grif->GetGriffinLowGainHit(g1)->GetArrayNumber());
264 obj.FillHistogram("GRIFFIN", "gE_lg", 5000, 0, 5000, grif->GetGriffinLowGainHit(g1)->GetEnergy(), 65, 0, 65, grif->GetGriffinLowGainHit(g1)->GetArrayNumber());
265 }
266 for(auto g1 = 0; g1 < grif->GetHighGainMultiplicity(); ++g1) {
267 obj.FillHistogram("HitPatterns", "gHP_hg", 65, 0, 65, grif->GetGriffinHighGainHit(g1)->GetArrayNumber());
268 obj.FillHistogram("GRIFFIN", "gE_hg", 5000, 0, 5000, grif->GetGriffinHighGainHit(g1)->GetEnergy(), 65, 0, 65, grif->GetGriffinHighGainHit(g1)->GetArrayNumber());
269 }
270 }
271}
const Double_t gbTLow
const Double_t keV
const Double_t us
const Double_t bELow
const Double_t ps
const Double_t ggTHigh
const Double_t gbTHigh
const Double_t ns
const Double_t gzTHigh
void MakeAnalysisHistograms(TRuntimeObjects &obj)
const Double_t ms
const Double_t gbgoTLow
const ULong64_t kCycleLength
const Double_t zELow
bool PromptCoincidence(TGriffinHit *hit1, TGriffinHit *hit2)
const Double_t s
const Double_t gbgoTHigh
const Double_t MeV
const Double_t gzTLow
const Double_t ggTLow
virtual Long64_t GetTimeStamp(Option_t *="") const
virtual double GetEnergy(Option_t *opt="") const
virtual UShort_t GetArrayNumber() const
! Simply returns the detector number, overwritten for detectors that have crystals/segments
virtual Int_t GetDetector() const
!
virtual Double_t GetTime(const ETimeFlag &correct_flag=ETimeFlag::kAll, Option_t *opt="") const
Returns a time value to the nearest nanosecond!
virtual Short_t GetMultiplicity() const
Definition TDetector.h:73
virtual TDetectorHit * GetHit(const int &index) const
Definition TDetector.cxx:61
TDetectorHit * GetHit(const int &idx)
Definition TGriffin.cxx:375
Short_t GetMultiplicity() const override
Definition TGriffin.h:52
TGriffinHit * GetGriffinHit(const int &i)
!
Definition TGriffin.h:47
static int RunNumber()
Definition TRunInfo.h:201
Object passed to the online histograms.
TH1 * FillHistogram(const char *name, int bins, double low, double high, double value, double weight=1)
std::shared_ptr< T > GetDetector()
Returns a pointer to the detector of type T.