4 bool triggered =
false;
5 for (
unsigned sliceIdx = 0; sliceIdx != slices->size(); ++sliceIdx) {
7 double _nxhits_slice = 0;
8 double _nyhits_slice = 0;
26 bool Intrusion =
false;
28 for (
unsigned i = 0;
i != slice.
NCell(); ++
i) {
30 if(chit->
IsMC()==
true) ++_trueHits;
33 if(_trueHits==0)
continue;
34 for (
unsigned i = 0;
i != slice.
NCell(); ++
i) {
36 if (chit->
ADC() < 501)
continue;
37 if (chit->
ADC() > 3100) ++_nSatHits;
38 _meanADC += chit->
ADC();
39 float tns = chit->
TNS();
43 if (PX_max < chit->Plane() ) PX_max = chit->
Plane();
44 if (PX_min > chit->
Plane() ) PX_min = chit->
Plane();
45 if (TX_max < tns ) TX_max = tns;
46 if (TX_min > tns ) TX_min = tns;
58 else if (!_penetrated) {
59 if (chit->
Cell() <
_xdelt && entry_id!=1) _penetrated = 1;
66 if (PY_max < chit->Plane() ) PY_max = chit->
Plane();
67 if (PY_min > chit->
Plane() ) PY_min = chit->
Plane();
68 if (TY_max < chit->TNS() ) TY_max = chit->
TNS();
69 if (TY_min > chit->
TNS() ) TY_min = chit->
TNS();
82 else if (!_penetrated) {
83 if (chit->
Cell() <
_ydelt && entry_id!=3) _penetrated = 1;
98 else if (!_penetrated) {
99 if (chit->
Plane() <
_zdelt && entry_id!=5) _penetrated=1;
104 _meanADC = _meanADC/(_nxhits_slice+_nyhits_slice+0.);
108 if(_trueHits>2&&_penetrated==1&&_nSatHits>4&&_meanADC>1200&&_deltaP>2&&_deltaTNS>0&&_nxhits_slice>2&&_nyhits_slice>2){
121 bool triggered =
false;
122 for (
unsigned sliceIdx = 0; sliceIdx != slices->size(); ++sliceIdx) {
133 double PX_min = 9999;
134 double PY_min = 9999;
137 double CX_min = 9999;
138 double CY_min = 9999;
139 double TX_max = -9e7;
140 double TY_max = -9e7;
143 double _trueHits = 0;
150 bool Intrusion =
false;
152 for (
unsigned i = 0;
i != slice.
NCell(); ++
i) {
154 if(chit->
IsMC()==
true) ++_trueHits;
157 if(_trueHits==0)
continue;
158 for (
unsigned i = 0;
i != slice.
NCell(); ++
i) {
160 if (chit->
ADC() < 501)
continue;
161 if (chit->
ADC() > 3100) ++_nSatHits;
162 sumADC += chit->
ADC();
163 sq_sum_ADC += chit->
ADC()*chit->
ADC();
164 _meanADC += chit->
ADC();
165 float tns = chit->
TNS();
169 if (PX_max < chit->Plane() ) PX_max = chit->
Plane();
170 if (PX_min > chit->
Plane() ) PX_min = chit->
Plane();
171 if (CX_max < chit->Cell() ) CX_max = chit->
Cell();
172 if (CX_min > chit->
Cell() ) CX_min = chit->
Cell();
173 if (TX_max < tns ) TX_max = tns;
174 if (TX_min > tns ) TX_min = tns;
186 else if (!_penetrated) {
187 if (chit->
Cell() <
_xdelt && entry_id!=1) _penetrated = 1;
194 if (PY_max < chit->Plane() ) PY_max = chit->
Plane();
195 if (PY_min > chit->
Plane() ) PY_min = chit->
Plane();
196 if (CY_max < chit->Cell() ) CY_max = chit->
Cell();
197 if (CY_min > chit->
Cell() ) CY_min = chit->
Cell();
198 if (TY_max < chit->TNS() ) TY_max = chit->
TNS();
199 if (TY_min > chit->
TNS() ) TY_min = chit->
TNS();
212 else if (!_penetrated) {
213 if (chit->
Cell() <
_ydelt && entry_id!=3) _penetrated = 1;
228 else if (!_penetrated) {
229 if (chit->
Plane() <
_zdelt && entry_id!=5) _penetrated=1;
234 _meanADC = _meanADC/(nx+ny+0.);
238 if(_trueHits>2&&_penetrated==1&&_deltaP>2&&_deltaTNS>0&&nx>2&&ny>2){
239 float tracklength =
Distance(CX_max, CY_max, PX_max, CX_min, CY_min, PX_min);
240 if( tracklength<80)
continue;
242 if( slice.
NCell()*1.0/tracklength < 0.44)
continue;
244 double tracklength_of_xview=
Distance(CX_max, 0, PX_max, CX_min, 0, PX_min);
245 double tracklength_of_yview=
Distance(0, CY_max, PY_max, 0, CY_min, PY_min);
246 double number_of_cells_per_length_xview=0;
247 double number_of_cells_per_length_yview=0;
248 double weighted_off_center_xx=0, weighted_off_center_xz=0, weighted_off_center_yy=0, weighted_off_center_yz=0;
249 double stdev_cells_per_plane_xview=0, stdev_cells_per_plane_yview=0;
251 NumberOfCellsPerLength(slice, tracklength_of_xview,tracklength_of_yview,number_of_cells_per_length_xview,number_of_cells_per_length_yview, PX_max, PY_max, PX_min, PY_min);
252 if( number_of_cells_per_length_xview<0.2 || number_of_cells_per_length_yview<0.2)
continue;
255 StdevCellsPerPlane(slice, stdev_cells_per_plane_xview, stdev_cells_per_plane_yview, PX_max, PX_min, PY_max, PY_min );
256 if (stdev_cells_per_plane_xview > 1.25 || stdev_cells_per_plane_yview > 1.2)
continue;
257 WeightedCenterCut(slice, PX_min, PX_max, PY_min, PY_max, CX_min, CX_max, CY_min, CY_max, weighted_off_center_xx, weighted_off_center_xz, weighted_off_center_yy, weighted_off_center_yz );
258 if (weighted_off_center_xz > 0.13 || weighted_off_center_yz > 0.13)
continue;
276 std::vector<std::vector<int>> x_hits_plane(PX_max-PX_min+1,std::vector<int>(0));
277 std::vector<std::vector<int>> y_hits_plane(PY_max-PY_min+1,std::vector<int>(0));
278 for (
auto const hit_ptr: hits.
AllCells()) {
282 x_hits_plane[hit.
Plane()-PX_min].push_back(hit.
Cell());
286 y_hits_plane[hit.
Plane()-PY_min].push_back(hit.
Cell());
290 std::vector<int> x_cells_plane(PX_max-PX_min+1,0);
291 std::vector<int> y_cells_plane(PY_max-PY_min+1,0);
292 for (
unsigned i=0;
i!=x_hits_plane.size();
i++){
293 std::sort(x_hits_plane[
i].begin(),x_hits_plane[
i].
end());
294 std::vector<int>::iterator
it;
295 it = std::unique (x_hits_plane[
i].begin(),x_hits_plane[
i].
end());
297 x_cells_plane[
i]=x_hits_plane[
i].size();
299 for (
unsigned i=0;
i!=y_hits_plane.size();
i++){
300 std::sort(y_hits_plane[
i].begin(),y_hits_plane[
i].
end());
301 std::vector<int>::iterator
it;
302 it = std::unique (y_hits_plane[
i].begin(),y_hits_plane[
i].
end());
304 y_cells_plane[
i]=y_hits_plane[
i].size();
308 for (
unsigned i=0;
i!=x_cells_plane.size();
i++){
309 sum_of_cells+=x_cells_plane[
i];
311 number_of_cells_per_length_xview=sum_of_cells/tracklength_of_xview;
313 for (
unsigned i=0;
i!=y_cells_plane.size();
i++){
314 sum_of_cells+=y_cells_plane[
i];
316 number_of_cells_per_length_yview=sum_of_cells/tracklength_of_yview;
321 int number_of_surface_hits=0;
322 for(
unsigned i=0;
i!=hits.
NCell();++
i){
325 return number_of_surface_hits;
341 std::vector<std::vector<int>> x_hits_plane(PX_max-PX_min+1,std::vector<int>(0));
342 std::vector<std::vector<int>> y_hits_plane(PY_max-PY_min+1,std::vector<int>(0));
343 for (
auto const hit_ptr: hits.
AllCells()) {
347 x_hits_plane[hit.
Plane()-PX_min].push_back(hit.
Cell());
351 y_hits_plane[hit.
Plane()-PY_min].push_back(hit.
Cell());
359 std::vector<int> x_cells_plane(PX_max-PX_min+1,0);
360 std::vector<int> y_cells_plane(PY_max-PY_min+1,0);
361 for (
unsigned i=0;
i!=x_hits_plane.size();
i++){
362 std::sort(x_hits_plane[
i].begin(),x_hits_plane[
i].
end());
363 std::vector<int>::iterator
it;
364 it = std::unique (x_hits_plane[
i].begin(),x_hits_plane[
i].
end());
366 x_cells_plane[
i]=x_hits_plane[
i].size();
368 for (
unsigned i=0;
i!=y_hits_plane.size();
i++){
369 std::sort(y_hits_plane[
i].begin(),y_hits_plane[
i].
end());
370 std::vector<int>::iterator
it;
371 it = std::unique (y_hits_plane[
i].begin(),y_hits_plane[
i].
end());
373 y_cells_plane[
i]=y_hits_plane[
i].size();
375 sum = std::accumulate(x_cells_plane.begin(),x_cells_plane.end(),0.0);
376 mean = sum / x_cells_plane.size();
377 sq_sum = std::inner_product(x_cells_plane.begin(),x_cells_plane.end(),x_cells_plane.begin(),0.0);
378 stdev =
std::sqrt(sq_sum / x_cells_plane.size()-mean *
mean);
381 stdev_cells_per_plane_xview=stdev/
mean;
382 sum = std::accumulate(y_cells_plane.begin(),y_cells_plane.end(),0.0);
383 mean = sum / y_cells_plane.size();
384 sq_sum = std::inner_product(y_cells_plane.begin(),y_cells_plane.end(),y_cells_plane.begin(),0.0);
385 stdev =
std::sqrt(sq_sum / y_cells_plane.size()-mean *
mean);
386 stdev_cells_per_plane_yview=stdev/
mean;
393 return std::sqrt((plane1-plane2)*(plane1-plane2)*2.82072+(xcell1-xcell2)*(xcell1-xcell2)+(ycell1-ycell2)*(ycell1-ycell2));
397 double zcl::FastMonopoleTriggers::WeightedCenterCut(
rb::Cluster const &
hits,
float PX_min,
float PX_max,
float PY_min,
float PY_max,
float CX_min,
float CX_max,
float CY_min,
float CY_max,
double &weighted_off_center_xx,
double &weighted_off_center_xz,
double &weighted_off_center_yy,
double &weighted_off_center_yz )
const{
398 float Weighted_PX_Center=0.5,Weighted_CX_Center=0.5,Weighted_PY_Center=0.5,Weighted_CY_Center=0.5,SumOfXHitsADC=0,SumOfYHitsADC=0;
399 for (
auto const hit_ptr: hits.
AllCells()) {
402 SumOfXHitsADC+=hit.
ADC();
403 Weighted_PX_Center+=hit.
ADC()*(hit.
Plane()-PX_min);
404 Weighted_CX_Center+=hit.
ADC()*(hit.
Cell()-CX_min);
407 SumOfYHitsADC+=hit.
ADC();
408 Weighted_PY_Center+=hit.
ADC()*(hit.
Plane()-PY_min);
409 Weighted_CY_Center+=hit.
ADC()*(hit.
Cell()-CY_min);
412 if(PX_max!=PX_min) Weighted_PX_Center=Weighted_PX_Center/SumOfXHitsADC/(PX_max-PX_min);
413 if(PY_max!=PY_min) Weighted_PY_Center=Weighted_PY_Center/SumOfYHitsADC/(PY_max-PY_min);
414 if(CX_max!=CX_min) Weighted_CX_Center=Weighted_CX_Center/SumOfXHitsADC/(CX_max-CX_min);
415 if(CY_max!=CY_min) Weighted_CY_Center=Weighted_CY_Center/SumOfYHitsADC/(CY_max-CY_min);
419 weighted_off_center_xx=
std::abs(Weighted_CX_Center-0.5);
420 weighted_off_center_xz=
std::abs(Weighted_PX_Center-0.5);
421 weighted_off_center_yy=
std::abs(Weighted_CY_Center-0.5);
422 weighted_off_center_yz=
std::abs(Weighted_PY_Center-0.5);
429 double number_hits_in_overlap_planes=0;
430 for (
auto const hit_ptr: hits.
AllCells()) {
432 if(hit.
Plane()<=max && hit.
Plane()>=
min) number_hits_in_overlap_planes++;
435 return number_hits_in_overlap_planes/hits.
NCell();
int NumberOfSurfaceHits(rb::Cluster const &hits) const
T max(const caf::Proxy< T > &a, T b)
void StdevCellsPerPlane(rb::Cluster const &hits, double &stdev_cells_per_plane_xview, double &stdev_cells_per_plane_yview, double PX_max, double PX_min, double PY_max, double PY_min) const
unsigned int NCell(geo::View_t view) const
Number of cells in view view.
unsigned short Plane() const
Vertical planes which measure X.
A collection of associated CellHits.
double NumberOfHitsInOverlapPlanesCut(rb::Cluster const &hits, float PX_min, float PX_max, float PY_min, float PY_max) const
float Distance(double const &xcell1, double const &ycell1, double const &plane1, double const &xcell2, double const &ycell2, double const &plane2) const
unsigned distance(const T &t1, const T &t2)
bool SurfAssign(rb::CellHit const &hit) const
Horizontal planes which measure Y.
art::PtrVector< rb::CellHit > AllCells() const
Get all cells from both views.
unsigned short Cell() const
double WeightedCenterCut(rb::Cluster const &hits, float PX_min, float PX_max, float PY_min, float PY_max, float CX_min, float CX_max, float CY_min, float CY_max, double &weighted_off_center_xx, double &weighted_off_center_xz, double &weighted_off_center_yy, double &weighted_off_center_yz) const
static float min(const float a, const float b, const float c)
art::Ptr< rb::CellHit > Cell(geo::View_t view, unsigned int viewIdx) const
Get the ith cell from view view.
A rawdata::RawDigit with channel information decoded.
int16_t ADC(uint32_t i) const
T min(const caf::Proxy< T > &a, T b)
bool is_trigger_by_epoch1_fmmtrigger(art::Handle< std::vector< rb::Cluster > > slices, int &_penetrated, double &_nSatHits, double &_meanADC)
T max(sqlite3 *const db, std::string const &table_name, std::string const &column_name)
bool passed_epoch2_fmmtrigger(art::Handle< std::vector< rb::Cluster > > slices, int &_penetrated, double &_nSatHits, double &_meanADC)
void NumberOfCellsPerLength(rb::Cluster const &hits, double tracklength_of_xview, double tracklength_of_yview, double &number_of_cells_per_length_xview, double &number_of_cells_per_length_yview, double const &PX_max, double const &PY_max, double const &PX_min, double const &PY_min) const