SND@LHC Software
Loading...
Searching...
No Matches
snd::analysis_tools::ScifiPlane Class Reference

#include <sndScifiPlane.h>

Collaboration diagram for snd::analysis_tools::ScifiPlane:

Classes

struct  ScifiHit
 
struct  xy_pair
 

Public Member Functions

 ScifiPlane (std::vector< sndScifiHit * > snd_hits, const Configuration &configuration, Scifi *scifi_geometry, int station)
 
const int GetStation () const
 
const std::vector< ScifiHitGetHits () const
 
const ROOT::Math::XYZPoint GetCentroid () const
 
const ROOT::Math::XYZPoint GetCentroidError () const
 
const xy_pair< double > GetTotQdc (bool only_positive=false) const
 
const xy_pair< double > GetTotEnergy (bool only_positive=false) const
 
const xy_pair< int > GetNHits () const
 
const ROOT::Math::XYZPoint GetCluster (int max_gap) const
 
void FindCentroid ()
 
bool HasShower () const
 
void TimeFilter (double min_timestamp, double max_timestamp)
 
xy_pair< double > GetPointQdc (const ROOT::Math::XYZPoint &point, double radius) const
 

Private Attributes

std::vector< ScifiHithits_
 
Configuration configuration_
 
ROOT::Math::XYZPoint centroid_
 
ROOT::Math::XYZPoint centroid_error_
 
int station_
 

Detailed Description

Definition at line 13 of file sndScifiPlane.h.

Constructor & Destructor Documentation

◆ ScifiPlane()

snd::analysis_tools::ScifiPlane::ScifiPlane ( std::vector< sndScifiHit * >  snd_hits,
const Configuration configuration,
Scifi scifi_geometry,
int  station 
)

Definition at line 16 of file sndScifiPlane.cxx.

16 : configuration_(configuration), centroid_(std::nan(""), std::nan(""), std::nan("")), centroid_error_(std::nan(""), std::nan(""), std::nan("")), station_(station)
17{
18 for ( auto snd_hit : snd_hits)
19 {
20 if (!snd_hit->isValid()) continue;
21 ScifiHit hit;
22 hit.channel_index = 512 * snd_hit->GetMat() + 128 * snd_hit->GetSiPM() + snd_hit->GetSiPMChan();
23 int detectorID = snd_hit->GetDetectorID();
25 hit.timestamp = snd_hit->GetTime(0) / ShipUnit::snd_TDC2ns; // timestamp is in clock cycles
26 }
27 else {
28 hit.timestamp = (scifi_geometry->GetCorrectedTime(detectorID, snd_hit->GetTime(0)*ShipUnit::snd_TDC2ns, 0) / ShipUnit::snd_TDC2ns); // timestamp is in clock cycles
29 }
30 hit.qdc = snd_hit->GetSignal(0);
31 hit.is_x = snd_hit->isVertical();
32
33 TVector3 A, B;
34 scifi_geometry->GetSiPMPosition(detectorID, A, B);
35 hit.z = A.Z();
36 if (hit.is_x)
37 {
38 hit.x = A.X();
39 hit.y = std::nan("");
40 }
41 else
42 {
43 hit.x = std::nan("");
44 hit.y = A.Y();
45 }
46 hits_.push_back(hit);
47 }
48}
void GetSiPMPosition(Int_t SiPMChan, TVector3 &A, TVector3 &B)
Definition Scifi.cxx:727
Double_t GetCorrectedTime(Int_t fDetectorID, Double_t rawTime, Double_t L)
Definition Scifi.cxx:614
ROOT::Math::XYZPoint centroid_error_
std::vector< ScifiHit > hits_
ROOT::Math::XYZPoint centroid_

Member Function Documentation

◆ FindCentroid()

void snd::analysis_tools::ScifiPlane::FindCentroid ( )

Definition at line 229 of file sndScifiPlane.cxx.

