GRSISort "v4.0.0.5"
An extension of the ROOT analysis Framework
Loading...
Searching...
No Matches
TGriffinAngles.cxx
Go to the documentation of this file.
1#include "TGriffinAngles.h"
2#include "TGriffin.h"
3#include "TGRSIOptions.h"
4
5double TGriffinAngles::fRounding = 0.001;
7
8TGriffinAngles::TGriffinAngles(double distance, bool folding, bool grouping, bool addback)
9 : fDistance(distance), fFolding(folding), fGrouping(grouping), fAddback(addback)
10{
11 SetName("GriffinAngles");
12 // get user settings for excluded detectors/crystals, and custom grouping
13 if(TGRSIOptions::Get() != nullptr) {
14 auto* settings = TGRSIOptions::UserSettings();
15 if(settings != nullptr) {
16 // try quietly to get the vectors of excluded crystals and detectors, catching (and disregarding) any exceptions
17 try {
18 fExcludedCrystals = settings->GetIntVector("ExcludedCrystal", true);
19 } catch(std::out_of_range&) {}
20 try {
21 fExcludedDetectors = settings->GetIntVector("ExcludedDetector", true);
22 } catch(std::out_of_range&) {}
23 try {
24 fCustomGrouping = settings->GetIntVector("CustomGrouping", true);
25 std::sort(fCustomGrouping.begin(), fCustomGrouping.end());
26 } catch(std::out_of_range&) {}
27 } else if(fVerbosity > EVerbosity::kQuiet) {
28 std::cout << "Failed to find user settings in TGRSIOptions, can't get user settings for excluded detectors/crystals" << std::endl;
29 }
30 } else if(fVerbosity > EVerbosity::kQuiet) {
31 std::cout << "Failed to find TGRSIOptions, can't get user settings for excluded detectors/crystals" << std::endl;
32 }
33
34 // check that we only use grouping if folding is also enabled, and we do not have any custom grouping information
35 if(fGrouping && !fFolding && fCustomGrouping.empty()) {
36 std::cout << DRED << "Warning, grouping is only possible if folding is also active. Setting folding to true!" << RESET_COLOR << std::endl;
37 fFolding = true;
38 }
39
41 std::cout << "Creating angles for detectors " << fDistance << " mm from center of array, " << (fAddback ? "" : "not ") << "using addback, " << (fFolding ? "" : "not ") << "using folding, and " << (fGrouping ? "" : "not ") << "using grouping. Any angles less than " << fRounding << " degrees apart will be considered the same." << std::endl;
42 }
43
44 // loop over all possible detector/crystal combinations
45 for(int firstDet = 1; firstDet <= 16; ++firstDet) {
46 if(ExcludeDetector(firstDet)) { continue; }
47 for(int firstCry = 0; firstCry < 4; ++firstCry) {
48 if(ExcludeCrystal(firstDet, firstCry)) { continue; }
49 for(int secondDet = 1; secondDet <= 16; ++secondDet) {
50 if(ExcludeDetector(secondDet)) { continue; }
51 for(int secondCry = 0; secondCry < 4; ++secondCry) {
52 if(ExcludeCrystal(secondDet, secondCry)) { continue; }
53 // exclude hits in the same crystal or, if addback is enabled, in the same detector
54 if(firstDet == secondDet && (firstCry == secondCry || fAddback)) { continue; }
55 double angle = TGriffin::GetPosition(firstDet, firstCry, fDistance).Angle(TGriffin::GetPosition(secondDet, secondCry, fDistance)) * 180. / TMath::Pi();
56 if(fVerbosity == kAll) {
57 std::cout << "det./cry. " << firstDet << "/" << firstCry << " with " << secondDet << "/" << secondCry << ", at " << fDistance << " mm = " << angle << std::endl;
58 }
59 // if folding is enable we fold the distribution at 90 degree and only use angles between 0 and 90 degree
60 if(fFolding && angle > 90.) {
61 angle = 180. - angle;
62 }
63 // round down to zero (mainly here for nicer looks)
64 if(angle < fRounding) {
65 angle = 0.;
66 }
67 // if the lower bound and the upper bound are the same we have a new angle
68 if(fAngles.lower_bound(angle - fRounding) == fAngles.upper_bound(angle + fRounding)) {
69 fAngles.insert(angle);
70 }
71 // this should always work, either this is a new angle, in which case it gets initialized to zero and then incremented to one,
72 // or we increment the existing counter
73 // the key is integer, so by dividing by rounding and then casting to integer we can avoid duplicates close to each other
74 // factor 2 to include that the "normal" rounding is +- fRounding
75 fAngleCount[static_cast<int>(std::round(angle / fRounding))]++;
76 if(fVerbosity == kAll) {
77 std::cout << "after folding and rounding: angle " << angle << ", counts " << fAngleCount[static_cast<int>(std::round(angle / fRounding))] << std::endl;
78 }
79 } // second crystal loop
80 } //second detector loop
81 } // first crystal loop
82 } // first detector loop
83
84 // create map of indices before we group so that we have an index for each (folded) angle
85 for(auto it = fAngles.begin(); it != fAngles.end(); ++it) {
86 fAngleMap.insert(std::make_pair(*it, std::distance(fAngles.begin(), it)));
87 }
88
89 if(fGrouping) {
90 // If grouping is enable we group angles that are close to each other together.
91 // Which angles are grouped is somewhat arbitrary, this grouping was chosen to get similar numbers of detector combinations and thus statistics for each angle group.
92 // Due to the way lower_bound works, we use the highest angle of each group as the angle of that group.
93 // This is just for the purpose of this algorithm, when plotting the correct average angle of the group should be used!
94
95 // check if we have custom grouping supplied and its format
96 if(!fCustomGrouping.empty()) {
97 if(fCustomGrouping.size() != fAngles.size()) {
98 std::ostringstream str;
99 str << DRED << "Custom grouping with " << fCustomGrouping.size() << " entries supplied, but we have " << fAngles.size() << " angles!" << RESET_COLOR;
100 throw std::runtime_error(str.str());
101 }
102 if(fCustomGrouping[0] != 0) {
103 std::ostringstream str;
104 str << DRED << "Custom grouping with " << fCustomGrouping.size() << " entries supplied, but the first entry is " << fCustomGrouping[0] << ", not 0!" << RESET_COLOR;
105 throw std::runtime_error(str.str());
106 }
107 if(!std::is_sorted(fCustomGrouping.begin(), fCustomGrouping.end())) {
108 std::ostringstream str;
109 str << DRED << "Custom grouping with " << fCustomGrouping.size() << " entries is not sorted!" << RESET_COLOR;
110 throw std::runtime_error(str.str());
111 }
112 auto tmp = fCustomGrouping;
113 auto end = std::unique(tmp.begin(), tmp.end());
114 if(fCustomGrouping.back() != end - tmp.begin() - 1) {
115 std::ostringstream str;
116 str << DRED << "Custom grouping with " << fCustomGrouping.size() << " entries supplied, but the last entry is " << fCustomGrouping[0] << " which does not match the number of unique groups minus one (" << end - tmp.begin() - 1 << ")!" << RESET_COLOR;
117 throw std::runtime_error(str.str());
118 }
119 } else if(fDistance == 110.) {
120 fCustomGrouping = {0, 1, 2, 2, 3, 3, 4, 4, 5, 5, 5, 6, 6, 6, 7, 7, 7, 8, 8, 8, 9, 9, 9, 10, 10, 10};
121 } else if(fDistance == 145.) {
122 fCustomGrouping = {0, 1, 2, 3, 3, 4, 4, 5, 5, 5, 6, 6, 6, 7, 7, 7, 8, 8, 9, 9, 9, 10, 10, 10, 11, 11};
123 } else {
124 std::ostringstream str;
125 str << DRED << "No custom grouping supplied and distance is not 110 mm or 145 mm, but " << fDistance << "!" << RESET_COLOR;
126 throw std::runtime_error(str.str());
127 }
128
129 // now we can loop over the grouping vector and the angles and store the largest angle per group in the new groupedAngles vector
130 // we also update our angle map with the new indices
131 std::set<double> groupedAngles;
132 auto angle = fAngles.begin();
133 double previousAngle = *angle;
134 for(size_t i = 0; i < fCustomGrouping.size(); ++i, ++angle) {
135 if(i != 0 && fCustomGrouping[i] != fCustomGrouping[i - 1]) { // we found a new group so we insert the previous angle
136 groupedAngles.insert(previousAngle);
137 }
138 fAngleMap[*angle] = fCustomGrouping[i];
139 previousAngle = *angle;
140 }
141 // add the last angle since we never get to the next group after the last
142 groupedAngles.insert(previousAngle);
143
144 // update the angles to the new grouped angles
145 fAngles = groupedAngles;
146 }
147}
148
149int TGriffinAngles::Index(double angle)
150{
151 auto matches = [angle, this](std::pair<double, int> val) { return std::abs(val.first - angle) < fRounding; };
152 auto it = std::find_if(fAngleMap.begin(), fAngleMap.end(), matches);
153
154 if(it == fAngleMap.end()) {
155 std::cout << "Failed to find angle " << angle << " in map!" << std::endl;
156 return -1;
157 }
158 if(it->second >= static_cast<int>(fAngles.size())) {
159 std::cout << "Found index " << it->second << " for angle " << angle << " which is outside range (" << fAngles.size() << ")" << std::endl;
160 return -1;
161 }
162 return it->second;
163}
164
165double TGriffinAngles::AverageAngle(int index) const
166{
167 double result = 0.;
168 int nofMatches = 0;
169 for(const auto& val : fAngleMap) {
170 if(val.second == index) {
171 result += val.first;
172 ++nofMatches;
173 }
174 }
175 return result / nofMatches;
176}
177
178void TGriffinAngles::FoldOrGroup(TGraphErrors* z0, TGraphErrors* z2, TGraphErrors* z4, bool verbose) const
179{
180 /// Apply folding and/or grouping to the theory graphs.
181 /// This assumes that the theory graphs all have the exact same length of 49 or 51 for singles or addback, respectively.
182
183 // these are simulated data, so if we add points together we take their average as the new value and change the uncertainties accordingly
184#if ROOT_VERSION_CODE >= ROOT_VERSION(6, 20, 0)
185 if(verbose) {
186 std::cout << "-------------------- original --------------------" << std::endl;
187 for(int i = 0; i < z0->GetN(); ++i) {
188 std::cout << std::setw(8) << z0->GetPointX(i) << ": " << std::setw(8) << z0->GetPointY(i) << " " << std::setw(8) << z2->GetPointY(i) << " " << std::setw(8) << z4->GetPointY(i) << std::endl;
189 }
190 std::cout << "--------------------------------------------------" << std::endl;
191 }
192#endif
193
194 // folding first
195 if(fFolding) {
196 // we first change the angles of the data points, then we re-sort the graphs, and finally we combine points for the same angle into one
197 auto* angle0 = z0->GetX();
198 auto* angle2 = z2->GetX();
199 auto* angle4 = z4->GetX();
200 for(int i = 0; i < z0->GetN(); ++i) {
201 if(angle0[i] > 90.) {
202#if ROOT_VERSION_CODE >= ROOT_VERSION(6, 20, 0)
203 if(verbose) {
204 std::cout << i << ": Folding from " << std::setw(8) << z0->GetPointX(i) << " " << std::setw(8) << z2->GetPointX(i) << " " << std::setw(8) << z4->GetPointX(i);
205 }
206#endif
207 angle0[i] = 180. - angle0[i];
208 angle2[i] = angle0[i];
209 angle4[i] = angle0[i];
210#if ROOT_VERSION_CODE >= ROOT_VERSION(6, 20, 0)
211 if(verbose) {
212 std::cout << " to " << std::setw(8) << z0->GetPointX(i) << " " << std::setw(8) << z2->GetPointX(i) << " " << std::setw(8) << z4->GetPointX(i) << std::endl;
213 }
214#endif
215 }
216 }
217 z0->Sort();
218 z2->Sort();
219 z4->Sort();
220#if ROOT_VERSION_CODE >= ROOT_VERSION(6, 20, 0)
221 if(verbose) {
222 std::cout << "-------------------- sorted --------------------" << std::endl;
223 for(int i = 0; i < z0->GetN(); ++i) {
224 std::cout << std::setw(8) << z0->GetPointX(i) << ": " << std::setw(8) << z0->GetPointY(i) << " " << std::setw(8) << z2->GetPointY(i) << " " << std::setw(8) << z4->GetPointY(i) << std::endl;
225 }
226 std::cout << "------------------------------------------------" << std::endl;
227 }
228#endif
229 angle0 = z0->GetX();
230 auto* counts0 = z0->GetY();
231 auto* errors0 = z0->GetEY();
232 auto* counts2 = z2->GetY();
233 auto* errors2 = z2->GetEY();
234 auto* counts4 = z4->GetY();
235 auto* errors4 = z4->GetEY();
236 for(int i = 1; i < z0->GetN(); ++i) {
237 if(std::abs(angle0[i] - angle0[i - 1]) < fRounding) {
238 if(verbose) {
239 std::cout << i << ": Adding angles " << angle0[i - 1] << " and " << angle0[i] << std::endl;
240 }
241 // add counts from this point to the previous point, then delete this point, and update our pointers to the data
242 counts0[i - 1] += counts0[i];
243 counts0[i - 1] /= 2.;
244 errors0[i - 1] = TMath::Sqrt(TMath::Power(errors0[i - 1], 2) + TMath::Power(errors0[i], 2)) / 2.;
245 counts2[i - 1] += counts2[i];
246 counts2[i - 1] /= 2.;
247 errors2[i - 1] = TMath::Sqrt(TMath::Power(errors2[i - 1], 2) + TMath::Power(errors2[i], 2)) / 2.;
248 counts4[i - 1] += counts4[i];
249 counts4[i - 1] /= 2.;
250 errors4[i - 1] = TMath::Sqrt(TMath::Power(errors4[i - 1], 2) + TMath::Power(errors4[i], 2)) / 2.;
251 z0->RemovePoint(i);
252 z2->RemovePoint(i);
253 z4->RemovePoint(i);
254 angle0 = z0->GetX();
255 counts0 = z0->GetY();
256 errors0 = z0->GetEY();
257 counts2 = z2->GetY();
258 errors2 = z2->GetEY();
259 counts4 = z4->GetY();
260 errors4 = z4->GetEY();
261 --i; // decrement to re-check the new ith point
262 } else {
263 if(verbose) {
264 std::cout << i << ": Not adding angles " << angle0[i - 1] << " and " << angle0[i] << std::endl;
265 }
266 }
267 }
268 }
269#if ROOT_VERSION_CODE >= ROOT_VERSION(6, 20, 0)
270 if(verbose) {
271 std::cout << "-------------------- after folding --------------------" << std::endl;
272 for(int i = 0; i < z0->GetN(); ++i) {
273 std::cout << std::setw(8) << z0->GetPointX(i) << ": " << std::setw(8) << z0->GetPointY(i) << " " << std::setw(8) << z2->GetPointY(i) << " " << std::setw(8) << z4->GetPointY(i) << std::endl;
274 }
275 std::cout << "-------------------------------------------------------" << std::endl;
276 }
277#endif
278
279 // grouping
280 if(fGrouping) {
281 // Due to the way lower_bound works, we use the highest angle of each group as the angle of that group.
282 // This is just for the purpose of this algorithm, when plotting the correct average angle of the group should be used!
283 auto* counts0 = z0->GetY();
284 auto* errors0 = z0->GetEY();
285 auto* counts2 = z2->GetY();
286 auto* errors2 = z2->GetEY();
287 auto* counts4 = z4->GetY();
288 auto* errors4 = z4->GetEY();
289 for(int i = 0; i < z0->GetN(); ++i) {
290 switch(i) {
291 case 0:
292 case 1: // first and second angle are not grouped
293#if ROOT_VERSION_CODE >= ROOT_VERSION(6, 20, 0)
294 if(verbose) {
295 std::cout << i << ": Leaving as is " << z0->GetPointX(i) << ": " << std::setw(8) << z0->GetPointY(i) << " " << std::setw(8) << z2->GetPointY(i) << " " << std::setw(8) << z4->GetPointY(i) << std::endl;
296 }
297#endif
298 break;
299 case 2:
300 case 3:
301 case 4:
302 // three groups of two angles each, so we add the this point to the next one, delete it, and update the pointers
303 // because we delete the point, we don't skip cases here like when we create the angles above
304#if ROOT_VERSION_CODE >= ROOT_VERSION(6, 20, 0)
305 if(verbose) {
306 std::cout << i << ": Grouping " << z0->GetPointX(i) << " and " << z0->GetPointX(i + 1) << " from " << std::setw(8) << z0->GetPointY(i) << " " << std::setw(8) << z2->GetPointY(i) << " " << std::setw(8) << z4->GetPointY(i);
307 }
308#endif
309 counts0[i + 1] += counts0[i];
310 counts0[i + 1] /= 2.;
311 errors0[i + 1] = TMath::Sqrt(TMath::Power(errors0[i], 2) + TMath::Power(errors0[i + 1], 2)) / 2.;
312 counts2[i + 1] += counts2[i];
313 counts2[i + 1] /= 2.;
314 errors2[i + 1] = TMath::Sqrt(TMath::Power(errors2[i], 2) + TMath::Power(errors2[i + 1], 2)) / 2.;
315 counts4[i + 1] += counts4[i];
316 counts4[i + 1] /= 2.;
317 errors4[i + 1] = TMath::Sqrt(TMath::Power(errors4[i], 2) + TMath::Power(errors4[i + 1], 2)) / 2.;
318 z0->RemovePoint(i);
319 z2->RemovePoint(i);
320 z4->RemovePoint(i);
321 counts0 = z0->GetY();
322 errors0 = z0->GetEY();
323 counts2 = z2->GetY();
324 errors2 = z2->GetEY();
325 counts4 = z4->GetY();
326 errors4 = z4->GetEY();
327#if ROOT_VERSION_CODE >= ROOT_VERSION(6, 20, 0)
328 if(verbose) {
329 std::cout << " to " << z0->GetPointX(i) << " with " << std::setw(8) << z0->GetPointY(i) << " " << std::setw(8) << z2->GetPointY(i) << " " << std::setw(8) << z4->GetPointY(i) << std::endl;
330 }
331#endif
332 // no need to decrement the index, by removing one point the old i+1 became i and we don't want to process that one (we just added to it)
333 break;
334 default:
335 // all others are groups of three, so we add this point and the next one to the one two ahead, delete them, and update the pointers
336#if ROOT_VERSION_CODE >= ROOT_VERSION(6, 20, 0)
337 if(verbose) {
338 std::cout << i << ": Grouping " << z0->GetPointX(i) << ", " << z0->GetPointX(i + 1) << ", and " << z0->GetPointX(i + 2) << " from " << std::setw(8) << z0->GetPointY(i) << " " << std::setw(8) << z2->GetPointY(i) << " " << std::setw(8) << z4->GetPointY(i);
339 }
340#endif
341 counts0[i + 2] += counts0[i] + counts0[i + 1];
342 counts0[i + 2] /= 3.;
343 errors0[i + 2] = TMath::Sqrt(TMath::Power(errors0[i], 2) + TMath::Power(errors0[i + 1], 2) + TMath::Power(errors0[i + 2], 2)) / 3.;
344 counts2[i + 2] += counts2[i] + counts2[i + 1];
345 counts2[i + 2] /= 3.;
346 errors2[i + 2] = TMath::Sqrt(TMath::Power(errors2[i], 2) + TMath::Power(errors2[i + 1], 2) + TMath::Power(errors0[i + 2], 2)) / 3.;
347 counts4[i + 2] += counts4[i] + counts4[i + 1];
348 counts4[i + 2] /= 3.;
349 errors4[i + 2] = TMath::Sqrt(TMath::Power(errors4[i], 2) + TMath::Power(errors4[i + 1], 2) + TMath::Power(errors0[i + 2], 2)) / 3.;
350 z0->RemovePoint(i);
351 z2->RemovePoint(i);
352 z4->RemovePoint(i);
353 z0->RemovePoint(i + 1);
354 z2->RemovePoint(i + 1);
355 z4->RemovePoint(i + 1);
356 counts0 = z0->GetY();
357 errors0 = z0->GetEY();
358 counts2 = z2->GetY();
359 errors2 = z2->GetEY();
360 counts4 = z4->GetY();
361 errors4 = z4->GetEY();
362#if ROOT_VERSION_CODE >= ROOT_VERSION(6, 20, 0)
363 if(verbose) {
364 std::cout << " to " << z0->GetPointX(i) << " with " << std::setw(8) << z0->GetPointY(i) << " " << std::setw(8) << z2->GetPointY(i) << " " << std::setw(8) << z4->GetPointY(i) << std::endl;
365 }
366#endif
367 // no need to decrement the index, by removing two points the old i+2 became i and we don't want to processs that one since we just added to it
368 break;
369 }
370 }
371 }
372#if ROOT_VERSION_CODE >= ROOT_VERSION(6, 20, 0)
373 if(verbose) {
374 std::cout << "-------------------- after folding and grouping --------------------" << std::endl;
375 for(int i = 0; i < z0->GetN(); ++i) {
376 std::cout << std::setw(8) << z0->GetPointX(i) << ": " << std::setw(8) << z0->GetPointY(i) << " " << std::setw(8) << z2->GetPointY(i) << " " << std::setw(8) << z4->GetPointY(i) << std::endl;
377 }
378 std::cout << "--------------------------------------------------------------------" << std::endl;
379 }
380#endif
381}
382
383bool TGriffinAngles::ExcludeDetector(int detector) const
384{
385 /// Returns true if any of the detectors in fExcludedDetectors matches the given detector.
386 return std::any_of(fExcludedDetectors.begin(), fExcludedDetectors.end(), [&detector](auto exclude) { return detector == exclude; });
387}
388
389bool TGriffinAngles::ExcludeCrystal(int detector, int crystal) const
390{
391 /// Returns true if any of the crystals in fExcludedCrystals matches the given detector and crystal (using 4*(detector-1)+crystal+1).
392 return std::any_of(fExcludedCrystals.begin(), fExcludedCrystals.end(), [&detector, &crystal](auto exclude) { return 4 * (detector - 1) + crystal + 1 == exclude; });
393}
394
395void TGriffinAngles::Print(Option_t*) const
396{
397 std::cout << "List of unique angles " << std::setw(2) << fAngles.size() << " Map from angles to indices " << std::setw(2) << fAngleMap.size() << " # of combinations " << std::setw(2) << fAngleCount.size() << std::endl;
398 //List of unique angles aa Map from angles to indices mm # of combinations cc"<<std::endl;
399 std::cout << "index angle angle index average angle angle counts" << std::endl;
400 //ii: aa.aaaa aa.aaaa ii aa.aaaa aa.aaaa ccc
401 auto it = fAngles.begin();
402 auto it2 = fAngleMap.begin();
403 auto it3 = fAngleCount.begin();
404 while(it != fAngles.end() || it2 != fAngleMap.end() || it3 != fAngleCount.end()) {
405 if(it != fAngles.end()) {
406 std::cout << std::setw(2) << std::distance(fAngles.begin(), it) << ": " << std::setw(7) << *it << " ";
407 ++it;
408 } else {
409 //ii: aa.aaaa aa.aaaa ii aa.aaaa aa.aaaa ccc
410 std::cout << " ";
411 }
412 if(it2 != fAngleMap.end()) {
413 std::cout << std::setw(7) << it2->first << " " << std::setw(2) << it2->second << " " << std::setw(7) << AverageAngle(it2->second) << " ";
414 ++it2;
415 } else {
416 //aa.aaaa ii aa.aaaa aa.aaaa ccc
417 std::cout << " ";
418 }
419 if(it3 != fAngleCount.end()) {
420 std::cout << std::setw(7) << it3->first * fRounding << " " << std::setw(3) << it3->second;
421 ++it3;
422 }
423 std::cout << std::endl;
424 }
425}
426
428{
429 if(fDistance != griffinAngles->fDistance) {
430 std::cerr << "Warning, merging files with different Griffin distances " << fDistance << " != " << griffinAngles->fDistance << std::endl;
431 }
432 if(fFolding != griffinAngles->fFolding) {
433 std::cerr << "Warning, merging files with different folding settings " << (fFolding ? "true" : "false") << " != " << (griffinAngles->fFolding ? "true" : "false") << std::endl;
434 }
435 if(fGrouping != griffinAngles->fGrouping) {
436 std::cerr << "Warning, merging files with different grouping settings " << (fGrouping ? "true" : "false") << " != " << (griffinAngles->fGrouping ? "true" : "false") << std::endl;
437 }
438 if(fAddback != griffinAngles->fAddback) {
439 std::cerr << "Warning, merging files with different addback settings " << (fAddback ? "true" : "false") << " != " << (griffinAngles->fAddback ? "true" : "false") << std::endl;
440 }
441 if(fRounding != griffinAngles->fRounding) {
442 std::cerr << "Warning, merging files with different rounding " << fRounding << " != " << griffinAngles->fRounding << std::endl;
443 }
444 if(fExcludedDetectors != griffinAngles->fExcludedDetectors) {
445 std::cerr << "Warning, merging files with different detectors excluded < ";
446 for(auto det : fExcludedDetectors) { std::cerr << det << " "; }
447 std::cerr << "> != < ";
448 for(auto det : griffinAngles->fExcludedDetectors) { std::cerr << det << " "; }
449 std::cerr << ">" << std::endl;
450 }
451 if(fExcludedCrystals != griffinAngles->fExcludedCrystals) {
452 std::cerr << "Warning, merging files with different crystals excluded < ";
453 for(auto cry : fExcludedCrystals) { std::cerr << cry << " "; }
454 std::cerr << "> != < ";
455 for(auto cry : griffinAngles->fExcludedCrystals) { std::cerr << cry << " "; }
456 std::cerr << ">" << std::endl;
457 }
458 if(fCustomGrouping != griffinAngles->fCustomGrouping) {
459 std::cerr << "Warning, merging files with different custom grouping < ";
460 for(auto group : fCustomGrouping) { std::cerr << group << " "; }
461 std::cerr << "> != < ";
462 for(auto group : griffinAngles->fCustomGrouping) { std::cerr << group << " "; }
463 std::cerr << ">" << std::endl;
464 }
465 if(fAngles != griffinAngles->fAngles) {
466 std::cerr << "Warning, merging files with different angles < ";
467 for(auto angle : fAngles) { std::cerr << angle << " "; }
468 std::cerr << "> != < ";
469 for(auto angle : griffinAngles->fAngles) { std::cerr << angle << " "; }
470 std::cerr << ">" << std::endl;
471 }
472 if(fAngleMap != griffinAngles->fAngleMap) {
473 std::cerr << "Warning, merging files with different angle maps < ";
474 for(auto angle : fAngleMap) { std::cerr << angle.first << "/" << angle.second << " "; }
475 std::cerr << "> != < ";
476 for(auto angle : griffinAngles->fAngleMap) { std::cerr << angle.first << "/" << angle.second << " "; }
477 std::cerr << ">" << std::endl;
478 }
479 if(fAngleCount != griffinAngles->fAngleCount) {
480 std::cerr << "Warning, merging files with different angle counts < ";
481 for(auto angle : fAngleCount) { std::cerr << angle.first << "/" << angle.second << " "; }
482 std::cerr << "> != < ";
483 for(auto angle : griffinAngles->fAngleCount) { std::cerr << angle.first << "/" << angle.second << " "; }
484 std::cerr << ">" << std::endl;
485 }
486}
EVerbosity
Definition Globals.h:143
@ kQuiet
Definition Globals.h:144
@ kAll
Definition Globals.h:148
#define DRED
Definition Globals.h:18
#define RESET_COLOR
Definition Globals.h:5
static TGRSIOptions * Get(int argc=0, char **argv=nullptr)
Do not use!
static TUserSettings * UserSettings()
void Add(TGriffinAngles *griffinAngles)
std::set< double > fAngles
set of unique angles, when grouping is used, the largest angle of the group is used!
static double fRounding
we consider any angles whose difference is less than this to be equal
void FoldOrGroup(TGraphErrors *z0, TGraphErrors *z2, TGraphErrors *z4, bool verbose=false) const
bool ExcludeDetector(int detector) const
bool fGrouping
flag indicating whether we group close angles together
std::set< double >::iterator end() const
std::vector< int > fCustomGrouping
list of custom groups
void Print(Option_t *="") const override
double fDistance
distance of detector from center of array in mmm
bool fFolding
flag indicating whether we fold our distribution around 90 degree
bool fAddback
flag indicating whether we use addback
TGriffinAngles(double distance=145., bool folding=false, bool grouping=false, bool addback=true)
int Index(double angle)
std::map< int, int > fAngleCount
Maps angles (divided by rounding and cast to integers) to number of combinations contributing to it.
std::vector< int > fExcludedDetectors
list of detectors that are excluded in calculating the angles
std::map< double, int > fAngleMap
Maps angles to indices. This is fairly straight forward without grouping, but if grouping is used mul...
static EVerbosity fVerbosity
verbosity level
std::vector< int > fExcludedCrystals
list of crystals that are excluded in calculating the angles, the crystals are numbered as 4*(det-1)+...
bool ExcludeCrystal(int detector, int crystal) const
double AverageAngle(int index) const
static TVector3 GetPosition(int DetNbr, int CryNbr=5, double dist=110.0)
!
Definition TGriffin.cxx:332