27 #include "TLinearFitter.h"
51 if (plane != NULL )plane->
print();
60 Info(
"process",
"Event %ld",
nEvents);
66 for(Int_t i = 0; i<
NTRACK; i++)
81 for(
int a = 0;a<4;a++){
82 for(
int b =0;b<4;b++){
87 Corig[a][b]=(2*inc)*(2*inc)+theta*theta;
94 for(Int_t plane = 0; plane<2; plane++){
97 Int_t counter = 0, color = 51, nClLeftOld = 0, nTry = 0;
104 Info(
"process",
"Look for seeds on plane %d %d",plane,counter);
107 Int_t iMin =
NCHAN, iMax = 0, jMin =
NADC, jMax = 0;
108 Int_t nClLeft =
GetMapEdges(clArray,plane,iMin,iMax,jMin,jMax);
111 Info(
"process",
"Only %d clusters left",nClLeft);
114 if(nClLeft==nClLeftOld){
116 Info(
"process",
"No new tracks, %d clusters left",nClLeft);
119 nClLeftOld = nClLeft;
124 for(Int_t side = 0; side<4; side++){
125 Int_t
type = -1, index0 = 0, index1 = 0, imax0 = 0, imax1 = 0, sign = 0;
127 case 0: type =
TCLUSTER; index0 = jMax; sign = -1;
break;
128 case 1: type =
CCLUSTER; index0 = iMax; sign = -1;
break;
129 case 2: type =
CCLUSTER; index0 = iMin; sign = +1;
break;
130 case 3: type =
TCLUSTER; index0 = jMin; sign = +1;
break;
134 index1 = index0 + sign;
136 if(index0<0 || index0>=
NADC)
continue;
137 while(
NTcl[index1]<=0 && TMath::Abs(index1-index0) < 50) index1 += sign;
138 if(index1<0 || index1>=
NADC)
continue;
139 imax0 =
NTcl[index0];
140 imax1 =
NTcl[index1];
142 if(index0<0 || index0>=
NCHAN)
continue;
143 while(
NCcl[index1]<=0 && TMath::Abs(index1-index0) < 50) index1 += sign;
144 if(index1<0 || index1>=
NCHAN)
continue;
145 imax0 =
NCcl[index0];
146 imax1 =
NCcl[index1];
150 if(TMath::Abs(index1-index0) >= 10){
151 for(Int_t i0 = 0; i0<imax0; i0++){
155 else icl0 =
Tcl[index0][i0];
159 if(
gHarpoDebug>1) Info(
"process",
"side %d: i0 = %d, i1 = %d",side,index0,index1);
165 Info(
"process",
"side %d => index = %d %d (%d %d)",side,index0,index1,imax0,imax1);
166 for(Int_t i0 = 0; i0<imax0; i0++){
170 else icl0 =
Tcl[index0][i0];
174 if(
used[icl0] != -1000)
continue;
176 for(Int_t i1 = 0; i1<imax1; i1++){
179 else icl1 =
Tcl[index1][i1];
182 if(
used[icl1] != -1000)
continue;
189 if(
gHarpoDebug>1) Info(
"process",
"nTry = %d",nTry);
190 FindTrack(clArray,icl0,icl1,Corig,plane,color);
202 fHist[plane]->GetXaxis()->SetRangeUser(0,511);
203 fHist[plane]->GetYaxis()->SetRangeUser(0,304);
212 for(Int_t itr = 0; itr<
NTRACK; itr++){
225 iMin =
NCHAN+1; iMax = -1; jMin =
NADC+1; jMax = -1;
226 Int_t nCl = clArray->GetEntries();
227 for(Int_t icl = 0; icl<nCl; icl++){
231 if(
used[icl] != -1000)
continue;
234 if(cluster->
GetPlane() != plane)
continue;
240 if(index>jMax) jMax = index;
241 if(index<jMin) jMin = index;
244 if(index>iMax) iMax = index;
245 if(index<iMin) iMin = index;
251 Info(
"GetMapEdges",
"%d %d %d %d (%d)",iMin, iMax, jMin, jMax, nClLeft);
260 Info(
"FindTrack",
"%d %d",icl0,icl1);
272 if(type != 0 && type != 1)
return;
274 Double_t x0 = cluster0->
GetMean();
275 Double_t x1 = cluster1->
GetMean();
276 Double_t index0 = cluster0->
GetIndex();
277 Double_t index1 = cluster1->
GetIndex();
278 Int_t q0 = cluster0->
GetQ();
280 if(TMath::Abs(x0-x1)>15)
return;
283 Xorig[1-
type][0] = index0;
284 Xorig[2+
type][0] = x1-x0;
285 Xorig[3-
type][0] = index1-index0;
287 Info(
"process",
"FindNext(%f,%f,%f,%f,%d)",Xorig[0][0],Xorig[1][0],Xorig[2][0],Xorig[3][0],
fNtr[plane]);
288 if(TMath::Abs(Xorig[2+type][0]/Xorig[3-type][0])>1.5)
return;
291 TArrow* arrow =
new TArrow();
292 arrow->SetLineColor(kBlack);
293 arrow->SetFillColor(color);
295 arrow->DrawArrow(Xorig[0][0],Xorig[1][0],Xorig[0][0]+Xorig[2][0],Xorig[1][0]+Xorig[3][0],0.01);
296 fHist[plane]->GetXaxis()->SetRangeUser(Xorig[0][0] - 10,Xorig[0][0] + 10);
297 fHist[plane]->GetYaxis()->SetRangeUser(Xorig[1][0] - 10,Xorig[1][0] + 10);
307 Double_t angleorig = TMath::ATan2(Xorig[2][0],Xorig[3][0]);
308 FindNext(Xorig,Corig,angleorig,
fNtr[plane],clArray,0,1,plane,color,3001,0,q0);
310 Info(
"process",
"Found new track %d",
fNtr[plane]);
312 Int_t nClTr =
AddTrack(clArray,plane);
317 Info(
"process",
"Track %d: %d clusters",
fNtr[plane]-1,nClTr);
333 Warning(
"AddTrack",
"Wrong track number %d",
fNtr[plane]);
346 Info(
"AddTrack",
"New seed: (%g, %g, %g, %g)",X[0][0],X[1][0],X[2][0],X[3][0]);
349 for(
int a = 0;a<4;a++){
350 for(
int b =0;b<4;b++){
361 TArrayI* arr =
new TArrayI(2000);
362 Double_t angleorig = TMath::ATan2(X[2][0],X[3][0]);
363 arr =
FindNextClosest(X, Corig, angleorig,fNtr[plane], clArray, arr, 1,plane, kGray, 3002, 2, -1);
365 for(Int_t i = 0; i<2000; i++){
366 if(arr->At(i)-1<0)
break;
375 if(
fGraphs[plane][fNtr[plane]])
376 fGraphs[plane][fNtr[plane]]->Delete();
377 fGraphs[plane][fNtr[plane]] =
new TGraph();
379 fGraphs[plane][fNtr[plane]]->SetName(Form(
"TrackX%d",fNtr[plane]));
381 fGraphs[plane][fNtr[plane]]->SetName(Form(
"TrackY%d",fNtr[plane]));
383 for(Int_t i = 0; i<n; i++){
384 Int_t icl = arr->At(i) - 1;
387 if(icl>=4000)
continue;
388 fTrack[icl] = fNtr[plane];
390 Double_t mean = cluster->
GetMean();
395 fGraphs[plane][fNtr[plane]]->SetPoint(
fGraphs[plane][fNtr[plane]]->GetN(),index,mean);
399 fGraphs[plane][fNtr[plane]]->SetPoint(
fGraphs[plane][fNtr[plane]]->GetN(),mean,index);
421 Double_t fThetaMin = -0.5;
422 Double_t fDminForward = 20, fDminTransverse = 2;
427 Info(
"SpliceTracks",
"plane %d: Ntracks = %d",plane,
fNtr[plane]);
429 Double_t x[
fNtr[plane]][2], z[fNtr[plane]][2], th[fNtr[plane]][2];
430 Double_t px[fNtr[plane]][2], pz[fNtr[plane]][2];
431 Int_t spliced[fNtr[plane]];
433 for(Int_t itr = 0; itr<fNtr[plane]; itr++){
435 if(
fGraphs[plane][itr]->GetN()<=fOffset+1)
continue;
436 fGraphs[plane][itr]->GetPoint(0 + fOffset,z[itr][0],x[itr][0]);
438 fGraphs[plane][itr]->GetPoint(1 + fOffset,xTmp,zTmp);
439 px[itr][0] = xTmp - x[itr][0];
440 pz[itr][0] = zTmp - z[itr][0];
441 fGraphs[plane][itr]->GetPoint(
fGraphs[plane][itr]->GetN()-1-fOffset,z[itr][1],x[itr][1]);
442 fGraphs[plane][itr]->GetPoint(
fGraphs[plane][itr]->GetN()-2-fOffset,xTmp,zTmp);
443 px[itr][1] = xTmp - x[itr][1];
444 pz[itr][1] = zTmp - z[itr][1];
452 for(Int_t itr1 = 0; itr1<fNtr[plane]; itr1++){
453 if(spliced[itr1] != itr1)
continue;
454 for(Int_t itr2 = itr1+1; itr2<fNtr[plane]; itr2++){
455 if(spliced[itr2] != itr2)
continue;
457 Int_t iMin = -1, jMin = -1;
458 for(Int_t i = 0; i<2; i++){
459 for(Int_t j = 0; j<2; j++){
460 Double_t d = TMath::Abs(x[itr1][i]-x[itr2][j]) + TMath::Abs(z[itr1][i]-z[itr2][j]);
469 Double_t dx1 = x[itr1][1-iMin] - x[itr1][iMin];
470 Double_t dz1 = z[itr1][1-iMin] - z[itr1][iMin];
471 Double_t dx2 = x[itr2][1-jMin] - x[itr2][jMin];
472 Double_t dz2 = z[itr2][1-jMin] - z[itr2][jMin];
475 Double_t cross = dx1*dx2 + dz1*dz2;
476 Double_t d1 = dx1*dx1 + dz1*dz1;
477 Double_t d2 = dx2*dx2 + dz2*dz2;
478 Double_t sd1 = TMath::Sqrt(d1);
479 Double_t sd2 = TMath::Sqrt(d2);
480 if(d1*d2 == 0)
continue;
481 Double_t costheta = cross/TMath::Sqrt(d1*d2);
485 Info(
"SpliceTracks",
"%d %d: (d=%g,theta=%g)",itr1,itr2,dMin,costheta);
487 if(costheta > fThetaMin)
continue;
493 Double_t dForward1 = TMath::Abs((x[itr1][iMin] - x[itr2][jMin])*dx1 +
494 (z[itr1][iMin] - z[itr2][jMin])*dz1)/sd1;
495 Double_t dTransverse1 = TMath::Abs((x[itr1][iMin] - x[itr2][jMin])*dz1 -
496 (z[itr1][iMin] - z[itr2][jMin])*dx1)/sd1;
497 Double_t dForward2 = TMath::Abs((x[itr1][iMin] - x[itr2][jMin])*dx2 +
498 (z[itr1][iMin] - z[itr2][jMin])*dz2)/sd2;
499 Double_t dTransverse2 = TMath::Abs((x[itr1][iMin] - x[itr2][jMin])*dz2 -
500 (z[itr1][iMin] - z[itr2][jMin])*dx2)/sd2;
503 Info(
"SpliceTracks",
"%d %d: (dF=%g,dT=%g,theta=%g)",itr1,itr2,dForward1,dTransverse1,costheta);
506 Double_t xRect[4] = {x[itr1][iMin] + ( fDminForward*dx1/sd1 - fDminTransverse*dz1/sd1),
507 x[itr1][iMin] + (- fDminForward*dx1/sd1 - fDminTransverse*dz1/sd1),
508 x[itr1][iMin] + (- fDminForward*dx1/sd1 + fDminTransverse*dz1/sd1),
509 x[itr1][iMin] + ( fDminForward*dx1/sd1 + fDminTransverse*dz1/sd1)};
510 Double_t zRect[4] = {z[itr1][iMin] + ( fDminForward*dz1/sd1 + fDminTransverse*dx1/sd1),
511 z[itr1][iMin] + (- fDminForward*dz1/sd1 + fDminTransverse*dx1/sd1),
512 z[itr1][iMin] + (- fDminForward*dz1/sd1 - fDminTransverse*dx1/sd1),
513 z[itr1][iMin] + ( fDminForward*dz1/sd1 - fDminTransverse*dx1/sd1)};
515 std::cout <<
"----------------" << std::endl;
516 for(Int_t i = 0; i<4; i++) std::cout << xRect[i] <<
" ";
517 std::cout << std::endl;
518 for(Int_t i = 0; i<4; i++) std::cout << zRect[i] <<
" ";
519 std::cout << std::endl;
520 TPolyLine* rect =
new TPolyLine(4,zRect,xRect,
"FL");
521 rect->SetFillColor(kRed);
522 rect->SetLineWidth(2);
525 fHist[plane]->GetYaxis()->SetRangeUser(x[itr1][iMin]-20,x[itr1][iMin]+20);
526 fHist[plane]->GetXaxis()->SetRangeUser(z[itr1][iMin]-20,z[itr1][iMin]+20);
534 if(((dForward1<fDminForward && dTransverse1<fDminTransverse) ||
535 (dForward2<fDminForward && dTransverse2<fDminTransverse)) &&
542 TArrow* arrow1 =
new TArrow();
543 arrow1->SetLineColor(kBlue);
544 arrow1->DrawArrow(z[itr1][iMin],x[itr1][iMin],z[itr1][iMin] + dz1,x[itr1][iMin]+dx1);
545 TArrow* arrow2 =
new TArrow();
546 arrow2->SetLineColor(kRed);
547 arrow2->DrawArrow(x[itr2][jMin],x[itr2][jMin],z[itr2][jMin]+dz2,x[itr2][jMin]+dx2);
549 gSystem->Sleep(1000);
550 fHist[plane]->GetXaxis()->SetRangeUser(90,420);
551 fHist[plane]->GetYaxis()->SetRangeUser(0,288);
556 spliced[itr2] = spliced[itr1];
557 x[itr1][iMin] = x[itr2][1-jMin];
558 z[itr1][iMin] = z[itr2][1-jMin];
559 px[itr1][iMin] = px[itr2][1-jMin];
560 pz[itr1][iMin] = pz[itr2][1-jMin];
561 th[itr1][iMin] = th[itr2][1-jMin];
564 Info(
"SpliceTrack",
"Splicing %d %d (dF=%g,dT=%g,theta=%g)",itr1,itr2,dForward1,dTransverse1,costheta);
568 Double_t a1 = ((x[itr1][iMin] - x[itr2][jMin])*dz2 - (z[itr1][iMin] - z[itr2][jMin])*dx2)/(dz1*dx2 - dz2*dx1);
569 Double_t a2 = ((x[itr1][iMin] - x[itr2][jMin])*dz1 - (z[itr1][iMin] - z[itr2][jMin])*dx1)/(dz1*dx2 - dz2*dx1);
572 Info(
"SpliceTracks",
"%d %d: %g/%g %g/%g",itr1,itr2,a1,sd1,a2,sd2);
576 if( TMath::Abs(a1 - 0.5) > 1)
continue;
577 if( TMath::Abs(a2 - 0.5) > 1)
continue;
578 if( (( a1*sd1 < -40 ) || (a1*sd1 > sd1 + 40)) &&
579 (( a2*sd2 < -40 ) || (a2*sd2 > sd2 + 40)))
585 TArrow* arrow1 =
new TArrow();
586 arrow1->SetLineColor(kBlue+3);
587 arrow1->DrawArrow(x[itr1][iMin],z[itr1][iMin],x[itr1][iMin] + a1*dx1,z[itr1][iMin]+a1*dz1);
588 TArrow* arrow2 =
new TArrow();
589 arrow2->SetLineColor(kRed+3);
590 arrow2->DrawArrow(x[itr2][jMin],z[itr2][jMin],x[itr2][jMin]+a2*dx2,z[itr2][jMin]+a2*dz2);
592 gSystem->Sleep(1000);
593 fHist[plane]->GetXaxis()->SetRangeUser(90,420);
594 fHist[plane]->GetYaxis()->SetRangeUser(0,288);
597 Double_t dMax = a1*sd1;
598 Int_t i0 = iMin, tr0 = itr1;
599 Int_t i1 = iMin, tr1 = itr1;
600 if((1-a1)*sd1 > dMax){
614 if((1-a2)*sd2 > dMax){
622 spliced[itr2] = spliced[itr1];
623 x[itr1][0] = x[tr0][i0];
624 z[itr1][0] = z[tr0][i0];
625 px[itr1][0] = px[tr0][i0];
626 pz[itr1][0] = pz[tr0][i0];
627 th[itr1][0] = th[tr0][i0];
628 x[itr1][1] = x[tr1][i1];
629 z[itr1][1] = z[tr1][i1];
630 px[itr1][1] = px[tr1][i1];
631 pz[itr1][1] = pz[tr1][i1];
632 th[itr1][1] = th[tr1][i1];
635 Info(
"SpliceTrack",
"** Splicing %d %d",itr1,itr2);
640 Info(
"SpliceTracks",
"** %d spliced (%d)",nSplice,nLoops);
641 if(nSplice == 0)
break;
649 Int_t newTrIndex[fNtr[plane]];
650 Int_t nClTr[fNtr[plane]];
651 Double_t qTr[fNtr[plane]];
653 for(Int_t itr = 0; itr<fNtr[plane]; itr++){
656 if(spliced[itr] != itr){
657 Int_t itrNew = spliced[itr];
658 while(spliced[itrNew] != itrNew)
659 itrNew = spliced[itrNew];
660 spliced[itr] = itrNew;
661 Int_t n =
fGraphs[plane][itr]->GetN();
662 for(Int_t k = 0; k<n; k++){
664 fGraphs[plane][itr]->GetPoint(k,xTmp,zTmp);
665 fGraphs[plane][itrNew]->SetPoint(
fGraphs[plane][itrNew]->GetN(),xTmp,zTmp);
673 newTrIndex[itr] = -1;
681 Int_t nCl = clArray->GetEntries();
683 Info(
"Splice",
"ntot = %d",nCl);
686 for(Int_t icl = 0; icl<nCl; icl++){
687 Int_t oldTrack =
fTrack[icl];
688 if(oldTrack<0 || oldTrack>=fNtr[plane])
continue;
689 Int_t splicedTrack = spliced[oldTrack];
692 nClTr[splicedTrack]++;
696 for(Int_t itr1 = 0; itr1<fNtr[plane]; itr1++){
698 for(Int_t itr2 = itr1+1; itr2<fNtr[plane]; itr2++){
700 Info(
"SpliceTracks",
"%d => %d (%g %g %g)",itr2,spliced[itr1],
fQcommon[itr1][itr1],
fQcommon[itr2][itr2],
fQcommon[itr1][itr2]);
704 Info(
"SpliceTracks",
"OK");
705 spliced[itr2] = spliced[itr1];
706 nClTr[spliced[itr1]] += nClTr[itr2];
707 qTr[spliced[itr1]] +=
fQcommon[itr2][itr2];
711 for(Int_t itr1 = 0; itr1<fNtr[plane]; itr1++){
712 if(itr1 == spliced[itr1] && nClTr[itr1]>
fNclMin2 && qTr[itr1]>20000){
713 newTrIndex[itr1] = nFinal;
716 newTrIndex[itr1] = -1;
720 for(Int_t icl = 0; icl<nCl; icl++){
722 if(cluster->
GetPlane() != plane)
continue;
725 Int_t oldTrack =
fTrack[icl];
726 if(oldTrack<0 || oldTrack>=fNtr[plane])
continue;
727 Int_t sp = spliced[oldTrack];
729 if(newTrIndex[sp]<0)
continue;
733 for(Int_t itr = 0; itr<fNtr[plane]; itr++){
735 Info(
"Splice",
"nCl(%d) = %d",itr,nClTr[itr]);
737 if(spliced[itr] != itr)
continue;
738 if(newTrIndex[itr] == -1)
continue;
740 Info(
"Splice",
"nCl(%d) = %d / %d",itr,nClTr[itr],
fGraphs[plane][itr]->GetN());
741 Double_t x0 = -1000, slope = -1000, chi2 = -1000;
742 if(
fGraphs[plane][itr]->GetN()>10){
743 TLinearFitter* linfit =
new TLinearFitter(1,
"pol1",
"D");
744 linfit->AssignData(
fGraphs[plane][itr]->GetN(), 1,
fGraphs[plane][itr]->GetX(),
fGraphs[plane][itr]->GetY());
745 linfit->EvalRobust();
746 x0 = linfit->GetParameter(0);
747 slope = linfit->GetParameter(1);
748 chi2 = linfit->GetChisquare();
761 Int_t typeTrack = 0, IdTrackMatching = 0;
762 Double_t xSt = x[itr][0], zSt = z[itr][0], thSt = th[itr][0], xEnd = x[itr][1], zEnd = z[itr][1], thEnd = th[itr][1];
763 Double_t pxSt = px[itr][0], pzSt = pz[itr][0], pxEnd = px[itr][1], pzEnd = pz[itr][1];
765 Info(
"SpliceTracks",
"new Track: %d (%d cl, q=%g) => (%g, %g, %g) (%g, %g, %g)", newTrIndex[itr],nClTr[itr], qTrTmp, xSt, zSt, thSt, xEnd, zEnd, thEnd);
766 HarpoRecoKalmanTracks* tr =
new HarpoRecoKalmanTracks(newTrIndex[itr],Int_t(qTrTmp),slope,x0,chi2,typeTrack,IdTrackMatching,nClTr[itr],plane,xSt,zSt,pxSt,pzSt,xEnd,zEnd,pxEnd,pzEnd);
769 fId[newTrIndex[itr]][plane] = itr;
773 Info(
"SpliceTracks",
"Done");
785 Int_t nCl = clArray->GetEntries();
787 for(Int_t icl = 0; icl<nCl; icl++){
789 if(cluster->
GetPlane() != plane)
continue;
795 for(Int_t k = ind2; k< ind2 + width; k++){
796 Int_t i = -1, j = -1;
797 if(type ==
CCLUSTER){ i = ind; j = k;}
798 if(type ==
TCLUSTER){ i = k; j = ind;}
799 if(i<0 || i>=
NCHAN)
continue;
800 if(j<0 || j>=
NADC)
continue;
801 Double_t q = map->
GetData(i,j);
804 ULong64_t test = 1 << itr;
806 if(
fMapTmp[i][j] & test)
continue;
813 Info(
"GetQtrack",
"Plane %d, Track %d: Q = %g/%g",plane,itr,qTr,
fQtot);
823 Int_t nCl = clArray->GetEntries();
827 for(Int_t i = 0; i<
NCHAN; i++){
828 for(Int_t j = 0; j<
NADC; j++){
834 for(Int_t icl = 0; icl<nCl; icl++){
836 if(cluster->
GetPlane() != plane)
continue;
841 Int_t nTr = cluster->
GetNtr();
842 for(Int_t iTr = 0; iTr<nTr;iTr++){
844 for(Int_t k = ind2; k< ind2 + width; k++){
845 Int_t i = -1, j = -1;
846 if(type ==
CCLUSTER){ i = ind; j = k;}
847 if(type ==
TCLUSTER){ i = k; j = ind;}
848 if(i<0 || i>=NCHAN)
continue;
849 if(j<0 || j>=
NADC)
continue;
850 Double_t q = map->
GetData(i,j);
853 ULong64_t test = 1 << itr;
855 if(used4Q[i][j] & test)
continue;
858 used4Q[i][j] |= test;
864 for(Int_t i = 0; i<
fNtr[plane]; i++){
865 for(Int_t j = 0; j<fNtr[plane]; j++){
869 for(Int_t i = 0; i<
NCHAN; i++){
870 for(Int_t j = 0; j<
NADC; j++){
871 Double_t q = map->
GetData(i,j);
874 for(Int_t itr1 = 0; itr1<fNtr[plane]; itr1++){
875 ULong64_t test1 = 1 << itr1;
876 if(!(used4Q[i][j] & test1))
continue;
877 for(Int_t itr2 = itr1; itr2<fNtr[plane]; itr2++){
878 ULong64_t test2 = 1 << itr2;
879 if(!(used4Q[i][j] & test2))
continue;
887 for(Int_t i = 0; i<fNtr[plane]; i++){
888 for(Int_t j = 0; j<fNtr[plane]; j++){
891 std::cout << std::endl;
928 hNcl =
new TH1F(
"hNcl",
"hNcl",500,0,500);
930 for(Int_t itr = 0; itr<
NTRACK; itr++){
936 const Int_t nbinsDist = 200;
937 Double_t xminDist = 1;
938 Double_t xmaxDist = 10000;
939 Double_t logxminDist = TMath::Log(xminDist);
940 Double_t logxmaxDist = TMath::Log(xmaxDist);
941 Double_t binwidthDist = (logxmaxDist-logxminDist)/nbinsDist;
942 Double_t xbinsDist[nbinsDist+1];
943 xbinsDist[0] = xminDist;
944 for (Int_t i=1;i<=nbinsDist;i++)
945 xbinsDist[i] = TMath::Exp(logxminDist+i*binwidthDist);
948 const Int_t nbinsTheta = 100;
949 Double_t xminTheta = 1e-5;
950 Double_t xmaxTheta = 2;
951 Double_t logxminTheta = TMath::Log(xminTheta);
952 Double_t logxmaxTheta = TMath::Log(xmaxTheta);
953 Double_t binwidthTheta = (logxmaxTheta-logxminTheta)/nbinsTheta;
954 Double_t xbinsTheta[nbinsTheta+1];
955 xbinsTheta[0] = xminTheta;
956 for (Int_t i=1;i<=nbinsTheta;i++)
957 xbinsTheta[i] = TMath::Exp(logxminTheta+i*binwidthTheta);
959 hDistTracks =
new TH1F(
"hDistTracks",
"",nbinsDist,xbinsDist);
960 hDistMinTracks =
new TH1F(
"hDistMinTracks",
"",nbinsDist,xbinsDist);
961 hThetaTracks =
new TH1F(
"hThetaTracks",
"",nbinsTheta,xbinsTheta);
963 hThetaTracksNtr =
new TH2F(
"hThetaTracksNtr",
"",nbinsDist,xbinsDist,nbinsTheta,xbinsTheta);
964 for(Int_t i = 0; i<10; i++)
965 hThetaDist[i] =
new TH2F(Form(
"hThetaDist%d",i),
"",nbinsDist,xbinsDist,nbinsTheta,xbinsTheta);
968 const Int_t nbinsQ = 100;
969 Double_t xminQ = 1e-3;
970 Double_t xmaxQ = 1e3;
971 Double_t logxminQ = TMath::Log(xminQ);
972 Double_t logxmaxQ = TMath::Log(xmaxQ);
973 Double_t binwidthQ = (logxmaxQ-logxminQ)/nbinsQ;
974 Double_t xbinsQ[nbinsQ+1];
976 for (Int_t i=1;i<=nbinsQ;i++)
977 xbinsQ[i] = TMath::Exp(logxminQ+i*binwidthQ);
978 fHistQQtestVsDist =
new TH2F(
"fHistQQtestVsDist",
"",nbinsQ,xbinsQ,nbinsDist,xbinsDist);
979 fHistQQmin =
new TH2F(
"fHistQQmin",
"",nbinsQ,xbinsQ,nbinsDist,xbinsDist);
991 for(Int_t i = 0; i<10; i++){
void print()
Ovreloaded medod whic do all job.
void SpliceTracks(Int_t plane)
HarpoRecoEvent * GetRecoEvent()
Double_t fStartPointZ[NTRACK][10]
Int_t InitPlane(TClonesArray *clArray, Int_t plane)
void AddIdClusterTrack(Int_t val)
Double_t GetData(Int_t i, Int_t j)
Set/Get Data Cell.
Double_t GetQtrack(Int_t itr, Int_t plane)
Double_t fQcommon[NTRACK][NTRACK]
Object containing the reconstruction information for one event (with array of HarpoRecoClusters Harpo...
Int_t AddTrack(TClonesArray *clArray, Int_t plane)
void AddTracks(HarpoRecoTracks *val)
Int_t InitPlane(TClonesArray *clArray, Int_t plane)
TArrayI * FindNextClosest(TMatrixD X, TMatrixD C, Double_t angle, Int_t Ntr, TClonesArray *clArray, TArrayI *arr, Int_t ncl, Int_t plane, Int_t color, Int_t fill, Int_t smooth, Int_t q)
Int_t GetMapEdges(TClonesArray *clArray, Int_t plane, Int_t &iMin, Int_t &iMax, Int_t &jMin, Int_t &jMax)
Track finder with Kalman filter.
TGraph * fGraphs[2][NTRACK]
void RemoveAllClusterTrack()
Double_t GetQtracks(Int_t plane)
HarpoRecoTracks object, Obtained with Kalman filter.
void Save(char *mode=NULL)
Cluster object, containing position, charge and quality information.
TFile * OpenHistFile(const char *ananame)
unpacked dcc data The class contains the data map for DCC or Feminos The data is stored as a 2D TMatr...
ULong64_t fMapTmp[NCHAN][NADC]
Int_t GetIdClusterTrack()
TH1F * hNcl
Redefine empty default.
Double_t fStartDirX[NTRACK][10]
Int_t fStartIndex[NTRACK]
HarpoEventHeader * GetHeader()
TArrayI * FindNext(TMatrixD X, TMatrixD C, Double_t angle, Int_t Ntr, TClonesArray *clArray, TArrayI *arr, Int_t ncl, Int_t plane, Int_t color, Int_t fill, Int_t smooth, Int_t q)
Double_t fStartPointX[NTRACK][10]
const ULong_t gkNDetectors
Double_t fStartDirZ[NTRACK][10]
HarpoDccMap * GetDccMap(Long_t plane=XDCC)
void FindTrack(TClonesArray *clArray, Int_t icl0, Int_t icl1, TMatrixD Corig, Int_t plane, Int_t &color)
TClonesArray * GetClustersArray()
void SetTrackType(Int_t val)
Bool_t CheckIdClusterTrack(Int_t val)