230{
231 centroid_.SetXYZ(0, 0, 0);
232 double tot_qdc_x{0}, sum_qdc2_x{0.0};
233 double tot_qdc_y{0}, sum_qdc2_y{0.0};
234 std::vector<ScifiHit> cleaned_hits = hits_;
235
236 cleaned_hits.erase(std::remove_if(cleaned_hits.begin(), cleaned_hits.end(),
237 [&](auto &hit)
238 { return hit.qdc <= 0; }),
239 cleaned_hits.end());
240 int counts_x = std::count_if(cleaned_hits.begin(), cleaned_hits.end(), [](auto &hit)
241 { return hit.is_x; });
242 int counts_y = static_cast<int>(cleaned_hits.size())-counts_x;
244 centroid_ = ROOT::Math::XYZPoint(std::nan(""), std::nan(""), std::nan(""));
245 return;
246 }
247
248 for (auto &hit : cleaned_hits)
249 {
250 if (hit.qdc > 0)
251 {
252 if (hit.is_x)
253 {
254 centroid_.SetX(centroid_.X() + hit.x * hit.qdc);
255 tot_qdc_x += hit.qdc;
256 sum_qdc2_x += hit.qdc*hit.qdc;
257 }
258 else
259 {
260 centroid_.SetY(centroid_.Y() + hit.y * hit.qdc);
261 tot_qdc_y += hit.qdc;
262 sum_qdc2_y += hit.qdc*hit.qdc;
263 }
264
265 centroid_.SetZ(centroid_.Z() + hit.z * hit.qdc);
266 }
267 }
268 centroid_ = ROOT::Math::XYZPoint((tot_qdc_x > 0) ? centroid_.X() / tot_qdc_x : std::nan(""), (tot_qdc_y > 0) ? centroid_.Y() / tot_qdc_y : std::nan(""), (tot_qdc_x+tot_qdc_y > 0) ? centroid_.Z() / (tot_qdc_x+tot_qdc_y) : std::nan(""));
269 auto qdc_x_error_scaler = sqrt(sum_qdc2_x)/tot_qdc_x;
270 auto qdc_y_error_scaler = sqrt(sum_qdc2_y)/tot_qdc_y;
271 auto qdc_z_error_scaler = sqrt(sum_qdc2_x+sum_qdc2_y)/(tot_qdc_x+tot_qdc_y);
272 centroid_error_ = ROOT::Math::XYZPoint(configuration_.scifi_centroid_error_x*qdc_x_error_scaler,
273 configuration_.scifi_centroid_error_y*qdc_y_error_scaler,
274 configuration_.scifi_centroid_error_z*qdc_z_error_scaler);
275}

◆ GetCentroid()

const ROOT::Math::XYZPoint snd::analysis_tools::ScifiPlane::GetCentroid ( ) const
inline

Definition at line 41 of file sndScifiPlane.h.

41{ return centroid_; }

◆ GetCentroidError()

const ROOT::Math::XYZPoint snd::analysis_tools::ScifiPlane::GetCentroidError ( ) const
inline

Definition at line 42 of file sndScifiPlane.h.

42{ return centroid_error_; }

◆ GetCluster()

const ROOT::Math::XYZPoint snd::analysis_tools::ScifiPlane::GetCluster ( int  max_gap) const

Definition at line 113 of file sndScifiPlane.cxx.

