11 #include "../../littrack/std/utils/CbmLitMemoryManagment.h"
67 , fUsedHitsAllCut(0.
f)
164 unsigned int size =
fData.size();
165 for (
unsigned int iHit = 0; iHit < size; iHit++) {
166 if (
fData[iHit].fIsUsed ==
true)
continue;
173 fData[iHit].fHit.fX,
fData[iHit].fHit.fY, &indmin, &indmax);
184 vector<CbmRichHoughHit>::iterator itmin, itmax;
196 *indmin = itmin -
fData.begin();
197 *indmax = itmax -
fData.begin();
199 int arSize = *indmax - *indmin + 1;
207 for (
int i = *indmin;
i <= *indmax;
i++) {
208 if (
fData[
i].fIsUsed ==
true)
continue;
209 float ry = y0 -
fData[
i].fHit.fY;
211 float rx = x0 -
fData[
i].fHit.fX;
212 float d = rx * rx + ry * ry;
219 for (
unsigned short j = 0; j <
fNofBinsXY; j++) {
223 for (
unsigned short k = 0; k <
fNofBinsR; k++) {
229 unsigned int indmax) {
230 for (
int iPart = 0; iPart <
fNofParts; iPart++) {
238 unsigned int nofHits =
fHitInd[iPart].size();
240 float dx = 1.0f /
fDx, dy = 1.0f /
fDy, dr = 1.0f /
fDr;
242 vector<CbmRichHoughHit> data;
243 data.resize(nofHits);
244 for (
unsigned int i = 0;
i < nofHits;
i++) {
248 typedef vector<CbmRichHoughHit>::iterator iH;
250 for (iH iHit1 = data.begin(); iHit1 != data.end(); iHit1++) {
251 float iH1X = iHit1->fHit.fX;
252 float iH1Y = iHit1->fHit.fY;
253 double time1 = iHit1->fTime;
255 for (iH iHit2 = iHit1 + 1; iHit2 != data.end(); iHit2++) {
256 float iH2X = iHit2->fHit.fX;
257 float iH2Y = iHit2->fHit.fY;
258 double time2 = iHit2->fTime;
261 float rx0 = iH1X - iH2X;
262 float ry0 = iH1Y - iH2Y;
263 float r12 = rx0 * rx0 + ry0 * ry0;
267 float t10 = iHit1->fX2plusY2 - iHit2->fX2plusY2;
268 for (iH iHit3 = iHit2 + 1; iHit3 != data.end(); iHit3++) {
269 float iH3X = iHit3->fHit.fX;
270 float iH3Y = iHit3->fHit.fY;
271 double time3 = iHit3->fTime;
276 float rx1 = iH1X - iH3X;
277 float ry1 = iH1Y - iH3Y;
278 float r13 = rx1 * rx1 + ry1 * ry1;
281 float rx2 = iH2X - iH3X;
282 float ry2 = iH2Y - iH3Y;
283 float r23 = rx2 * rx2 + ry2 * ry2;
286 float det = rx2 * ry0 - rx0 * ry2;
287 if (det == 0.0
f)
continue;
288 float t19 = 0.5f / det;
289 float t5 = iHit2->fX2plusY2 - iHit3->fX2plusY2;
291 float xc = (t5 * ry0 - t10 * ry2) * t19;
293 int intX = int(xcs * dx);
294 if (intX < 0 || intX >=
fNofBinsX)
continue;
296 float yc = (t10 * rx2 - t5 * rx0) * t19;
298 int intY = int(ycs * dy);
299 if (intY < 0 || intY >=
fNofBinsY)
continue;
302 float t6 = iH1X - xc;
303 float t7 = iH1Y - yc;
305 float r =
sqrt(t6 * t6 + t7 * t7);
307 int intR = int(r * dr);
308 if (intR < 0 || intR >=
fNofBinsR)
continue;
320 int maxBinR = -1, maxR = -1;
321 unsigned int size =
fHistR.size();
322 for (
unsigned int k = 0; k < size; k++) {
323 if (
fHistR[k] > maxBinR) {
332 int maxBinXY = -1, maxXY = -1;
334 for (
unsigned int k = 0; k < size; k++) {
335 if (
fHist[k] > maxBinXY) {
341 if (maxBinXY <
fHTCut)
return;
350 r = (maxR + 0.5f) *
fDr;
351 for (
int j = indmin; j < indmax + 1; j++) {
353 rx =
fData[j].fHit.fX - xc;
354 ry =
fData[j].fHit.fY - yc;
356 dr =
fabs(
sqrt(rx * rx + ry * ry) - r);
357 if (dr > 0.9
f)
continue;
377 for (
int j = indmin; j < indmax + 1; j++) {
380 rx =
fData[j].fHit.fX - xc;
381 ry =
fData[j].fHit.fY - yc;
384 dr =
fabs(
sqrt(rx * rx + ry * ry) - r);
385 if (dr > drCOPCut)
continue;
418 for (
int j = indmin; j < indmax + 1; j++) {
424 if (dr < dCut) {
fData[j].fIsUsed =
true; }
431 set<unsigned int> usedHitsAll;
432 vector<unsigned int> goodRingIndex;
433 goodRingIndex.reserve(nofRings);
436 for (
int iRing = 0; iRing < nofRings; iRing++) {
440 bool isGoodRingAll =
true;
441 int nofUsedHitsAll = 0;
442 for (
int iHit = 0; iHit < nofHits; iHit++) {
443 set<unsigned int>::iterator it = usedHitsAll.find(ring->
GetHitId(iHit));
444 if (it != usedHitsAll.end()) { nofUsedHitsAll++; }
447 isGoodRingAll =
false;
452 for (
unsigned int iRSet = 0; iRSet < goodRingIndex.size(); iRSet++) {
455 goodRingIndex.push_back(iRing);
456 for (
int iHit = 0; iHit < nofHits; iHit++) {
457 usedHitsAll.insert(ring->
GetHitId(iHit));
464 goodRingIndex.clear();
474 for (
int iHit1 = 0; iHit1 < nofHits1; iHit1++) {
475 unsigned int hitId1 = ring1->
GetHitId(iHit1);
476 for (
int iHit2 = 0; iHit2 < nofHits2; iHit2++) {
477 unsigned int hitId2 = ring2->
GetHitId(iHit2);
478 if (hitId1 != hitId2)
continue;
480 float hitX =
fData[hitIndData].fHit.fX;
481 float hitY =
fData[hitIndData].fHit.fY;
500 unsigned int size =
fData.size();
501 for (
unsigned int i = 0;
i < size;
i++) {
502 if (
fData[
i].fHit.fId == hitId)
return i;
512 float t1, t2, t3, t4, t5, t6, t8, t9, t10, t11, t14, t16, t19, t21, t41;
518 t5 = t1 - t2 + t3 - t4;
522 t10 = t8 - t1 + t9 - t3;
526 t19 = 1.0f / (t14 * t6 - t16 * t11);
528 *xc = 0.5e0 * (t5 * t6 - t10 * t11) * t19;
529 *yc = 0.5e0 * (t10 * t14 - t5 * t16) * t19;
531 t21 = (
x[0] - *xc) * (
x[0] - *xc);
532 t41 = (
y[0] - *yc) * (
y[0] - *yc);
533 *r =
sqrt(t21 + t41);