346 Double_t px, Double_t py, Double_t pz, Double_t signal, Int_t refID) {
353 if(signal <= 0.0) signal = 1e-16;
355 Double_t AvalancheRadius = 0.01;
357 Double_t radius = AvalancheRadius;
358 if(radius <= 0.0)
return false;
364 radius = gRandom->Gaus(AvalancheRadius, 0.005);
365 if(radius > AvalancheRadius/2.0 && radius < AvalancheRadius*2.0 && radius > 0.0)
break;
369 radius = AvalancheRadius;
375 Int_t n_strip_layers = StripLayers.size();
376 vector<StripCluster> cluster_layers(n_strip_layers,
StripCluster());
378 for(Int_t ilayer = 0; ilayer < n_strip_layers; ++ilayer) {
379 if( !StripLayers[ilayer].IsPointInsideStripLayer(xloc, yloc) )
continue;
380 cluster_layers[ilayer] =
MakeCluster(ilayer, xloc, yloc, signal, radius);
385 for (Int_t ilayer = 0; ilayer < n_strip_layers; ++ilayer) {
386 for (
size_t ielement = 0; ielement < cluster_layers[ilayer].Strips.size(); ++ielement) {
387 Int_t strip_num = cluster_layers[ilayer].Strips.at(ielement);
388 Double_t strip_signal = cluster_layers[ilayer].Signals.at(ielement);
389 Int_t zone_id = StripLayers[ilayer].GetZoneID();
396 for (Int_t ilayer = 0; ilayer < n_strip_layers; ++ilayer) {
397 for (
size_t ielement = 0; ielement < cluster_layers[ilayer].Strips.size(); ++ielement) {
398 Int_t strip_num = cluster_layers[ilayer].Strips.at(ielement);
399 Double_t strip_signal = cluster_layers[ilayer].Signals.at(ielement);
400 Int_t zone_id = StripLayers[ilayer].GetZoneID();
406 RealPointsX.push_back(x);
407 RealPointsY.push_back(y);
408 RealPointsMC.push_back(refID);
413 if(Verbosity) cerr <<
"WARNING: BmnSiliconModule::AddRealPointFullOne(): Point (" << x <<
" : " << y <<
" : " << z <<
") is out of the readout plane or inside a dead zone\n";
420 Double_t px, Double_t py, Double_t pz, Double_t signal, Int_t refID) {
428 Double_t pz_loc = pz;
433 Double_t z_dist_mid = ModuleThickness*0.5;
434 Double_t vector_coeff_mid =
fabs(z_dist_mid/pz_loc);
436 Double_t x_vec_from_in_to_mid = vector_coeff_mid*px_loc;
437 Double_t y_vec_from_in_to_mid = vector_coeff_mid*py_loc;
441 Double_t x_mid = x_loc + x_vec_from_in_to_mid;
442 Double_t y_mid = y_loc + y_vec_from_in_to_mid;
448 if(signal <= 0.0) signal = 1e-16;
450 Double_t AvalancheRadius = 0.015;
452 Double_t radius = AvalancheRadius;
453 if(radius <= 0.0)
return false;
459 radius = gRandom->Gaus(AvalancheRadius, 0.005);
460 if(radius > AvalancheRadius/2.0 && radius < AvalancheRadius*2.0 && radius > 0.0)
break;
464 radius = AvalancheRadius;
470 Int_t n_strip_layers = StripLayers.size();
471 vector<StripCluster> cluster_layers(n_strip_layers,
StripCluster());
473 for(Int_t ilayer = 0; ilayer < n_strip_layers; ++ilayer) {
474 if( !StripLayers[ilayer].IsPointInsideStripLayer(x_mid, y_mid) )
continue;
475 cluster_layers[ilayer] =
MakeCluster(ilayer, x_mid, y_mid, signal, radius);
480 for (Int_t ilayer = 0; ilayer < n_strip_layers; ++ilayer) {
481 for (
size_t ielement = 0; ielement < cluster_layers[ilayer].Strips.size(); ++ielement) {
482 Int_t strip_num = cluster_layers[ilayer].Strips.at(ielement);
483 Double_t strip_signal = cluster_layers[ilayer].Signals.at(ielement);
484 Int_t zone_id = StripLayers[ilayer].GetZoneID();
491 for (Int_t ilayer = 0; ilayer < n_strip_layers; ++ilayer) {
492 for (
size_t ielement = 0; ielement < cluster_layers[ilayer].Strips.size(); ++ielement) {
493 Int_t strip_num = cluster_layers[ilayer].Strips.at(ielement);
494 Double_t strip_signal = cluster_layers[ilayer].Signals.at(ielement);
495 Int_t zone_id = StripLayers[ilayer].GetZoneID();
501 RealPointsX.push_back(x);
502 RealPointsY.push_back(y);
503 RealPointsMC.push_back(refID);
508 if(Verbosity) cerr <<
"WARNING: BmnSiliconModule::AddRealPointFullOne(): Point (" << x <<
" : " << y <<
" : " << z <<
") is out of the readout plane or inside a dead zone\n";
514 Double_t px, Double_t py, Double_t pz, Double_t signal, Int_t refID) {
523 Double_t pz_loc = pz;
528 Double_t chargedCarrierCloudRadius = gRandom->Uniform(0.0001, 0.0005);
529 Double_t step = 0.0010;
535 Double_t z_distance_from_in_to_out;
539 z_distance_from_in_to_out = ModuleThickness - (z_loc-ZStartModulePosition);
542 z_distance_from_in_to_out = z_loc - ZStartModulePosition;
546 if(z_distance_from_in_to_out <= 0.0) {
547 if(Verbosity) cerr <<
"WARNING: BmnSiliconModule::AddRealPointFull(): Point (" << x <<
" : " << y <<
" : " << z <<
") has a track with an incorrect length\n";
552 Double_t vector_coeff =
fabs(z_distance_from_in_to_out/pz_loc);
555 Double_t x_vec_from_in_to_out = vector_coeff*px_loc;
556 Double_t y_vec_from_in_to_out = vector_coeff*py_loc;
557 Double_t z_vec_from_in_to_out = vector_coeff*pz_loc;
560 Double_t x_out = x_loc + x_vec_from_in_to_out;
561 Double_t y_out = y_loc + y_vec_from_in_to_out;
566 if(x_out < XMinModule || x_out > XMaxModule ) {
567 if(x_out < XMinModule ) x_out = XMinModule;
568 if(x_out > XMaxModule ) x_out = XMaxModule;
570 x_vec_from_in_to_out = x_out - x_loc;
571 vector_coeff = x_vec_from_in_to_out/px_loc;
573 y_vec_from_in_to_out= vector_coeff*py_loc;
574 z_vec_from_in_to_out = vector_coeff*pz_loc;
576 y_out = y_loc + y_vec_from_in_to_out;
580 if(y_out < YMinModule || y_out > YMaxModule ) {
581 if(y_out < YMinModule ) y_out = YMinModule;
582 if(y_out > YMaxModule ) y_out = YMaxModule;
584 y_vec_from_in_to_out = y_out - y_loc;
585 vector_coeff = y_vec_from_in_to_out/py_loc;
587 x_vec_from_in_to_out= vector_coeff*px_loc;
588 z_vec_from_in_to_out = vector_coeff*pz_loc;
590 x_out = x_loc + x_vec_from_in_to_out;
596 Double_t track_length = std::sqrt( (x_vec_from_in_to_out)*(x_vec_from_in_to_out) + (y_vec_from_in_to_out)*(y_vec_from_in_to_out) + (z_vec_from_in_to_out)*(z_vec_from_in_to_out) );
598 Double_t current_length = 0.0;
599 Double_t current_length_ratio = 0.0;
601 Int_t collisions_cnt = 0;
602 Double_t current_step = 0.0;
605 std::vector<CollPoint> collision_points;
607 while(current_length < track_length) {
611 current_length += current_step;
613 if(current_length > track_length)
break;
615 current_length_ratio = current_length/track_length;
617 Double_t current_x = x_loc + current_length_ratio*x_vec_from_in_to_out;
618 Double_t current_y = y_loc + current_length_ratio*y_vec_from_in_to_out;
619 Double_t current_z = z_loc + current_length_ratio*z_vec_from_in_to_out;
621 collision_points.push_back(CollPoint(current_x, current_y, current_z));
626 if(signal <= 0.0) signal = 0;
628 Int_t n_collisions = collision_points.size();
629 Double_t signal_per_step = signal / n_collisions;
631 for(Int_t icollpoint = 0; icollpoint < n_collisions; ++icollpoint) {
634 Double_t current_radius = chargedCarrierCloudRadius;
635 if(chargedCarrierCloudRadius <= 0.0) {
636 if(Verbosity) cerr <<
"WARNING: BmnSiliconModule::AddRealPointFull(): incorrect chargedCarrierCloudRadius!\n";
640 Double_t current_x = collision_points[icollpoint].x;
641 Double_t current_y = collision_points[icollpoint].y;
645 for(Int_t ilayer = 0; ilayer < n_strip_layers; ++ilayer) {
646 if( !StripLayers[ilayer].IsPointInsideStripLayer(current_x, current_y) )
continue;
648 Double_t pitch = StripLayers[ilayer].GetPitch();
649 Double_t radius_in_strips = current_radius/pitch;
650 Int_t zone_id = StripLayers[ilayer].GetZoneID();
652 Int_t first_strip_in_layer = StripLayers[ilayer].GetFirstStripNumber();
653 Int_t last_strip_in_layer = StripLayers[ilayer].GetLastStripNumber();
656 Double_t central_strip_pos = StripLayers[ilayer].ConvertPointToStripPosition(current_x, current_y);
657 Double_t left_strip_pos = central_strip_pos - radius_in_strips;
658 Double_t right_strip_pos = central_strip_pos + radius_in_strips;
665 Double_t central_strip_signal = signal_per_step;
666 Double_t left_strip_signal = 0.0;
667 Double_t right_strip_signal = 0.0;
669 Double_t central_strip_fraction = 1.0;
670 Double_t left_strip_fraction = 0.0;
671 Double_t right_strip_fraction = 0.0;
674 if( (Int_t)left_strip_pos != (Int_t)central_strip_pos ) {
675 if( left_strip_pos >= first_strip_in_layer ) {
676 left_strip_fraction = (Int_t)central_strip_pos - left_strip_pos;
677 if(left_strip_fraction > 1.0) left_strip_fraction = 1.0;
681 if( (Int_t)right_strip_pos != (Int_t)central_strip_pos ) {
682 if( right_strip_pos < (last_strip_in_layer+1) ) {
683 right_strip_fraction = right_strip_pos - (Int_t)(central_strip_pos+1);
684 if(right_strip_fraction > 1.0) right_strip_fraction = 1.0;
689 central_strip_signal = signal_per_step*central_strip_fraction / (central_strip_fraction + left_strip_fraction + right_strip_fraction);
690 left_strip_signal = signal_per_step*left_strip_fraction / (central_strip_fraction + left_strip_fraction + right_strip_fraction);
691 right_strip_signal = signal_per_step*right_strip_fraction / (central_strip_fraction + left_strip_fraction + right_strip_fraction);
704 RealPointsX.push_back(x);
705 RealPointsY.push_back(y);
706 RealPointsMC.push_back(refID);
711 if(Verbosity) cerr <<
"WARNING: BmnSiliconModule::AddRealPointFull(): Point (" << x <<
" : " << y <<
" : " << z <<
") is out of the readout plane or inside a dead zone\n";
718 Double_t ClusterDistortion = 0.0;
724 if(radius <= 0.0)
return cluster;
726 Double_t Pitch = StripLayers[layer_num].GetPitch();
728 Double_t RadiusInZones = radius/Pitch;
729 Double_t
Sigma = RadiusInZones/3.33;
731 TF1 gausF(
"gausF",
"gaus", -4*
Sigma, 4*
Sigma);
732 gausF.SetParameter(0, 1.0);
733 gausF.SetParameter(1, 0.0);
734 gausF.SetParameter(2,
Sigma);
736 Double_t SCluster = gausF.Integral(-4*
Sigma, 4*
Sigma);
739 Double_t var_level = ClusterDistortion;
741 Int_t NStripsInLayer = StripLayers[layer_num].GetNStrips();
742 Int_t FirstStripInLayer = StripLayers[layer_num].GetFirstStripNumber();
743 Int_t LastStripInLayer = StripLayers[layer_num].GetLastStripNumber();
745 if(NStripsInLayer == 0)
return cluster;
747 Double_t CenterZonePos = StripLayers[layer_num].ConvertPointToStripPosition(xcoord, ycoord);
749 if( CenterZonePos < FirstStripInLayer || CenterZonePos >= LastStripInLayer+1 )
return cluster;
752 gausF.SetParameter(1, CenterZonePos);
753 Double_t total_signal = 0.0;
755 Double_t LeftZonePos = CenterZonePos - RadiusInZones;
756 Double_t RightZonePos = CenterZonePos + RadiusInZones;
757 Double_t OutLeftBorder = 0.0;
758 Double_t OutRightBorder = 0.0;
760 if(LeftZonePos < FirstStripInLayer) {
761 OutLeftBorder = LeftZonePos;
762 LeftZonePos = FirstStripInLayer;
764 if(RightZonePos < FirstStripInLayer) {
765 OutRightBorder = RightZonePos;
766 RightZonePos = FirstStripInLayer;
768 if(LeftZonePos >= LastStripInLayer+1) {
769 OutLeftBorder = LeftZonePos;
770 LeftZonePos = LastStripInLayer+1 - 0.001;
772 if(RightZonePos >= LastStripInLayer+1) {
773 OutRightBorder = RightZonePos;
774 RightZonePos = LastStripInLayer+1 - 0.001;
781 if((Int_t)LeftZonePos == (Int_t)RightZonePos) {
783 Int_t NumCurrentZone = (Int_t) LeftZonePos;
785 h = RightZonePos - LeftZonePos;
788 Double_t xp = LeftZonePos +
dist;
789 Double_t S = gausF.Integral(xp, xp+h);
790 Double_t Energy = (signal*Gain)*(S/SCluster);
791 Energy += rand.Gaus(0.0, var_level*Energy);
792 if(Energy < 0.0) Energy = 0.0;
794 if(NumCurrentZone >= FirstStripInLayer && NumCurrentZone < LastStripInLayer+1 && Energy > 0.0) {
795 cluster.
AddStrip(NumCurrentZone, Energy);
796 total_signal += Energy;
804 Int_t NumCurrentZone = (Int_t) LeftZonePos;
806 h = ((Int_t)LeftZonePos + 1) - LeftZonePos;
809 Double_t xp = LeftZonePos +
dist;
810 Double_t S = gausF.Integral(xp, xp+h);
811 Double_t Energy = (signal*Gain)*(S/SCluster);
812 Energy += rand.Gaus(0.0, var_level*Energy);
813 if(Energy < 0.0) Energy = 0.0;
815 if(NumCurrentZone >= FirstStripInLayer && NumCurrentZone < LastStripInLayer+1 && Energy > 0.0) {
816 cluster.
AddStrip(NumCurrentZone, Energy);
817 total_signal += Energy;
823 for(Int_t
i = (Int_t)LeftZonePos + 1;
i < (Int_t)RightZonePos; ++
i) {
829 xp = LeftZonePos +
dist;
830 S = gausF.Integral(xp, xp+h);
831 Energy = (signal*Gain)*(S/SCluster);
832 Energy += rand.Gaus(0.0, var_level*Energy);
833 if(Energy < 0.0) Energy = 0.0;
835 if(NumCurrentZone >= FirstStripInLayer && NumCurrentZone < LastStripInLayer+1 && Energy > 0.0) {
836 cluster.
AddStrip(NumCurrentZone, Energy);
837 total_signal += Energy;
843 NumCurrentZone = (Int_t)RightZonePos;
845 h = RightZonePos - (Int_t)RightZonePos;
848 xp = LeftZonePos +
dist;
849 S = gausF.Integral(xp, xp+h);
850 Energy = (signal*Gain)*(S/SCluster);
851 Energy += rand.Gaus(0.0, var_level*Energy);
852 if(Energy < 0.0) Energy = 0.0;
854 if(NumCurrentZone >= FirstStripInLayer && NumCurrentZone < LastStripInLayer+1 && Energy > 0.0) {
855 cluster.
AddStrip(NumCurrentZone, Energy);
856 total_signal += Energy;
867 Double_t mean_fit_pos = 0.0;
869 Int_t NOutLeftBins = 0;
870 Int_t NOutRightBins = 0;
871 if(OutLeftBorder != 0.0) {
872 NOutLeftBins = (Int_t)(
fabs(OutLeftBorder)) + 1;
874 if(OutRightBorder != 0.0) {
875 NOutRightBins = (Int_t)(OutRightBorder - RightZonePos) + 1;
878 Int_t begin_strip_num = cluster.
Strips.at(0);
880 Int_t nstrips = last_strip_num - begin_strip_num + 1;
882 TH1F hist(
"hist_for_fit",
"hist_for_fit", nstrips+NOutLeftBins+NOutRightBins, begin_strip_num-NOutLeftBins, last_strip_num+1+NOutRightBins);
883 Int_t hist_index = 0;
886 Double_t value = cluster.
Signals.at(
i);
887 hist.SetBinContent(hist_index+1+NOutLeftBins, value);
892 if(NOutLeftBins > 0.0) {
893 Double_t first = OutLeftBorder;
894 Double_t last = (Int_t)OutLeftBorder;
895 Double_t S = gausF.Integral(first, last);
896 Double_t Energy = (signal*Gain)*(S/SCluster);
897 Energy += rand.Gaus(0.0, var_level*Energy);
898 if(Energy < 0.0) Energy = 0.0;
899 Double_t value = Energy;
900 hist.SetBinContent(1, value);
903 for(Int_t
i = 1;
i < NOutLeftBins; ++
i) {
905 first = (Int_t)OutLeftBorder+
dist;
907 S = gausF.Integral(first, last);
908 Energy = (signal*Gain)*(S/SCluster);
909 Energy += rand.Gaus(0.0, var_level*Energy);
910 if(Energy < 0.0) Energy = 0.0;
912 hist.SetBinContent(1+
i, value);
918 if(NOutRightBins > 0.0) {
920 for(Int_t
i = hist_index;
i < hist_index+NOutRightBins-1; ++
i) {
922 Double_t first = NStripsInLayer+
dist;
923 Double_t last = first + h;
924 Double_t S = gausF.Integral(first, last);
925 Double_t Energy = (signal*Gain)*(S/SCluster);
926 Energy += rand.Gaus(0.0, var_level*Energy);
927 if(Energy < 0.0) Energy = 0.0;
928 Double_t value = Energy;
929 hist.SetBinContent(1+
i, value);
933 Double_t first = (Int_t)OutRightBorder;
934 Double_t last = first + (OutRightBorder - (Int_t)OutRightBorder);
935 Double_t S = gausF.Integral(first, last);
936 Double_t Energy = (signal*Gain)*(S/SCluster);
937 Energy += rand.Gaus(0.0, var_level*Energy);
938 if(Energy < 0.0) Energy = 0.0;
939 Double_t value = Energy;
940 hist.SetBinContent(hist_index+NOutRightBins, value);
943 TF1* gausFitFunction = 0;
944 TString fit_params =
"WQ0";
946 #ifdef DRAW_REAL_CLUSTER_HISTOGRAMS
951 hist.Fit(
"gaus", fit_params);
952 gausFitFunction = hist.GetFunction(
"gaus");
953 if(gausFitFunction) {
954 mean_fit_pos = gausFitFunction->GetParameter(1);
957 mean_fit_pos = hist.GetMean();
961 mean_fit_pos = hist.GetMean();
977 Int_t n_strip_layers = StripLayers.size();
980 for(Int_t ilayer = 0; ilayer < n_strip_layers; ++ilayer) {
981 StripLayers[ilayer].FindClustersAndStripHits();
985 map<Int_t, Bool_t> UsageStatus_LowerClusters;
986 map<Int_t, Bool_t> UsageStatus_UpperClusters;
989 for(Int_t ilayer = 0; ilayer < n_strip_layers; ++ilayer) {
990 for(Int_t i_strip_hit = 0; i_strip_hit < StripLayers[ilayer].GetNStripHits(); ++i_strip_hit) {
991 StripCluster cluster = StripLayers[ilayer].GetStripClusters()[i_strip_hit];
993 UsageStatus_LowerClusters[cluster.GetUniqueID()] =
false;
996 UsageStatus_UpperClusters[cluster.GetUniqueID()] =
false;
1001 for(Int_t ilayer = 0; ilayer < n_strip_layers; ++ilayer) {
1002 for(Int_t jlayer = ilayer+1; jlayer < n_strip_layers; ++jlayer) {
1005 if(StripLayers[ilayer].GetType() == StripLayers[jlayer].GetType())
continue;
1007 Double_t iAngleRad = StripLayers[ilayer].GetAngleRad();
1008 Double_t jAngleRad = StripLayers[jlayer].GetAngleRad();
1010 if( Abs(iAngleRad - jAngleRad) < 0.01 )
continue;
1012 Double_t i_PiOver2MinusAngleRad = PiOver2() - iAngleRad;
1013 Double_t j_PiOver2MinusAngleRad = PiOver2() - jAngleRad;
1015 for(Int_t i_strip_hit = 0; i_strip_hit < StripLayers[ilayer].GetNStripHits(); ++i_strip_hit) {
1016 for(Int_t j_strip_hit = 0; j_strip_hit < StripLayers[jlayer].GetNStripHits(); ++j_strip_hit) {
1018 Double_t iStripHitPos = StripLayers[ilayer].GetStripHitPos(i_strip_hit);
1019 Double_t jStripHitPos = StripLayers[jlayer].GetStripHitPos(j_strip_hit);
1021 Double_t iStripHitError = StripLayers[ilayer].GetStripHitError(i_strip_hit);
1022 Double_t jStripHitError = StripLayers[jlayer].GetStripHitError(j_strip_hit);
1027 if( Abs(iAngleRad) >= 0.01 && Abs(jAngleRad) >= 0.01 ) {
1029 Double_t iB = StripLayers[ilayer].CalculateStripEquationB(iStripHitPos);
1030 Double_t jB = StripLayers[jlayer].CalculateStripEquationB(jStripHitPos);
1032 xcoord = (jB - iB)/(Tan(i_PiOver2MinusAngleRad) - Tan(j_PiOver2MinusAngleRad));
1033 ycoord = Tan(i_PiOver2MinusAngleRad)*xcoord + iB;
1037 if( Abs(iAngleRad) < 0.01 ) {
1039 if(StripLayers[ilayer].GetStripNumberingOrder() ==
LeftToRight ) {
1040 xcoord = StripLayers[ilayer].GetXLeftStripBorderPoint() + (iStripHitPos-StripLayers[ilayer].GetFirstStripNumber())*StripLayers[ilayer].GetPitch();
1043 xcoord = StripLayers[ilayer].GetXRightStripBorderPoint() - (iStripHitPos-StripLayers[ilayer].GetFirstStripNumber())*StripLayers[ilayer].GetPitch();
1046 Double_t jB = StripLayers[jlayer].CalculateStripEquationB(jStripHitPos);
1047 ycoord = Tan(j_PiOver2MinusAngleRad)*xcoord + jB;
1050 if( Abs(jAngleRad) < 0.01 ) {
1052 if(StripLayers[jlayer].GetStripNumberingOrder() ==
LeftToRight ) {
1053 xcoord = StripLayers[jlayer].GetXLeftStripBorderPoint() + (jStripHitPos-StripLayers[jlayer].GetFirstStripNumber())*StripLayers[jlayer].GetPitch();
1056 xcoord = StripLayers[jlayer].GetXRightStripBorderPoint() - (jStripHitPos-StripLayers[jlayer].GetFirstStripNumber())*StripLayers[jlayer].GetPitch();
1059 Double_t iB = StripLayers[ilayer].CalculateStripEquationB(iStripHitPos);
1060 ycoord = Tan(i_PiOver2MinusAngleRad)*xcoord + iB;
1065 if( StripLayers[ilayer].IsPointInsideStripLayer(xcoord, ycoord) && StripLayers[jlayer].IsPointInsideStripLayer(xcoord, ycoord) ) {
1070 IntersectionPointsX.push_back(xrot);
1071 IntersectionPointsY.push_back(yrot);
1074 Double_t AngleIntersecRad = Abs(iAngleRad - jAngleRad);
1077 Double_t RhombSide1 = (iStripHitError*StripLayers[ilayer].GetPitch())/Sin(AngleIntersecRad);
1078 Double_t RhombSide2 = (jStripHitError*StripLayers[jlayer].GetPitch())/Sin(AngleIntersecRad);
1081 Double_t xerr1 = Sin(Abs(jAngleRad))*RhombSide1;
1082 Double_t xerr2 = Sin(Abs(iAngleRad))*RhombSide2;
1083 Double_t xcoord_err = xerr1 + xerr2;
1086 Double_t yerr1 = Cos(Abs(jAngleRad))*RhombSide1;
1087 Double_t yerr2 = Cos(Abs(iAngleRad))*RhombSide2;
1088 Double_t ycoord_err = yerr1 + yerr2;
1090 IntersectionPointsXErrors.push_back(xcoord_err);
1091 IntersectionPointsYErrors.push_back(ycoord_err);
1095 BmnMatch iStripMatch = StripLayers[ilayer].GetStripMatch((Int_t)iStripHitPos);
1096 BmnMatch jStripMatch = StripLayers[jlayer].GetStripMatch((Int_t)jStripHitPos);
1099 intersection_match.
AddLink(iStripMatch);
1100 intersection_match.
AddLink(jStripMatch);
1102 IntersectionPointMatches.push_back(intersection_match);
1106 BmnMatch iStripDigitNumberMatch = StripLayers[ilayer].GetStripDigitNumberMatch((Int_t)iStripHitPos);
1107 BmnMatch jStripDigitNumberMatch = StripLayers[jlayer].GetStripDigitNumberMatch((Int_t)jStripHitPos);
1109 BmnMatch intersection_digit_num_match;
1110 intersection_digit_num_match.
AddLink(iStripDigitNumberMatch);
1111 intersection_digit_num_match.
AddLink(jStripDigitNumberMatch);
1113 IntersectionPointDigitNumberMatches.push_back(intersection_digit_num_match);
1118 Int_t iLayerClusterSize = StripLayers[ilayer].GetStripHitClusterSize(i_strip_hit);
1119 Int_t jLayerClusterSize = StripLayers[jlayer].GetStripHitClusterSize(j_strip_hit);
1122 Double_t iLayerStripPosition = StripLayers[ilayer].ConvertPointToStripPosition(xcoord, ycoord);
1123 Double_t jLayerStripPosition = StripLayers[jlayer].ConvertPointToStripPosition(xcoord, ycoord);
1126 IntersectionPoints_LowerLayerClusterSize.push_back(iLayerClusterSize);
1127 IntersectionPoints_UpperLayerClusterSize.push_back(jLayerClusterSize);
1128 IntersectionPoints_LowerLayerStripPosition.push_back(iLayerStripPosition);
1129 IntersectionPoints_UpperLayerStripPosition.push_back(jLayerStripPosition);
1130 IntersectionPoints_LowerLayerStripTotalSignal.push_back(StripLayers[ilayer].GetStripHitTotalSignal(i_strip_hit));
1131 IntersectionPoints_UpperLayerStripTotalSignal.push_back(StripLayers[jlayer].GetStripHitTotalSignal(j_strip_hit));
1132 LowerClusters.push_back(StripLayers[ilayer].GetStripClusters().at(i_strip_hit));
1133 UpperClusters.push_back(StripLayers[jlayer].GetStripClusters().at(j_strip_hit));
1134 UsageStatus_LowerClusters[StripLayers[ilayer].GetStripClusters()[i_strip_hit].GetUniqueID()] =
true;
1135 UsageStatus_UpperClusters[StripLayers[jlayer].GetStripClusters()[j_strip_hit].GetUniqueID()] =
true;
1138 IntersectionPoints_LowerLayerClusterSize.push_back(jLayerClusterSize);
1139 IntersectionPoints_UpperLayerClusterSize.push_back(iLayerClusterSize);
1140 IntersectionPoints_LowerLayerStripPosition.push_back(jLayerStripPosition);
1141 IntersectionPoints_UpperLayerStripPosition.push_back(iLayerStripPosition);
1142 IntersectionPoints_LowerLayerStripTotalSignal.push_back(StripLayers[jlayer].GetStripHitTotalSignal(j_strip_hit));
1143 IntersectionPoints_UpperLayerStripTotalSignal.push_back(StripLayers[ilayer].GetStripHitTotalSignal(i_strip_hit));
1144 LowerClusters.push_back(StripLayers[jlayer].GetStripClusters().at(j_strip_hit));
1145 UpperClusters.push_back(StripLayers[ilayer].GetStripClusters().at(i_strip_hit));
1146 UsageStatus_LowerClusters[StripLayers[jlayer].GetStripClusters()[j_strip_hit].GetUniqueID()] =
true;
1147 UsageStatus_UpperClusters[StripLayers[ilayer].GetStripClusters()[i_strip_hit].GetUniqueID()] =
true;
1160 for(
auto it : UsageStatus_LowerClusters)
1163 if(it.second ==
false)
1165 for(Int_t ilayer = 0; ilayer < n_strip_layers; ++ilayer)
1169 for (
size_t icluster = 0; icluster < StripLayers[ilayer].GetStripClusters().size(); ++icluster) {
1170 StripCluster cluster = StripLayers[ilayer].GetStripClusters()[icluster];
1172 if((
int)cluster.GetUniqueID() != it.first)
continue;
1177 if(StripLayers[ilayer].GetStripNumberingOrder() ==
LeftToRight ) {
1178 xcoord = StripLayers[ilayer].GetXLeftStripBorderPoint() + cluster.
MeanPosition*StripLayers[ilayer].GetPitch();
1181 xcoord = StripLayers[ilayer].GetXRightStripBorderPoint() - cluster.
MeanPosition*StripLayers[ilayer].GetPitch();
1184 ycoord = (StripLayers[ilayer].GetYMinLayer() + StripLayers[ilayer].GetYMaxLayer())*0.5;
1186 Double_t xerr = cluster.
Error*StripLayers[ilayer].GetPitch();
1187 Double_t yerr = StripLayers[ilayer].GetYSize()*0.5;
1199 PseudoIntersectionsX.push_back(xcoord);
1200 PseudoIntersectionsY.push_back(ycoord);
1202 PseudoIntersectionsXErrors.push_back(xerr);
1203 PseudoIntersectionsYErrors.push_back(yerr);
1205 PseudoIntersections_LowerLayerClusterSize.push_back(cluster.
Width);
1206 PseudoIntersections_UpperLayerClusterSize.push_back(0);
1208 PseudoIntersections_LowerLayerStripPosition.push_back(StripLayers[ilayer].ConvertPointToStripPosition(xcoord, ycoord));
1209 PseudoIntersections_UpperLayerStripPosition.push_back(0);
1211 PseudoIntersections_LowerLayerStripTotalSignal.push_back(cluster.
TotalSignal);
1212 PseudoIntersections_UpperLayerStripTotalSignal.push_back(0);
1214 PseudoIntersectionMatches.push_back(StripLayers[ilayer].GetStripMatch((Int_t)cluster.
MeanPosition));
1215 PseudoIntersectionDigitNumberMatches.push_back(StripLayers[ilayer].GetStripDigitNumberMatch((Int_t)cluster.
MeanPosition));
1217 LowerClusters_PseudoIntersections.push_back(cluster);
1218 UpperClusters_PseudoIntersections.push_back(
StripCluster());
1228 Int_t ilayer = layerA_index;
1229 Int_t jlayer = layerB_index;
1231 if(StripLayers[ilayer].GetType() == StripLayers[jlayer].GetType())
return false;
1233 Double_t iAngleRad = StripLayers[ilayer].GetAngleRad();
1234 Double_t jAngleRad = StripLayers[jlayer].GetAngleRad();
1236 if( Abs(iAngleRad - jAngleRad) < 0.01 )
return false;
1238 Double_t i_PiOver2MinusAngleRad = PiOver2() - iAngleRad;
1239 Double_t j_PiOver2MinusAngleRad = PiOver2() - jAngleRad;
1241 Double_t iStripHitPos = strip_pos_layerA;
1242 Double_t jStripHitPos = strip_pos_layerB;
1247 if( Abs(iAngleRad) >= 0.01 && Abs(jAngleRad) >= 0.01 ) {
1249 Double_t iB = StripLayers[ilayer].CalculateStripEquationB(iStripHitPos);
1250 Double_t jB = StripLayers[jlayer].CalculateStripEquationB(jStripHitPos);
1252 xcoord = (jB - iB)/(Tan(i_PiOver2MinusAngleRad) - Tan(j_PiOver2MinusAngleRad));
1253 ycoord = Tan(i_PiOver2MinusAngleRad)*xcoord + iB;
1257 if( Abs(iAngleRad) < 0.01 ) {
1259 if(StripLayers[ilayer].GetStripNumberingOrder() ==
LeftToRight ) {
1260 xcoord = StripLayers[ilayer].GetXLeftStripBorderPoint() + (iStripHitPos-StripLayers[ilayer].GetFirstStripNumber())*StripLayers[ilayer].GetPitch();
1263 xcoord = StripLayers[ilayer].GetXRightStripBorderPoint() - (iStripHitPos-StripLayers[ilayer].GetFirstStripNumber())*StripLayers[ilayer].GetPitch();
1266 Double_t jB = StripLayers[jlayer].CalculateStripEquationB(jStripHitPos);
1267 ycoord = Tan(j_PiOver2MinusAngleRad)*xcoord + jB;
1270 if( Abs(jAngleRad) < 0.01 ) {
1272 if(StripLayers[jlayer].GetStripNumberingOrder() ==
LeftToRight ) {
1273 xcoord = StripLayers[jlayer].GetXLeftStripBorderPoint() + (jStripHitPos-StripLayers[jlayer].GetFirstStripNumber())*StripLayers[jlayer].GetPitch();
1276 xcoord = StripLayers[jlayer].GetXRightStripBorderPoint() - (jStripHitPos-StripLayers[jlayer].GetFirstStripNumber())*StripLayers[jlayer].GetPitch();
1279 Double_t iB = StripLayers[ilayer].CalculateStripEquationB(iStripHitPos);
1280 ycoord = Tan(i_PiOver2MinusAngleRad)*xcoord + iB;
1291 if( StripLayers[ilayer].IsPointInsideStripLayer(xcoord, ycoord) && StripLayers[jlayer].IsPointInsideStripLayer(xcoord, ycoord) ) {