114{
115 std::vector<double> pos_x(configuration_.scifi_n_channels_per_plane, std::nan(""));
116 std::vector<double> pos_y(configuration_.scifi_n_channels_per_plane, std::nan(""));
117
118 for (auto &hit : hits_)
119 {
120 if (hit.qdc > 0)
121 {
122 if (hit.is_x)
123 {
124 pos_x[hit.channel_index] = hit.x;
125 }
126 else
127 {
128 pos_y[hit.channel_index] = hit.y;
129 }
130 }
131 }
132
133 auto largest_cluster = [&](const std::vector<double> &positions)
134 {
135 int n = static_cast<int>(positions.size());
136
137 int best_start = -1, best_end = -1, best_size = 0;
138 int start = -1, gap_count = 0, size = 0;
139
140 // Find the largest cluster
141 for (int i = 0; i < n; ++i)
142 {
143 if (!std::isnan(positions[i]))
144 {
145 if (start == -1)
146 start = i; // Start a new cluster
147 size++;
148 gap_count = 0; // Reset consecutive gap counter
149 }
150 else
151 {
152 gap_count++;
153 if (gap_count > max_gap)
154 { // Too many consecutive gaps, finalize previous cluster
155 if (size > best_size)
156 {
157 best_start = start;
158 best_end = i - gap_count; // End before the excessive gaps
159 best_size = size;
160 }
161 // Reset for a new potential cluster
162 start = -1;
163 gap_count = 0;
164 size = 0;
165 }
166 else
167 {
168 size++; // Gaps within max_gap still count in cluster size
169 }
170 }
171 }
172
173 // Check last cluster
174 if (size > best_size)
175 {
176 best_start = start;
177 best_end = n - 1;
178 best_size = size;
179 }
180
181 // Compute the average of non-gap values in the best cluster
182 if (best_start == -1 || best_end == -1)
183 return std::nan(""); // No valid cluster found
184
185 double sum = 0.0;
186 int count = 0;
187 for (int i = best_start; i <= best_end; ++i)
188 {
189 if (!std::isnan(positions[i]))
190 {
191 sum += positions[i];
192 count++;
193 }
194 }
195
196 return (count > 0) ? sum / count : std::nan("");
197 };
198
199 double cluster_x = largest_cluster(pos_x);
200 double cluster_y = largest_cluster(pos_y);
201 if (!(std::isnan(cluster_x) || std::isnan(cluster_y))) {
202 return ROOT::Math::XYZPoint(cluster_x, cluster_y, hits_[0].z);
203 }
204 return ROOT::Math::XYZPoint(std::nan(""), std::nan(""), std::nan(""));
205}
int i
Definition ShipAna.py:86

◆ GetHits()

const std::vector< ScifiHit > snd::analysis_tools::ScifiPlane::GetHits ( ) const
inline

Definition at line 40 of file sndScifiPlane.h.

40{ return hits_; }

◆ GetNHits()

const snd::analysis_tools::ScifiPlane::xy_pair< int > snd::analysis_tools::ScifiPlane::GetNHits ( ) const

Definition at line 50 of file sndScifiPlane.cxx.

51{
52 xy_pair<int> counts{0, 0};
53 counts.x = std::count_if(hits_.begin(), hits_.end(), [](auto &hit)
54 { return hit.is_x; });
55 counts.y = static_cast<int>(hits_.size()) - counts.x;
56
57 return counts;
58}

◆ GetPointQdc()

snd::analysis_tools::ScifiPlane::xy_pair< double > snd::analysis_tools::ScifiPlane::GetPointQdc ( const ROOT::Math::XYZPoint &  point,
double  radius 
) const

Definition at line 215 of file sndScifiPlane.cxx.

216{
217 xy_pair<double> qdc{0.0, 0.0};
218 for (const auto &hit : hits_) {
219 if (hit.is_x && std::abs(hit.x - point.X()) <= radius) {
220 qdc.x += hit.qdc;
221 }
222 else if (!hit.is_x && std::abs(hit.y - point.Y()) <= radius) {
223 qdc.y += hit.qdc;
224 }
225 }
226 return qdc;
227}

◆ GetStation()

const int snd::analysis_tools::ScifiPlane::GetStation ( ) const
inline

Definition at line 39 of file sndScifiPlane.h.

39{ return station_; }

◆ GetTotEnergy()

const snd::analysis_tools::ScifiPlane::xy_pair< double > snd::analysis_tools::ScifiPlane::GetTotEnergy ( bool  only_positive = false) const

