00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00035
00036
00037
00038
00039
00040
00041
00042
00043
00044
00045
00046
00047
00048
00049
00050
00051
00052
00053
00054
00055
00056
00057
00058
00059
00060
00061
00062
00063
00064
00065
00066
00067
00068
00069
00070
00071
00072
00073
00074
00075
00076
00077
00078
00079
00080
00081
00082
00083
00084
00085
00086
00087
00088
00089
00090
00091
00092
00093
00094
00095
00096
00097
00098
00099
00100
00101
00102
00103
00104
00105
00106
00107
00108
00109
00110
00111
00112
00113
00114
00115
00116
00117
00118
00119
00120
00121
00122
00123
00124
00125
00126
00127
00128
00129
00130
00131
00132
00133
00134
00135
00136
00137
00138
00139
00140
00141
00142
00143
00144
00145
00146
00147
00148
00149
00150
00151
00152
00153
00154
00155
00156
00157
00158
00159
00160
00161
00162
00163
00164
00165
00166
00167
00168
00169
00170
00171
00172
00173
00174
00175
00176
00177
00178
00179
00180
00181
00182
00183
00184
00185
00186
00187
00188
00189
00190
00191
00192
00193
00194
00195
00196
00197
00198
00199
00200
00201
00202
00203
00204
00205
00206
00207
00208
00209
00210
00211
00212
00213
00214
00215
00216
00217
00218
00219
00220
00221
00222
00223
00224
00225
00226
00227
00228
00229
00230
00231
00232 #include "Stiostream.h"
00233 #include <iterator>
00234 #include <stdlib.h>
00235 #include <string>
00236 #include <vector>
00237 #include <algorithm>
00238 #include <cmath>
00239 #if !defined(ST_NO_NAMESPACES)
00240 using std::string;
00241 using std::vector;
00242 using std::cout;
00243 using std::ostream;
00244 #endif
00245 #include "StAssociationMaker.h"
00246 #include "StMcParameterDB.h"
00247 #include "StTrackPairInfo.hh"
00248
00249 #include "StGlobals.hh"
00250 #include "PhysicalConstants.h"
00251 #include "SystemOfUnits.h"
00252 #include "StThreeVectorF.hh"
00253
00254 #include "StMessMgr.h"
00255
00256 #include "St_DataSet.h"
00257 #include "St_DataSetIter.h"
00258 #include "TH2.h"
00259
00260 #include "StEventTypes.h"
00261
00262 #include "StMcEventTypes.hh"
00263
00264 #include "StMcEvent.hh"
00265
00266 #include "StMemoryInfo.hh"
00267
00268
00269
00270
00271
00272
00273
00274
00275
00276
00277
00278
00279
00280
00281
00282
00283
00284
00285
00286
00287
00288
00289
00290
00291
00292
00293
00294
00295
00296
00297
00298
00299
00300
00301
00302
00303
00304
00305
00306
00307
00308
00309
00310
00311
00312
00313
00314
00315
00316
00317
00318
00319
00320
00321
00322
00323
00324
00325
00326
00327
00328
00329
00330 class compFuncMcTpcHit{
00331 public:
00332 bool operator()(const StMcTpcHit*) const;
00333 void setReferenceZ(float z) { mRefZ = z; }
00334 float mRefZ;
00335 };
00336 bool compFuncMcTpcHit::operator()(const StMcTpcHit* h) const {
00337
00338 return (h->position().z()) > mRefZ;
00339 }
00340
00341 class compFuncMcFtpcHit{
00342 public:
00343 bool operator()(const StMcFtpcHit*) const;
00344 void setReferencePhi(float phi) { mRefPhi = phi; }
00345 float mRefPhi;
00346 };
00347 bool compFuncMcFtpcHit::operator()(const StMcFtpcHit* h) const {
00348
00349
00350 return (h->position().phi()/degree) > mRefPhi;
00351 }
00352
00353 bool compTrack::operator()(const StGlobalTrack* t1, const StGlobalTrack* t2) const {
00354 return t1 < t2;
00355 }
00356
00357 bool compMcTrack::operator()(const StMcTrack* t1, const StMcTrack* t2) const {
00358 return t1 < t2;
00359 }
00360 bool compKinkVertex::operator()(const StKinkVertex* v1, const StKinkVertex* v2) const {
00361 return v1 < v2;
00362 }
00363 bool compV0Vertex::operator()(const StV0Vertex* v1, const StV0Vertex* v2) const {
00364 return v1 < v2;
00365 }
00366 bool compXiVertex::operator()(const StXiVertex* v1, const StXiVertex* v2) const {
00367 return v1 < v2;
00368 }
00369
00370 bool compMcVertex::operator()(const StMcVertex* v1, const StMcVertex* v2) const {
00371 return v1 < v2;
00372 }
00373
00374
00375 ostream& operator<<(ostream& out,
00376 const pair< const StGlobalTrack* const, StTrackPairInfo*>& p)
00377 {
00378 out << "const StGlobalTrack at : " << (void*)p.first << endl;
00379 out << "const StMcTrack at : " << (void*)(p.second->partnerMcTrack()) << endl;
00380 out << "Common TPC Hits : " << p.second->commonTpcHits() << endl;
00381 out << "Common SVT Hits : " << p.second->commonSvtHits() << endl;
00382 out << "Common SSD Hits : " << p.second->commonSsdHits() << endl;
00383 out << "Common FTPC Hits : " << p.second->commonFtpcHits() << endl;
00384
00385 return out;
00386 }
00387 ostream& operator<<(ostream& out,
00388 const pair< const StMcTrack* const, StTrackPairInfo*>& p)
00389 {
00390 out << "const StMcTrack at : " << (void*)p.first << endl;
00391 out << "const StGlobalTrack at : " << (void*)(p.second->partnerTrack()) << endl;
00392 out << "Common TPC Hits : " << p.second->commonTpcHits() << endl;
00393 out << "Common SVT Hits : " << p.second->commonSvtHits() << endl;
00394 out << "Common SSD Hits : " << p.second->commonSsdHits() << endl;
00395 out << "Common FTPC Hits : " << p.second->commonFtpcHits() << endl;
00396
00397 return out;
00398 }
00399
00400
00401 ostream& operator<<(ostream& out, const rcTrackMapType& tm)
00402 {
00403 #if __SUNPRO_CC != 0x500
00404 copy(tm.begin(),tm.end(), ostream_iterator<rcTrackMapValType>(out,"\n"));
00405
00406 #else
00407 out << "Sorry, can't use ostream_iterator with this version of Stiostream.h !" << endl ;
00408 #endif
00409 return out;
00410 }
00411
00412 ostream& operator<<(ostream& out, const mcTrackMapType& tm)
00413 {
00414 #if __SUNPRO_CC != 0x500
00415 copy(tm.begin(),tm.end(), ostream_iterator<mcTrackMapValType>(out,"\n"));
00416 #else
00417 out << "Sorry, can't use ostream_iterator with this version of Stiostream.h !" << endl ;
00418 #endif
00419 return out;
00420 }
00421
00422
00423 ClassImp(StAssociationMaker)
00424
00425
00426
00427 StAssociationMaker::StAssociationMaker(const char *name, const char *title):StMaker(name,title)
00428 {
00429
00430
00431
00432 mRcTpcHitMap = 0;
00433 mMcTpcHitMap = 0;
00434 mRcSvtHitMap = 0;
00435 mMcSvtHitMap = 0;
00436 mRcSsdHitMap = 0;
00437 mMcSsdHitMap = 0;
00438 mRcFtpcHitMap = 0;
00439 mMcFtpcHitMap = 0;
00440 mRcTrackMap = 0;
00441 mMcTrackMap = 0;
00442 mRcKinkMap = 0;
00443 mMcKinkMap = 0;
00444 mRcV0Map = 0;
00445 mMcV0Map = 0;
00446 mRcXiMap = 0;
00447 mMcXiMap = 0;
00448
00449 mTpcLocalHitResolution = 0;
00450 mSvtHitResolution = 0;
00451 mSsdHitResolution = 0;
00452 mFtpcHitResolution = 0;
00453
00454 doPrintMemoryInfo = kFALSE;
00455 mL3TriggerOn = false;
00456 mInTrackerOn = kTRUE;
00457 mEstTracksOn = false;
00458 mDistanceAssoc = kFALSE;
00459 }
00460
00461
00462 StAssociationMaker::~StAssociationMaker()
00463 {
00464
00465 if(Debug()) gMessMgr->Info() << "Inside StAssociationMaker Destructor" << endm;
00466
00467 if (mRcTpcHitMap) {
00468 mRcTpcHitMap->clear();
00469 SafeDelete(mRcTpcHitMap);
00470 }
00471 if (mMcTpcHitMap) {
00472 mMcTpcHitMap->clear();
00473 SafeDelete(mMcTpcHitMap);
00474 }
00475 if (mRcSvtHitMap) {
00476 mRcSvtHitMap->clear();
00477 SafeDelete(mRcSvtHitMap);
00478 }
00479 if (mMcSvtHitMap) {
00480 mMcSvtHitMap->clear();
00481 SafeDelete(mMcSvtHitMap);
00482 }
00483 if (mRcSsdHitMap) {
00484 mRcSsdHitMap->clear();
00485 SafeDelete(mRcSsdHitMap);
00486 }
00487 if (mMcSsdHitMap) {
00488 mMcSsdHitMap->clear();
00489 SafeDelete(mMcSsdHitMap);
00490 }
00491 if (mRcFtpcHitMap) {
00492 mRcFtpcHitMap->clear();
00493 SafeDelete(mRcFtpcHitMap);
00494 }
00495 if (mMcFtpcHitMap) {
00496 mMcFtpcHitMap->clear();
00497 SafeDelete(mMcFtpcHitMap);
00498 }
00499
00500 if (mRcTrackMap) {
00501
00502
00503 for (rcTrackMapIter i=mRcTrackMap->begin(); i!=mRcTrackMap->end(); i++){
00504 delete (*i).second;
00505 }
00506
00507 mRcTrackMap->clear();
00508 SafeDelete(mRcTrackMap);
00509 }
00510 if (mMcTrackMap) {
00511 mMcTrackMap->clear();
00512 SafeDelete(mMcTrackMap);
00513 }
00514 if (mRcKinkMap) {
00515 mRcKinkMap->clear();
00516 SafeDelete(mRcKinkMap);
00517 }
00518 if (mMcKinkMap) {
00519 mMcKinkMap->clear();
00520 SafeDelete(mMcKinkMap);
00521 }
00522 if (mRcV0Map) {
00523 mRcV0Map->clear();
00524 SafeDelete(mRcV0Map);
00525 }
00526 if (mMcV0Map) {
00527 mMcV0Map->clear();
00528 SafeDelete(mMcV0Map);
00529 }
00530 if (mRcXiMap) {
00531 mRcXiMap->clear();
00532 SafeDelete(mRcXiMap);
00533 }
00534 if (mMcXiMap) {
00535 mMcXiMap->clear();
00536 SafeDelete(mMcXiMap);
00537 }
00538 }
00539
00540
00541
00542 void StAssociationMaker::Clear(const char* c)
00543 {
00544
00545 if (doPrintMemoryInfo)
00546 StMemoryInfo::instance()->snapshot();
00547
00548
00549 if (mRcTpcHitMap) {
00550 mRcTpcHitMap->clear();
00551 SafeDelete(mRcTpcHitMap);
00552 if (doPrintMemoryInfo)
00553 cout << "Deleted Rec. Tpc Hit Map" << endl;
00554 }
00555 if (mMcTpcHitMap) {
00556 mMcTpcHitMap->clear();
00557 SafeDelete(mMcTpcHitMap);
00558 if (doPrintMemoryInfo)
00559 cout << "Deleted M.C. Tpc Hit Map" << endl;
00560 }
00561 if (mRcSvtHitMap) {
00562 mRcSvtHitMap->clear();
00563 SafeDelete(mRcSvtHitMap);
00564 if (doPrintMemoryInfo)
00565 cout << "Deleted Rec. Svt Hit Map" << endl;
00566 }
00567 if (mMcSvtHitMap) {
00568 mMcSvtHitMap->clear();
00569 SafeDelete(mMcSvtHitMap);
00570 if (doPrintMemoryInfo)
00571 cout << "Deleted M.C. Svt Hit Map" << endl;
00572 }
00573 if (mRcSsdHitMap) {
00574 mRcSsdHitMap->clear();
00575 SafeDelete(mRcSsdHitMap);
00576 if (doPrintMemoryInfo)
00577 cout << "Deleted Rec. Ssd Hit Map" << endl;
00578 }
00579 if (mMcSsdHitMap) {
00580 mMcSsdHitMap->clear();
00581 SafeDelete(mMcSsdHitMap);
00582 if (doPrintMemoryInfo)
00583 cout << "Deleted M.C. Ssd Hit Map" << endl;
00584 }
00585 if (mRcFtpcHitMap) {
00586 mRcFtpcHitMap->clear();
00587 SafeDelete(mRcFtpcHitMap);
00588 if (doPrintMemoryInfo)
00589 cout << "Deleted Rec. Ftpc Hit Map" << endl;
00590 }
00591 if (mMcFtpcHitMap) {
00592 mMcFtpcHitMap->clear();
00593 SafeDelete(mMcFtpcHitMap);
00594 if (doPrintMemoryInfo)
00595 cout << "Deleted M.C. Ftpc Hit Map" << endl;
00596 }
00597
00598 if (mRcTrackMap) {
00599
00600
00601 for (rcTrackMapIter i=mRcTrackMap->begin(); i!=mRcTrackMap->end(); i++){
00602 SafeDelete((*i).second);
00603 }
00604
00605 mRcTrackMap->clear();
00606 SafeDelete(mRcTrackMap);
00607 if (doPrintMemoryInfo)
00608 cout << "Deleted Rec. Track Map" << endl;
00609 }
00610 if (mMcTrackMap) {
00611 mMcTrackMap->clear();
00612 SafeDelete(mMcTrackMap);
00613 if (doPrintMemoryInfo)
00614 cout << "Deleted M.C. Track Map" << endl;
00615
00616 }
00617 if (mRcKinkMap) {
00618 mRcKinkMap->clear();
00619 SafeDelete(mRcKinkMap);
00620 if (doPrintMemoryInfo)
00621 cout << "Deleted Rec. Kink Map" << endl;
00622 }
00623 if (mMcKinkMap) {
00624 mMcKinkMap->clear();
00625 SafeDelete(mMcKinkMap);
00626 if (doPrintMemoryInfo)
00627 cout << "Deleted M.C. Kink Map" << endl;
00628 }
00629 if (mRcV0Map) {
00630 mRcV0Map->clear();
00631 SafeDelete(mRcV0Map);
00632 if (doPrintMemoryInfo)
00633 cout << "Deleted Rec. V0 Map" << endl;
00634 }
00635 if (mMcV0Map) {
00636 mMcV0Map->clear();
00637 SafeDelete(mMcV0Map);
00638 if (doPrintMemoryInfo)
00639 cout << "Deleted M.C. V0 Map" << endl;
00640 }
00641 if (mRcXiMap) {
00642 mRcXiMap->clear();
00643 SafeDelete(mRcXiMap);
00644 if (doPrintMemoryInfo)
00645 cout << "Deleted Rec. Xi Map" << endl;
00646 }
00647 if (mMcXiMap) {
00648 mMcXiMap->clear();
00649 SafeDelete(mMcXiMap);
00650 if (doPrintMemoryInfo)
00651 cout << "Deleted M.C. Xi Map" << endl;
00652 }
00653
00654 if (doPrintMemoryInfo) {
00655 StMemoryInfo::instance()->snapshot();
00656 StMemoryInfo::instance()->print();
00657 }
00658 StMaker::Clear(c);
00659 }
00660
00661
00662 Int_t StAssociationMaker::Finish()
00663 {
00664 return StMaker::Finish();
00665 }
00666
00667
00668 Int_t StAssociationMaker::Init()
00669 {
00670
00671
00672
00673
00674
00675
00676
00677 mTpcLocalHitResolution = new TH2F("TpcLocalHitResolution",
00678 "Delta Z Vs Delta X for Nearby Hits",
00679 50, -0.52, 0.52,
00680 50, -0.24, 0.24);
00681 mTpcLocalHitResolution->SetXTitle("Local (Xmc - Xrec) (cm)");
00682 mTpcLocalHitResolution->SetYTitle("Zmc - Zrec (cm)");
00683
00684
00685
00686
00687 mSvtHitResolution = new TH2F("SvtHitResolution",
00688 "Delta Z Vs Delta X for Nearby Hits",
00689 50, -0.12, 0.12,
00690 50, -0.12, 0.12);
00691 mSvtHitResolution->SetXTitle("Xmc - Xrec (cm)");
00692 mSvtHitResolution->SetYTitle("Zmc - Zrec (cm)");
00693
00694
00695
00696 mSsdHitResolution = new TH2F("SsdHitResolution",
00697 "Delta Z Vs Delta X for Nearby Hits",
00698 50, -0.12, 0.12,
00699 50, -0.12, 0.12);
00700 mSsdHitResolution->SetXTitle("Xmc - Xrec (cm)");
00701 mSsdHitResolution->SetYTitle("Zmc - Zrec (cm)");
00702
00703
00704
00705
00706 mFtpcHitResolution = new TH2F("FtpcHitResolution",
00707 "Delta Z Vs Delta X for Nearby Hits",
00708 50, -0.32, 0.32,
00709 50, -8, 0.8);
00710 mFtpcHitResolution->SetXTitle("Rmc - Rrec (cm)");
00711 mFtpcHitResolution->SetYTitle("PHImc - PHIrec (deg)");
00712
00713 gMessMgr->Info() << "Cuts used in association for this run: " << endm;
00714 gMessMgr->Info() << "\n" << *(StMcParameterDB::instance()) << "\n" << endm;
00715 if (mInTrackerOn) {
00716 gMessMgr->Info() << "Using IT Tracks in Association" << endm;
00717 }
00718 if (mEstTracksOn) {
00719 gMessMgr->Info() << "Using EST Tracks in Association" << endm;
00720 }
00721 return StMaker::Init();
00722 }
00723
00724
00725
00726 Int_t StAssociationMaker::Make()
00727 {
00728 gMessMgr->Info() << "AssociationMaker::Make()" << endm;
00729 if (doPrintMemoryInfo)
00730 StMemoryInfo::instance()->snapshot();
00731
00732 if (mL3TriggerOn && mInTrackerOn) {
00733 gMessMgr->Warning() << "AssociationMaker::Make(): L3 and IT Tracks are both on!" << endm;
00734 return kStWarn;
00735 }
00736 if (mEstTracksOn && mInTrackerOn) {
00737 gMessMgr->Warning() << "AssociationMaker::Make(): EST and IT Tracks are both on!" << endm;
00738 return kStWarn;
00739 }
00740 if (Debug()) {
00741 if (mDistanceAssoc==true) {
00742 gMessMgr->Info() << "Using Distance Association" << endm;
00743 }
00744 else {
00745 gMessMgr->Info() << "Using Id Association" << endm;
00746 }
00747
00748 }
00749
00750
00751
00752 StEvent* rEvent = 0;
00753 rEvent = (StEvent*) GetInputDS("StEvent");
00754 if (!rEvent) {
00755 gMessMgr->Warning() << "No StEvent!!! " << endm;
00756 gMessMgr->Warning() << "Bailing out ..." << endm;
00757 return kStWarn;
00758 }
00759
00760
00761 StTpcHitCollection* rcTpcHitColl;
00762 StMcTpcHitCollection* mcTpcHitColl;
00763
00764 StSvtHitCollection* rcSvtHitColl;
00765 StMcSvtHitCollection* mcSvtHitColl;
00766
00767 StSsdHitCollection* rcSsdHitColl;
00768 StMcSsdHitCollection* mcSsdHitColl;
00769
00770 StFtpcHitCollection* rcFtpcHitColl;
00771 StMcFtpcHitCollection* mcFtpcHitColl;
00772
00773
00774
00775
00776
00777 StL3Trigger* rL3Event = 0;
00778 if (!mL3TriggerOn) {
00779
00780 rcTpcHitColl = rEvent->tpcHitCollection();
00781 rcSvtHitColl = rEvent->svtHitCollection();
00782 rcSsdHitColl = rEvent->ssdHitCollection();
00783 rcFtpcHitColl = rEvent->ftpcHitCollection();
00784 }
00785 else {
00786
00787 if (rEvent){
00788 rL3Event = rEvent->l3Trigger();
00789 }
00790
00791 if (rL3Event) {
00792 rcTpcHitColl = rL3Event->tpcHitCollection();
00793 }
00794 else {
00795 return kStWarn ;
00796 }
00797
00798 rcSvtHitColl = 0;
00799 rcSsdHitColl = 0;
00800 rcFtpcHitColl = 0;
00801
00802 }
00803
00804
00805 StSPtrVecTrackNode& rcTrackNodes = (!mL3TriggerOn) ? rEvent->trackNodes() : rL3Event->trackNodes();
00806
00807
00808
00809
00810 StMcEvent* mEvent = 0;
00811 mEvent = (StMcEvent*) GetDataSet("StMcEvent");
00812 if (!mEvent) {
00813 gMessMgr->Error() << "No StMcEvent!!! " << endm;
00814 gMessMgr->Error() << "Bailing out ..." << endm;
00815 return kStWarn;
00816 }
00817
00818
00819
00820
00821
00822 mcTpcHitColl = mEvent->tpcHitCollection();
00823 mcSvtHitColl = mEvent->svtHitCollection();
00824 mcSsdHitColl = mEvent->ssdHitCollection();
00825 mcFtpcHitColl = mEvent->ftpcHitCollection();
00826
00827
00828
00829
00830
00831 StMcParameterDB* parDB = StMcParameterDB::instance();
00832
00833
00834
00835
00836
00837 if (rcTpcHitColl && mcTpcHitColl) {
00838 if (Debug()) gMessMgr->Info() << "Making TPC Hit Associations..." << endm;
00839
00840 StTpcHit* rcTpcHit;
00841 StMcTpcHit* mcTpcHit;
00842
00843
00844 mRcTpcHitMap = new rcTpcHitMapType;
00845 mMcTpcHitMap = new mcTpcHitMapType;
00846
00847 int matchedR = 0;
00848
00849 float tpcHitDistance = 9999;
00850 if (Debug()) cout << "In Sector : ";
00851
00852 for (unsigned int iSector=0;
00853 iSector<rcTpcHitColl->numberOfSectors(); iSector++) {
00854
00855 if (Debug()) cout << iSector + 1 << " "; flush(cout);
00856 StTpcSectorHitCollection* tpcSectHitColl = rcTpcHitColl->sector(iSector);
00857 for (unsigned int iPadrow=0;
00858 iPadrow<tpcSectHitColl->numberOfPadrows();
00859 iPadrow++) {
00860 StTpcPadrowHitCollection* tpcPadRowHitColl = tpcSectHitColl->padrow(iPadrow);
00861
00862 int padrowMatchesId = 0;
00863 int padrowMatchesDi = 0;
00864 if (Debug()>=2 && iPadrow>=0) {
00865 cout << endl;
00866 cout << "Padrow " << iPadrow+1 << endl;
00867 cout << "Reco hit index \tX pos\tY pos\tZ pos\tmIdTruth" << endl;
00868
00869 for (unsigned int iHit=0;
00870 iHit<tpcPadRowHitColl->hits().size();
00871 iHit++) {
00872 rcTpcHit = tpcPadRowHitColl->hits()[iHit];
00873 int qatru=rcTpcHit->qaTruth(); int idtru= rcTpcHit->idTruth();
00874 cout << iHit << "\t" << rcTpcHit->position() << "\t" << idtru << "\t" << qatru << endl;
00875 }
00876 cout << "MC Hit Key\tX pos\tY pos\tZ pos\tparent key" << endl;
00877 for (StMcTpcHitIterator jHit = mcTpcHitColl->sector(iSector)->padrow(iPadrow)->hits().begin();
00878 jHit != mcTpcHitColl->sector(iSector)->padrow(iPadrow)->hits().end();
00879 jHit++){
00880 mcTpcHit = *jHit;
00881 cout << mcTpcHit->key() << "\t" << mcTpcHit->position() << "\t" << mcTpcHit->parentTrack()->key() << endl;
00882
00883 }
00884 }
00885 for (unsigned int iHit=0; iHit<tpcPadRowHitColl->hits().size(); iHit++){
00886
00887
00888 rcTpcHit = tpcPadRowHitColl->hits()[iHit];
00889
00890
00891
00892
00893
00894
00895
00896
00897
00898
00899
00900
00901
00902
00903 StMcHit* closestTpcHit = 0;
00904
00905 float xDiff, yDiff, zDiff;
00906 xDiff = yDiff = zDiff = -999;
00907 for (StMcTpcHitIterator jHit = mcTpcHitColl->sector(iSector)->padrow(iPadrow)->hits().begin();
00908 jHit != mcTpcHitColl->sector(iSector)->padrow(iPadrow)->hits().end();
00909 jHit++){
00910
00911 mcTpcHit = *jHit;
00912
00913
00914
00915
00916 if (mDistanceAssoc || !rcTpcHit->idTruth()) {
00917 xDiff = mcTpcHit->position().x()-rcTpcHit->position().x();
00918 yDiff = mcTpcHit->position().y()-rcTpcHit->position().y();
00919 zDiff = mcTpcHit->position().z()-rcTpcHit->position().z();
00920 if (!closestTpcHit) {
00921 tpcHitDistance=xDiff*xDiff+zDiff*zDiff;
00922 closestTpcHit = mcTpcHit;
00923 }
00924 if (xDiff*xDiff+zDiff*zDiff<tpcHitDistance) {
00925 tpcHitDistance = xDiff*xDiff+zDiff*zDiff;
00926 closestTpcHit = mcTpcHit;
00927 }
00928
00929 if ( fabs(xDiff)>= parDB->xCutTpc() ||
00930 fabs(yDiff)>= parDB->yCutTpc() ||
00931 fabs(zDiff)>= parDB->zCutTpc(mcTpcHit->position().z())) continue;
00932 } else
00933 if (mcTpcHit->parentTrack()->key() != rcTpcHit->idTruth()) continue;
00934
00935 mRcTpcHitMap->insert(rcTpcHitMapValType (rcTpcHit, mcTpcHit) );
00936 mMcTpcHitMap->insert(mcTpcHitMapValType (mcTpcHit, rcTpcHit) );
00937 rcTpcHit->SetBit(StMcHit::kMatched,1);
00938 mcTpcHit->SetBit(StMcHit::kMatched,1);
00939 ++matchedR;
00940 ++padrowMatchesDi;
00941 }
00942 if (Debug()==2)
00943 if (closestTpcHit)
00944 mTpcLocalHitResolution->Fill(closestTpcHit->position().x()-
00945 rcTpcHit->position().x(),
00946 closestTpcHit->position().z()-
00947 rcTpcHit->position().z() );
00948 }
00949 if (Debug()>=2 && iPadrow>=0) {
00950 cout << "Reco Hits in Padrow " << tpcPadRowHitColl->hits().size() << endl;
00951 cout << "Id Matches in Padrow " << padrowMatchesId << endl;
00952 cout << "Distance Matches in pr " << padrowMatchesDi << endl;
00953 }
00954 }
00955 }
00956
00957 if (Debug()) {
00958 cout << "\n";
00959 gMessMgr->Info() << "Number of Entries in TPC Hit Maps: " << mRcTpcHitMap->size()
00960 << " counter " << matchedR << endm;
00961 }
00962 if (doPrintMemoryInfo) {
00963 if (Debug()) gMessMgr->Info() << "End of TPC Hit Associations\n" << endm;
00964 StMemoryInfo::instance()->snapshot();
00965 StMemoryInfo::instance()->print();
00966 }
00967 }
00968
00969
00970
00971 if (rcSvtHitColl && mcSvtHitColl) {
00972
00973 if (Debug()) gMessMgr->Info() << "Making SVT Hit Associations..." << endm;
00974
00975 StSvtHit* rcSvtHit;
00976 StMcSvtHit* mcSvtHit;
00977
00978
00979 mRcSvtHitMap = new rcSvtHitMapType;
00980 mMcSvtHitMap = new mcSvtHitMapType;
00981
00982 float svtHitDistance = 9999;
00983 unsigned int nSvtHits = rcSvtHitColl->numberOfHits();
00984 if (Debug()) cout << "In Barrel : ";
00985 for (unsigned int iBarrel=0; nSvtHits &&
00986 iBarrel<rcSvtHitColl->numberOfBarrels(); iBarrel++) {
00987
00988 if (Debug()) cout << iBarrel + 1 << " "; flush(cout);
00989
00990 for (unsigned int iLadder=0;
00991 iLadder<rcSvtHitColl->barrel(iBarrel)->numberOfLadders();
00992 iLadder++) {
00993
00994 for (unsigned int iWafer=0;
00995 iWafer<rcSvtHitColl->barrel(iBarrel)->ladder(iLadder)->numberOfWafers();
00996 iWafer++) {
00997
00998
00999 for (unsigned int iHit=0;
01000 iHit<rcSvtHitColl->barrel(iBarrel)->ladder(iLadder)->wafer(iWafer)->hits().size();
01001 iHit++){
01002
01003 rcSvtHit = rcSvtHitColl->barrel(iBarrel)->ladder(iLadder)->wafer(iWafer)->hits()[iHit];
01004
01005 StMcSvtHit* closestSvtHit = 0;
01006 float newDist = 0;
01007 for (unsigned int jHit=0;
01008 jHit<mcSvtHitColl->barrel(iBarrel)->ladder(iLadder)->wafer(iWafer)->hits().size();
01009 jHit++){
01010
01011 mcSvtHit = mcSvtHitColl->barrel(iBarrel)->ladder(iLadder)->wafer(iWafer)->hits()[jHit];
01012 if (mDistanceAssoc || !rcSvtHit->idTruth()) {
01013 float xDiff = mcSvtHit->position().x() - rcSvtHit->position().x();
01014 float yDiff = mcSvtHit->position().y() - rcSvtHit->position().y();
01015 float zDiff = mcSvtHit->position().z() - rcSvtHit->position().z();
01016 if ( zDiff > parDB->zCutSvt() ) break;
01017
01018 if (jHit==0) {
01019 svtHitDistance=xDiff*xDiff+yDiff*yDiff+zDiff*zDiff;
01020 closestSvtHit = mcSvtHit;
01021 }
01022 if ( (newDist = xDiff*xDiff+yDiff*yDiff+zDiff*zDiff) < svtHitDistance) {
01023 svtHitDistance = newDist;
01024 closestSvtHit = mcSvtHit;
01025 }
01026
01027 if ( fabs(xDiff)>= parDB->xCutSvt() ||
01028 fabs(yDiff)>= parDB->yCutSvt() ||
01029 fabs(zDiff)>= parDB->zCutSvt()) continue;
01030 } else
01031 if (mcSvtHit->parentTrack()->key() != rcSvtHit->idTruth()) continue;
01032
01033 mRcSvtHitMap->insert(rcSvtHitMapValType (rcSvtHit, mcSvtHit) );
01034 mMcSvtHitMap->insert(mcSvtHitMapValType (mcSvtHit, rcSvtHit) );
01035 rcSvtHit->SetBit(StMcHit::kMatched,1);
01036 mcSvtHit->SetBit(StMcHit::kMatched,1);
01037 }
01038 if (closestSvtHit)
01039 mSvtHitResolution->Fill(closestSvtHit->position().x()-
01040 rcSvtHit->position().x(),
01041 closestSvtHit->position().z()-
01042 rcSvtHit->position().z() );
01043 }
01044 }
01045 }
01046 }
01047
01048 if (Debug()) {
01049 cout << "\n";
01050 gMessMgr->Info() << "Number of Entries in SVT Hit Maps: " << mRcSvtHitMap->size() << endm;
01051 }
01052 if (doPrintMemoryInfo) {
01053 if (Debug()) gMessMgr->Info() << "End of SVT Hit Associations\n" << endm;
01054 StMemoryInfo::instance()->snapshot();
01055 StMemoryInfo::instance()->print();
01056 }
01057 }
01058
01059
01060
01061 if (rcSsdHitColl && mcSsdHitColl) {
01062
01063 if (Debug()) gMessMgr->Info() << "Making SSD Hit Associations..." << endm;
01064
01065 StSsdHit* rcSsdHit;
01066 StMcSsdHit* mcSsdHit;
01067
01068
01069 mRcSsdHitMap = new rcSsdHitMapType;
01070 mMcSsdHitMap = new mcSsdHitMapType;
01071
01072 float ssdHitDistance = 9999;
01073 unsigned int nSsdHits = rcSsdHitColl->numberOfHits();
01074 if (Debug()) cout << "In Ladder : ";
01075 for (unsigned int iLadder=0; nSsdHits &&
01076 iLadder<rcSsdHitColl->numberOfLadders(); iLadder++) {
01077
01078 if (Debug()) cout << iLadder + 1 << " "; flush(cout);
01079
01080 for (unsigned int iWafer=0;
01081 iWafer<rcSsdHitColl->ladder(iLadder)->numberOfWafers();
01082 iWafer++) {
01083
01084
01085 for (unsigned int iHit=0;
01086 iHit<rcSsdHitColl->ladder(iLadder)->wafer(iWafer)->hits().size();
01087 iHit++){
01088
01089 rcSsdHit = rcSsdHitColl->ladder(iLadder)->wafer(iWafer)->hits()[iHit];
01090
01091 StMcSsdHit* closestSsdHit = 0;
01092 float newDist = 0;
01093 for (unsigned int jHit=0;
01094 jHit<mcSsdHitColl->ladder(iLadder)->wafer(iWafer)->hits().size();
01095 jHit++){
01096
01097 mcSsdHit = mcSsdHitColl->ladder(iLadder)->wafer(iWafer)->hits()[jHit];
01098 float xDiff = mcSsdHit->position().x() - rcSsdHit->position().x();
01099 float yDiff = mcSsdHit->position().y() - rcSsdHit->position().y();
01100 float zDiff = mcSsdHit->position().z() - rcSsdHit->position().z();
01101 if ( zDiff > parDB->zCutSsd() ) break;
01102
01103 if (jHit==0) {
01104 ssdHitDistance=xDiff*xDiff+yDiff*yDiff+zDiff*zDiff;
01105 closestSsdHit = mcSsdHit;
01106 }
01107 if ( (newDist = xDiff*xDiff+yDiff*yDiff+zDiff*zDiff) < ssdHitDistance) {
01108 ssdHitDistance = newDist;
01109 closestSsdHit = mcSsdHit;
01110 }
01111 if (mDistanceAssoc || !rcSsdHit->idTruth()) {
01112 if ( fabs(xDiff)>= parDB->xCutSsd() ||
01113 fabs(yDiff)>= parDB->yCutSsd() ||
01114 fabs(zDiff)>= parDB->zCutSsd()) continue;
01115 } else
01116 if (mcSsdHit->parentTrack()->key() != rcSsdHit->idTruth()) continue;
01117
01118 mRcSsdHitMap->insert(rcSsdHitMapValType (rcSsdHit, mcSsdHit) );
01119 mMcSsdHitMap->insert(mcSsdHitMapValType (mcSsdHit, rcSsdHit) );
01120 rcSsdHit->SetBit(StMcHit::kMatched,1);
01121 mcSsdHit->SetBit(StMcHit::kMatched,1);
01122 }
01123 if (closestSsdHit)
01124 mSsdHitResolution->Fill(closestSsdHit->position().x()-
01125 rcSsdHit->position().x(),
01126 closestSsdHit->position().z()-
01127 rcSsdHit->position().z() );
01128 }
01129 }
01130 }
01131
01132 if (Debug()) {
01133 cout << "\n";
01134 gMessMgr->Info() << "Number of Entries in SSD Hit Maps: " << mRcSsdHitMap->size() << endm;
01135 }
01136 if (doPrintMemoryInfo) {
01137 if (Debug()) gMessMgr->Info() << "End of SSD Hit Associations\n" << endm;
01138 StMemoryInfo::instance()->snapshot();
01139 StMemoryInfo::instance()->print();
01140 }
01141 }
01142
01143
01144
01145 if (rcFtpcHitColl && mcFtpcHitColl) {
01146 if (Debug()) gMessMgr->Info() << "Making FTPC Hit Associations..." << endm;
01147
01148 StFtpcHit* rcFtpcHit = 0;
01149 StMcFtpcHit* mcFtpcHit = 0;
01150
01151
01152 mRcFtpcHitMap = new rcFtpcHitMapType;
01153 mMcFtpcHitMap = new mcFtpcHitMapType;
01154
01155 float ftpcHitDistance = 0;
01156 float minHitDistance = 0;
01157 if (Debug()) cout << "In Plane : ";
01158 for (unsigned int iPlane=0;
01159 iPlane<rcFtpcHitColl->numberOfPlanes(); iPlane++) {
01160
01161 if (Debug()) cout << iPlane + 1 << " "; flush(cout);
01162
01163 for (unsigned int iSector=0;
01164 iSector<rcFtpcHitColl->plane(iPlane)->numberOfSectors();
01165 iSector++) {
01166
01167 compFuncMcFtpcHit ftpcComp;
01168
01169
01170 for (unsigned int iHit=0;
01171 iHit<rcFtpcHitColl->plane(iPlane)->sector(iSector)->hits().size();
01172 iHit++){
01173
01174
01175 rcFtpcHit = rcFtpcHitColl->plane(iPlane)->sector(iSector)->hits()[iHit];
01176
01177 ftpcComp.setReferencePhi((rcFtpcHit->position().phi()/degree) - parDB->phiCutFtpc());
01178
01179 float rDiff, phiDiff, rDiffMin, phiDiffMin;
01180 rDiff = phiDiff = rDiffMin = phiDiffMin =0;
01181
01182
01183
01184 StMcFtpcHitIterator ftpcHitSeed = mcFtpcHitColl->plane(iPlane)->hits().begin();
01185 bool isFirst = true;
01186 for (StMcFtpcHitIterator jHit = ftpcHitSeed;
01187 jHit != mcFtpcHitColl->plane(iPlane)->hits().end();
01188 jHit++){
01189
01190 mcFtpcHit = *jHit;
01191 rDiff = mcFtpcHit->position().perp() - rcFtpcHit->position().perp();
01192 phiDiff = (mcFtpcHit->position().phi() - rcFtpcHit->position().phi())/degree;
01193
01194 if ( phiDiff > parDB->phiCutFtpc() ) break;
01195
01196 ftpcHitDistance = (mcFtpcHit->position() - rcFtpcHit->position()).mag2();
01197
01198 if (isFirst) {
01199 minHitDistance=ftpcHitDistance;
01200 rDiffMin=rDiff;
01201 phiDiffMin=phiDiff;
01202 isFirst = false;
01203 }
01204 if (ftpcHitDistance < minHitDistance) {
01205 minHitDistance = ftpcHitDistance;
01206 rDiffMin=rDiff;
01207 phiDiffMin=phiDiff;
01208 }
01209
01210 if ( fabs(rDiff)< parDB->rCutFtpc() && fabs(phiDiff) < parDB->phiCutFtpc()) {
01211
01212 mRcFtpcHitMap->insert(rcFtpcHitMapValType (rcFtpcHit, mcFtpcHit) );
01213 mMcFtpcHitMap->insert(mcFtpcHitMapValType (mcFtpcHit, rcFtpcHit) );
01214 rcFtpcHit->SetBit(StMcHit::kMatched,1);
01215 rcFtpcHit->setIdTruth(mcFtpcHit->parentTrack()->key(),100);
01216 mcFtpcHit->SetBit(StMcHit::kMatched,1);
01217 }
01218
01219 }
01220 if (!isFirst)
01221 mFtpcHitResolution->Fill(rDiffMin,phiDiffMin);
01222 }
01223 }
01224 }
01225
01226 if (Debug()) {
01227 cout << "\n";
01228 gMessMgr->Info() << "Number of Entries in Ftpc Hit Maps: " << mRcFtpcHitMap->size() << endm;
01229 }
01230 if (doPrintMemoryInfo) {
01231 if (Debug()) gMessMgr->Info() << "End of FTPC Hit Associations\n" << endm;
01232 StMemoryInfo::instance()->snapshot();
01233 StMemoryInfo::instance()->print();
01234 }
01235
01236 }
01237
01238
01239
01240 bool smallTpcHitMap, smallSvtHitMap, smallSsdHitMap, smallFtpcHitMap;
01241 smallTpcHitMap = smallSvtHitMap = smallSsdHitMap = smallFtpcHitMap = false;
01242
01243 if (mRcTpcHitMap && mRcTpcHitMap->size()>0 && mRcTpcHitMap->size() < parDB->reqCommonHitsTpc()) {
01244 gMessMgr->Warning() << "\n----------- WARNING ---------------\n"
01245 << " The Tpc Hit Map is too small for \n"
01246 << " any meaningful track association. \n"
01247 << " ------------------------------------ \n"
01248 << "Entries in Hit Map : " << mRcTpcHitMap->size() << "\n"
01249 << "Required Common Hits: " << parDB->reqCommonHitsTpc() << "\n"
01250 << "Suggest increase distance cuts." << endm;
01251 smallTpcHitMap = true;
01252 }
01253
01254 if (mRcSvtHitMap && mRcSvtHitMap->size()>0 && mRcSvtHitMap->size() < parDB->reqCommonHitsSvt()) {
01255 gMessMgr->Warning() << "\n----------- WARNING ---------------\n"
01256 << " The Svt Hit Map is too small for \n"
01257 << " any meaningful track association. \n"
01258 << " ------------------------------------ \n"
01259 << "Entries in Hit Map : " << mRcSvtHitMap->size() << "\n"
01260 << "Required Common Hits: " << parDB->reqCommonHitsSvt() << "\n"
01261 << "Suggest increase distance cuts." << endm;
01262 smallSvtHitMap = true;
01263 }
01264 if (mRcSsdHitMap && mRcSsdHitMap->size()>0 && mRcSsdHitMap->size() < parDB->reqCommonHitsSsd()) {
01265 gMessMgr->Warning() << "\n----------- WARNING ---------------\n"
01266 << " The Ssd Hit Map is too small for \n"
01267 << " any meaningful track association. \n"
01268 << " ------------------------------------ \n"
01269 << "Entries in Hit Map : " << mRcSsdHitMap->size() << "\n"
01270 << "Required Common Hits: " << parDB->reqCommonHitsSsd() << "\n"
01271 << "Suggest increase distance cuts." << endm;
01272 smallSsdHitMap = true;
01273 }
01274 if (mRcFtpcHitMap && mRcFtpcHitMap->size()>0 && mRcFtpcHitMap->size() < parDB->reqCommonHitsFtpc()) {
01275 gMessMgr->Warning() << "\n----------- WARNING ---------------\n"
01276 << " The Ftpc Hit Map is too small for \n"
01277 << " any meaningful track association. \n"
01278 << " ------------------------------------ \n"
01279 << "Entries in Hit Map : " << mRcFtpcHitMap->size() << "\n"
01280 << "Required Common Hits: " << parDB->reqCommonHitsFtpc() << "\n"
01281 << "Suggest increase distance cuts." << endm;
01282 smallFtpcHitMap = true;
01283 }
01284
01285 if ((smallTpcHitMap && smallSvtHitMap && smallSsdHitMap && smallFtpcHitMap) ||
01286 (!mRcTpcHitMap && !mRcSvtHitMap && mRcSsdHitMap && !mRcFtpcHitMap)) {
01287 gMessMgr->Error() << "No Useful Hit Map to make Track Associations" << endm;
01288 return kStWarn;
01289 }
01290
01291
01292
01293
01294
01295 StTrackNode* trkNode;
01296 const StGlobalTrack* rcTrack;
01297
01298 StHit* rcHit;
01299 StTpcHit* rcKeyTpcHit;
01300 StSvtHit* rcKeySvtHit;
01301 StSsdHit* rcKeySsdHit;
01302 StFtpcHit* rcKeyFtpcHit;
01303
01304 pair<rcTpcHitMapIter,rcTpcHitMapIter> boundsTpc;
01305 pair<rcSvtHitMapIter,rcSvtHitMapIter> boundsSvt;
01306 pair<rcSsdHitMapIter,rcSsdHitMapIter> boundsSsd;
01307 pair<rcFtpcHitMapIter,rcFtpcHitMapIter> boundsFtpc;
01308
01309 rcTpcHitMapIter tpcHMIter;
01310 rcSvtHitMapIter svtHMIter;
01311 rcSsdHitMapIter ssdHMIter;
01312 rcFtpcHitMapIter ftpcHMIter;
01313
01314 const StMcTpcHit* mcValueTpcHit;
01315 const StMcSvtHit* mcValueSvtHit;
01316 const StMcSsdHit* mcValueSsdHit;
01317 const StMcFtpcHit* mcValueFtpcHit;
01318
01319 const StMcTrack* trackCand;
01320 StTrackPairInfo* trkPair;
01321
01322 trackPing initializedTrackPing;
01323 initializedTrackPing.mcTrack = 0;
01324 initializedTrackPing.nPingsTpc = 0;
01325 initializedTrackPing.nPingsSvt = 0;
01326 initializedTrackPing.nPingsSsd = 0;
01327 initializedTrackPing.nPingsFtpc = 0;
01328
01329
01330
01331 mRcTrackMap = new rcTrackMapType;
01332 mMcTrackMap = new mcTrackMapType;
01333
01334 if(Debug()) gMessMgr->Info() << "Making Track Associations..." << endm;
01335
01336
01337
01338
01339 unsigned int trkNodeI;
01340 for (trkNodeI = 0; trkNodeI < rcTrackNodes.size(); trkNodeI++){
01341
01342 #ifndef ST_NO_TEMPLATE_DEF_ARGS
01343 vector<trackPing> candidates(20, initializedTrackPing);
01344 #else
01345 vector<trackPing, allocator<trackPing> > candidates(20, initializedTrackPing);
01346 #endif
01347 trkNode = rcTrackNodes[trkNodeI];
01348 if (!mEstTracksOn)
01349 rcTrack = dynamic_cast<const StGlobalTrack*>(trkNode->track(global));
01350 else {
01351 rcTrack = dynamic_cast<const StGlobalTrack*>(trkNode->track(global));
01352 }
01353 if (!rcTrack || !(rcTrack->detectorInfo()->hits().size()))
01354 continue;
01355 if (mInTrackerOn && rcTrack->encodedMethod()!=263) continue;
01356 if (!mInTrackerOn && rcTrack->encodedMethod()==263) continue;
01357 unsigned int nCandidates = 0;
01358
01359
01360
01361
01362
01363
01364 if (mRcTpcHitMap) {
01365 StPtrVecHit recTpcHits = rcTrack->detectorInfo()->hits(kTpcId);
01366 unsigned int recTpcHitI;
01367 for (recTpcHitI = 0; recTpcHitI < recTpcHits.size(); recTpcHitI++) {
01368
01369 rcHit = recTpcHits[recTpcHitI];
01370 rcKeyTpcHit = dynamic_cast<StTpcHit*>(rcHit);
01371
01372 if (!rcKeyTpcHit) continue;
01373 boundsTpc = mRcTpcHitMap->equal_range(rcKeyTpcHit);
01374
01375 for (tpcHMIter=boundsTpc.first; tpcHMIter!=boundsTpc.second; ++tpcHMIter) {
01376
01377 mcValueTpcHit = (*tpcHMIter).second;
01378 if (! mcValueTpcHit) continue;
01379 trackCand = mcValueTpcHit->parentTrack();
01380
01381
01382
01383
01384
01385
01386 if (nCandidates == 0) {
01387 candidates[0].mcTrack = trackCand;
01388 candidates[0].nPingsTpc = 1;
01389 nCandidates++;
01390
01391 }
01392
01393 else {
01394 for (unsigned int iCandidate=0; iCandidate<nCandidates; iCandidate++){
01395 if (trackCand==candidates[iCandidate].mcTrack){
01396 candidates[iCandidate].nPingsTpc++;
01397 break;
01398 }
01399 if (iCandidate == (nCandidates-1)){
01400 candidates[nCandidates].mcTrack = trackCand;
01401 candidates[nCandidates].nPingsTpc = 1;
01402 nCandidates++;
01403
01404
01405 if (nCandidates>=candidates.size()) {
01406 candidates.resize(nCandidates+20);
01407 if (Debug()) cout << "Resizing in the TPC hits of the track " << endl;
01408 }
01409 break;
01410 }
01411 }
01412
01413 }
01414 }
01415
01416 }
01417 }
01418
01419
01420
01421 if (mRcSvtHitMap) {
01422 StPtrVecHit recSvtHits = rcTrack->detectorInfo()->hits(kSvtId);
01423 unsigned int recSvtHitI;
01424 for (recSvtHitI = 0; recSvtHitI < recSvtHits.size(); recSvtHitI++) {
01425
01426
01427 rcHit = recSvtHits[recSvtHitI];
01428 rcKeySvtHit = dynamic_cast<StSvtHit*>(rcHit);
01429
01430 if (!rcKeySvtHit) continue;
01431 boundsSvt = mRcSvtHitMap->equal_range(rcKeySvtHit);
01432
01433 for (svtHMIter=boundsSvt.first; svtHMIter!=boundsSvt.second; ++svtHMIter) {
01434
01435 mcValueSvtHit = (*svtHMIter).second;
01436 trackCand = mcValueSvtHit->parentTrack();
01437
01438
01439
01440
01441
01442
01443 if (nCandidates == 0) {
01444 candidates[0].mcTrack = trackCand;
01445 candidates[0].nPingsSvt = 1;
01446 nCandidates++;
01447
01448 }
01449
01450 else {
01451 for (unsigned int iCandidate=0; iCandidate<nCandidates; iCandidate++){
01452 if (trackCand==candidates[iCandidate].mcTrack){
01453 candidates[iCandidate].nPingsSvt++;
01454 break;
01455 }
01456 if (iCandidate == (nCandidates-1)){
01457 candidates[nCandidates].mcTrack = trackCand;
01458 candidates[nCandidates].nPingsSvt = 1;
01459 nCandidates++;
01460
01461
01462 if (nCandidates>=candidates.size()) {
01463 candidates.resize(nCandidates+20);
01464 if (Debug()) cout << "Resizing in the SVT hits of the track " << endl;
01465 }
01466 break;
01467 }
01468 }
01469
01470 }
01471 }
01472
01473 }
01474 }
01475
01476
01477
01478 if (mRcSsdHitMap) {
01479 StPtrVecHit recSsdHits = rcTrack->detectorInfo()->hits(kSsdId);
01480 unsigned int recSsdHitI;
01481 for (recSsdHitI = 0; recSsdHitI < recSsdHits.size(); recSsdHitI++) {
01482
01483
01484 rcHit = recSsdHits[recSsdHitI];
01485 rcKeySsdHit = dynamic_cast<StSsdHit*>(rcHit);
01486
01487 if (!rcKeySsdHit) continue;
01488 boundsSsd = mRcSsdHitMap->equal_range(rcKeySsdHit);
01489
01490 for (ssdHMIter=boundsSsd.first; ssdHMIter!=boundsSsd.second; ++ssdHMIter) {
01491
01492 mcValueSsdHit = (*ssdHMIter).second;
01493 trackCand = mcValueSsdHit->parentTrack();
01494
01495
01496
01497
01498
01499
01500 if (nCandidates == 0) {
01501 candidates[0].mcTrack = trackCand;
01502 candidates[0].nPingsSsd = 1;
01503 nCandidates++;
01504
01505 }
01506
01507 else {
01508 for (unsigned int iCandidate=0; iCandidate<nCandidates; iCandidate++){
01509 if (trackCand==candidates[iCandidate].mcTrack){
01510 candidates[iCandidate].nPingsSsd++;
01511 break;
01512 }
01513 if (iCandidate == (nCandidates-1)){
01514 candidates[nCandidates].mcTrack = trackCand;
01515 candidates[nCandidates].nPingsSsd = 1;
01516 nCandidates++;
01517
01518
01519 if (nCandidates>=candidates.size()) {
01520 candidates.resize(nCandidates+20);
01521 if (Debug()) cout << "Resizing in the SSD hits of the track " << endl;
01522 }
01523 break;
01524 }
01525 }
01526
01527 }
01528 }
01529
01530 }
01531 }
01532
01533
01534
01535
01536
01537
01538 if (mRcFtpcHitMap) {
01539 StPtrVecHit recFtpcHitsW = rcTrack->detectorInfo()->hits(kFtpcWestId);
01540 StPtrVecHit recFtpcHitsE = rcTrack->detectorInfo()->hits(kFtpcEastId);
01541 unsigned int recFtpcHitI;
01542
01543
01544 for (recFtpcHitI = 0; recFtpcHitI < recFtpcHitsW.size(); recFtpcHitI++) {
01545
01546 rcHit = recFtpcHitsW[recFtpcHitI];
01547 rcKeyFtpcHit = dynamic_cast<StFtpcHit*>(rcHit);
01548
01549 if (!rcKeyFtpcHit) continue;
01550 boundsFtpc = mRcFtpcHitMap->equal_range(rcKeyFtpcHit);
01551
01552 for (ftpcHMIter=boundsFtpc.first; ftpcHMIter!=boundsFtpc.second; ++ftpcHMIter) {
01553
01554 mcValueFtpcHit = (*ftpcHMIter).second;
01555 trackCand = mcValueFtpcHit->parentTrack();
01556
01557
01558
01559
01560
01561
01562 if (nCandidates == 0) {
01563 candidates[0].mcTrack = trackCand;
01564 candidates[0].nPingsFtpc = 1;
01565 nCandidates++;
01566
01567 }
01568
01569 else {
01570 for (unsigned int iCandidate=0; iCandidate<nCandidates; iCandidate++){
01571 if (trackCand==candidates[iCandidate].mcTrack){
01572 candidates[iCandidate].nPingsFtpc++;
01573 break;
01574 }
01575 if (iCandidate == (nCandidates-1)){
01576 candidates[nCandidates].mcTrack = trackCand;
01577 candidates[nCandidates].nPingsFtpc = 1;
01578 nCandidates++;
01579
01580
01581 if (nCandidates>=candidates.size()) {
01582 candidates.resize(nCandidates+20);
01583 if (Debug()) cout << "Resizing in the East FTPC hits of the track " << endl;
01584 }
01585 break;
01586 }
01587 }
01588
01589 }
01590 }
01591
01592 }
01593
01594
01595 for (recFtpcHitI = 0; recFtpcHitI < recFtpcHitsE.size(); recFtpcHitI++) {
01596
01597
01598 rcHit = recFtpcHitsE[recFtpcHitI];
01599 rcKeyFtpcHit = dynamic_cast<StFtpcHit*>(rcHit);
01600
01601 if (!rcKeyFtpcHit) continue;
01602 boundsFtpc = mRcFtpcHitMap->equal_range(rcKeyFtpcHit);
01603
01604 for (ftpcHMIter=boundsFtpc.first; ftpcHMIter!=boundsFtpc.second; ++ftpcHMIter) {
01605
01606 mcValueFtpcHit = (*ftpcHMIter).second;
01607 trackCand = mcValueFtpcHit->parentTrack();
01608
01609
01610
01611
01612
01613
01614 if (nCandidates == 0) {
01615 candidates[0].mcTrack = trackCand;
01616 candidates[0].nPingsFtpc = 1;
01617 nCandidates++;
01618
01619 }
01620
01621 else {
01622 for (unsigned int iCandidate=0; iCandidate<nCandidates; iCandidate++){
01623 if (trackCand==candidates[iCandidate].mcTrack){
01624 candidates[iCandidate].nPingsFtpc++;
01625 break;
01626 }
01627 if (iCandidate == (nCandidates-1)){
01628 candidates[nCandidates].mcTrack = trackCand;
01629 candidates[nCandidates].nPingsFtpc = 1;
01630 nCandidates++;
01631
01632
01633 if (nCandidates>=candidates.size()) {
01634 candidates.resize(nCandidates+20);
01635 if (Debug()) cout << "Resizing in the West FTPC hits of the track " << endl;
01636 }
01637 break;
01638 }
01639 }
01640
01641 }
01642 }
01643
01644 }
01645 }
01646
01647
01648
01649
01650 if (nCandidates>20 || candidates.size()>20)
01651 cout << "We Have " << candidates.size() << " candidates!!! " << endl;
01652 if (candidates.size()<nCandidates) {
01653 cout << "The candidate track vector has grown more than expected!! " << endl;
01654 cout << "Something is wrong! We probably went out-of-bounds! " << endl;
01655 }
01656 for (unsigned int iCandidate=0; iCandidate<nCandidates; iCandidate++){
01657
01658
01659
01660 if (candidates[iCandidate].nPingsTpc >= parDB->reqCommonHitsTpc() ||
01661 candidates[iCandidate].nPingsSvt >= parDB->reqCommonHitsSvt() ||
01662 candidates[iCandidate].nPingsSsd >= parDB->reqCommonHitsSsd() ||
01663 candidates[iCandidate].nPingsFtpc >= parDB->reqCommonHitsFtpc()){
01664
01665
01666
01667 trkPair = new StTrackPairInfo(rcTrack,
01668 candidates[iCandidate].mcTrack,
01669 candidates[iCandidate].nPingsTpc,
01670 candidates[iCandidate].nPingsSvt,
01671 candidates[iCandidate].nPingsSsd,
01672 candidates[iCandidate].nPingsFtpc);
01673 mRcTrackMap->insert(rcTrackMapValType (rcTrack, trkPair));
01674 mMcTrackMap->insert(mcTrackMapValType (candidates[iCandidate].mcTrack, trkPair));
01675
01676 }
01677 }
01678
01679
01680 candidates.clear();
01681
01682 }
01683
01684 if (Debug() > 1) {
01685 cout << "The RcTrack map is now" << endl << *mRcTrackMap << endl;
01686 cout << "The McTrack map is now" << endl << *mMcTrackMap << endl;
01687 }
01688
01689 if(Debug()){
01690 gMessMgr->Info() << "Number of Entries in Track Maps: " << mRcTrackMap->size() << endm;
01691 }
01692 if (doPrintMemoryInfo) {
01693 cout << "End of Track Associations\n";
01694 StMemoryInfo::instance()->snapshot();
01695 StMemoryInfo::instance()->print();
01696 }
01697
01698
01699
01700
01701 if (!mL3TriggerOn) {
01702
01703 mRcKinkMap = new rcKinkMapType;
01704 mMcKinkMap = new mcKinkMapType;
01705 mRcV0Map = new rcV0MapType;
01706 mMcV0Map = new mcV0MapType;
01707 mRcXiMap = new rcXiMapType;
01708 mMcXiMap = new mcXiMapType;
01709
01710 if(Debug()) gMessMgr->Info() << "Making Vertex Associations" << endm;
01711
01712 StSPtrVecKinkVertex& kinks = rEvent->kinkVertices();
01713
01714
01715 if(Debug()) gMessMgr->Info() << "Kinks..." << endm;
01716
01717
01718
01719 pair<rcTrackMapIter, rcTrackMapIter> kinkBoundsDaughter, kinkBoundsParent;
01720
01721 const StMcVertex* primary = mEvent->primaryVertex();
01722 for (StKinkVertexIterator kvi = kinks.begin(); kvi!=kinks.end(); kvi++) {
01723
01724 StKinkVertex* rcKink = *kvi;
01725 StTrack* kinkDaughter = rcKink->daughter(0);
01726 const StGlobalTrack* gKinkDaughter = dynamic_cast<const StGlobalTrack*>(kinkDaughter);
01727 if (!gKinkDaughter) continue;
01728
01729 StTrack* kinkParent = rcKink->parent();
01730 const StGlobalTrack* gKinkParent = dynamic_cast<const StGlobalTrack*>(kinkParent);
01731 if (!gKinkParent) continue;
01732
01733
01734 kinkBoundsDaughter = mRcTrackMap->equal_range(gKinkDaughter);
01735
01736 for (rcTrackMapIter trkIter = kinkBoundsDaughter.first; trkIter!=kinkBoundsDaughter.second; trkIter++) {
01737 const StMcTrack* mcDaughter = (*trkIter).second->partnerMcTrack();
01738
01739 const StMcVertex* mcKink = mcDaughter->startVertex();
01740 if (mcKink == primary || mcKink == 0) continue;
01741 const StMcTrack* mcParent = mcKink->parent();
01742
01743
01744 kinkBoundsParent = mRcTrackMap->equal_range(gKinkParent);
01745
01746 for (rcTrackMapIter trkIter2 = kinkBoundsParent.first;
01747 trkIter2!=kinkBoundsParent.second; trkIter2++) {
01748
01749 if (mcParent == (*trkIter2).second->partnerMcTrack() ) {
01750
01751
01752 mRcKinkMap->insert(rcKinkMapValType (rcKink, mcKink));
01753 mMcKinkMap->insert(mcKinkMapValType (mcKink, rcKink));
01754 }
01755 }
01756
01757 }
01758 }
01759 if(Debug()){
01760 gMessMgr->Info() << "Number of Entries in kink Maps: " << mRcKinkMap->size() << endm;
01761 }
01762 if (doPrintMemoryInfo) {
01763 cout << "End of kink Associations\n";
01764 StMemoryInfo::instance()->snapshot();
01765 StMemoryInfo::instance()->print();
01766 }
01767
01768 if(Debug()) gMessMgr->Info() << "V0s..." << endm;
01769
01770 StSPtrVecV0Vertex& v0s = rEvent->v0Vertices();
01771
01772
01773 for (StV0VertexIterator v0vi = v0s.begin(); v0vi!=v0s.end(); v0vi++) {
01774 StV0Vertex* rcV0 = *v0vi;
01775 StTrack* v0Daughter1 = rcV0->daughter(0);
01776 if (mEstTracksOn && rcV0->chiSquared()==-16) continue;
01777 if (!mEstTracksOn && rcV0->chiSquared()<-16) continue;
01778 const StGlobalTrack* gV0Daughter1 = dynamic_cast<const StGlobalTrack*>(v0Daughter1);
01779 if (!gV0Daughter1) continue;
01780
01781 StTrack* v0Daughter2 = rcV0->daughter(1);
01782 const StGlobalTrack* gV0Daughter2 = dynamic_cast<const StGlobalTrack*>(v0Daughter2);
01783 if (!gV0Daughter2) continue;
01784
01785 pair<rcTrackMapIter, rcTrackMapIter> v0Bounds1 = mRcTrackMap->equal_range(gV0Daughter1);
01786 pair<rcTrackMapIter, rcTrackMapIter> v0Bounds2 = mRcTrackMap->equal_range(gV0Daughter2);
01787 for (rcTrackMapIter trkIter1 = v0Bounds1.first; trkIter1!=v0Bounds1.second; trkIter1++) {
01788 const StMcTrack* mcDaughter1 = (*trkIter1).second->partnerMcTrack();
01789 for (rcTrackMapIter trkIter2 = v0Bounds2.first; trkIter2!=v0Bounds2.second; trkIter2++) {
01790 const StMcTrack* mcDaughter2 = (*trkIter2).second->partnerMcTrack();
01791 if (mcDaughter1->startVertex() == mcDaughter2->startVertex() &&
01792 mcDaughter1->startVertex() != primary &&
01793 mcDaughter1->startVertex() != 0) {
01794
01795 mRcV0Map->insert(rcV0MapValType (rcV0, mcDaughter1->startVertex()));
01796 mMcV0Map->insert(mcV0MapValType (mcDaughter1->startVertex(), rcV0));
01797
01798 }
01799 }
01800 }
01801
01802 }
01803 if(Debug()) {
01804 gMessMgr->Info() << "Number of Entries in V0 Maps: " << mRcV0Map->size() << endm;
01805 }
01806 if (doPrintMemoryInfo) {
01807 cout << "End of V0 Associations\n";
01808 StMemoryInfo::instance()->snapshot();
01809 StMemoryInfo::instance()->print();
01810 }
01811
01812 if(Debug()) gMessMgr->Info() << "Xis..." << endm;
01813
01814 StSPtrVecXiVertex& xis = rEvent->xiVertices();
01815
01816
01817 for (StXiVertexIterator xvi = xis.begin(); xvi!=xis.end(); xvi++) {
01818 StXiVertex* rcXi = *xvi;
01819 if (mEstTracksOn && rcXi->chiSquared()==-16) continue;
01820 if (!mEstTracksOn && rcXi->chiSquared()<-16) continue;
01821 StV0Vertex* rcV0ofXi = rcXi->v0Vertex();
01822 StTrack* rcBachelor = rcXi->bachelor();
01823 const StGlobalTrack* gRcBachelor = dynamic_cast<const StGlobalTrack*>(rcBachelor);
01824 if (!gRcBachelor) continue;
01825 pair<rcTrackMapIter, rcTrackMapIter> xiBounds = mRcTrackMap->equal_range(gRcBachelor);
01826 for (rcTrackMapIter trkIter3 = xiBounds.first; trkIter3!= xiBounds.second; trkIter3++){
01827 const StMcTrack* mcBachelor = (*trkIter3).second->partnerMcTrack();
01828 const StMcVertex* mcXi = mcBachelor->startVertex();
01829 if (mcXi == primary || mcXi == 0) continue;
01830 pair<rcV0MapIter, rcV0MapIter> xiBoundsV0 = mRcV0Map->equal_range(rcV0ofXi);
01831 for (rcV0MapIter v0Iter = xiBoundsV0.first; v0Iter!= xiBoundsV0.second; v0Iter++){
01832 const StMcVertex* mcV0 = (*v0Iter).second;
01833 if (mcV0->parent() != 0 && mcXi == mcV0->parent()->startVertex()) {
01834
01835 mRcXiMap->insert(rcXiMapValType (rcXi, mcXi));
01836 mMcXiMap->insert(mcXiMapValType (mcXi, rcXi));
01837
01838 }
01839 }
01840 }
01841 }
01842 if(Debug()) {
01843 gMessMgr->Info() << "Number of Entries in Xi Maps: " << mRcXiMap->size() << endm;
01844 }
01845
01846 }
01847 if (doPrintMemoryInfo) {
01848 cout << "End of Make()\n";
01849 StMemoryInfo::instance()->snapshot();
01850 StMemoryInfo::instance()->print();
01851 }
01852 return kStOK;
01853 }