Definition at line 288 of file sndScifiPlane.cxx.

289{
290 xy_pair<double> energy{0.0, 0.0};
291 auto qdc = GetTotQdc(only_positive);
292
295
296 return energy;
297}
const xy_pair< double > GetTotQdc(bool only_positive=false) const

◆ GetTotQdc()

const snd::analysis_tools::ScifiPlane::xy_pair< double > snd::analysis_tools::ScifiPlane::GetTotQdc ( bool  only_positive = false) const

Definition at line 277 of file sndScifiPlane.cxx.

278{
279 xy_pair<double> qdc_sum{0.0, 0.0};
280 qdc_sum.x = std::accumulate(hits_.begin(), hits_.end(), 0.0, [&](double current_sum, auto &hit)
281 { return (hit.is_x && (hit.qdc > 0 || !only_positive)) ? (current_sum + hit.qdc) : current_sum; });
282 qdc_sum.y = std::accumulate(hits_.begin(), hits_.end(), 0.0, [&](double current_sum, auto &hit)
283 { return (!hit.is_x && (hit.qdc > 0 || !only_positive)) ? (current_sum + hit.qdc) : current_sum; });
284
285 return qdc_sum;
286}

◆ HasShower()

bool snd::analysis_tools::ScifiPlane::HasShower ( ) const

Definition at line 60 of file sndScifiPlane.cxx.

61{
63 {
64 throw std::runtime_error{"min_hits > window_width"};
65 }
66
68 {
69 throw std::runtime_error{"window_width > n_channels_per_plane"};
70 }
71
72 xy_pair<std::vector<int>> is_hit;
75
76 for (auto &hit : hits_)
77 {
78 if (hit.channel_index < 0 ||
79 hit.channel_index >= configuration_.scifi_n_channels_per_plane)
80 {
81 throw std::out_of_range{"hit.channel_index out of range"};
82 }
83 (hit.is_x ? is_hit.x : is_hit.y)[hit.channel_index] = 1;
84 }
85
86 auto density = [&](const std::vector<int> &hit_arr)
87 {
88 int count{0};
89
90 // Initial count for the first window
91 for (int i{0}; i < configuration_.scifi_shower_window_width; ++i)
92 {
93 count += hit_arr[i];
94 }
95
97 return true;
98
99 // Slide the window across the array
101 {
102 count += hit_arr[i] - hit_arr[i - configuration_.scifi_shower_window_width];
104 return true;
105 }
106
107 return false;
108 };
109
110 return density(is_hit.x) && density(is_hit.y);
111}

◆ TimeFilter()

void snd::analysis_tools::ScifiPlane::TimeFilter ( double  min_timestamp,
double  max_timestamp 
)

Definition at line 207 of file sndScifiPlane.cxx.

208{
209 hits_.erase(std::remove_if(hits_.begin(), hits_.end(),
210 [&](auto &hit)
211 { return hit.timestamp < min_timestamp || hit.timestamp > max_timestamp; }),
212 hits_.end());
213}

Member Data Documentation

◆ centroid_

ROOT::Math::XYZPoint snd::analysis_tools::ScifiPlane::centroid_
private

Definition at line 58 of file sndScifiPlane.h.

◆ centroid_error_

ROOT::Math::XYZPoint snd::analysis_tools::ScifiPlane::centroid_error_
private

Definition at line 59 of file sndScifiPlane.h.

◆ configuration_

Configuration snd::analysis_tools::ScifiPlane::configuration_
private

Definition at line 57 of file sndScifiPlane.h.

◆ hits_

std::vector<ScifiHit> snd::analysis_tools::ScifiPlane::hits_
private

Definition at line 56 of file sndScifiPlane.h.

◆ station_

int snd::analysis_tools::ScifiPlane::station_
private

Definition at line 61 of file sndScifiPlane.h.


The documentation for this class was generated from the following files: