3 #include "StEEmcPool/StEEmcA2EMaker/StEEmcA2EMaker.h"
4 #include "StEEmcUtil/EEmcGeom/EEmcGeomDefs.h"
5 #include "StEEmcUtil/EEmcGeom/EEmcGeomSimple.h"
7 #include "StEventTypes.h"
9 #include "StMuDSTMaker/COMMON/StMuDstMaker.h"
10 #include "StMuDSTMaker/COMMON/StMuDst.h"
11 #include "StMuDSTMaker/COMMON/StMuPrimaryVertex.h"
12 #include "StMuDSTMaker/COMMON/StMuTrack.h"
14 #include "StEmcUtil/geometry/StEmcGeom.h"
15 #include "StEmcUtil/others/emcDetectorName.h"
16 #include "StEmcADCtoEMaker/StBemcData.h"
17 #include "StEmcADCtoEMaker/StEmcADCtoEMaker.h"
18 #include "StEmcRawMaker/defines.h"
19 #include "StEmcRawMaker/StBemcRaw.h"
20 #include "StEmcRawMaker/StBemcTables.h"
21 #include "StEmcRawMaker/StEmcRawMaker.h"
23 #include "StGammaRawMaker.h"
30 StGammaRawMaker::StGammaRawMaker(
const char *name):
StMaker(name)
46 StGammaRawMaker::~StGammaRawMaker()
52 Int_t StGammaRawMaker::Init()
59 return StMaker::Init();
78 for (Int_t index = 0; index < kEEmcNumSectors * kEEmcNumSubSectors * kEEmcNumEtas; index++)
81 for(Int_t layer = 0; layer < 4; layer++)
84 mEEtowers[index][layer] = 0;
86 for(Int_t sec = 0; sec < kEEmcNumSectors; sec++)
89 for(Int_t plane = 0; plane < kEEmcNumSmdUVs; plane++)
92 for(Int_t strip = 0; strip < kEEmcNumStrips; strip++)
95 mEEstrips[sec][plane][strip] = 0;
108 memset(mBarrelEmcTower, 0,
sizeof(mBarrelEmcTower));
109 memset(mBarrelEmcPreshower, 0,
sizeof(mBarrelEmcPreshower));
112 mBarrelSmdEtaStrip.clear();
113 mBarrelSmdPhiStrip.clear();
123 mMuDstMaker =
dynamic_cast<StMuDstMaker*
>(GetMakerInheritsFrom(
"StMuDstMaker"));
126 LOG_WARN <<
"Make - No MuDst found!" << endm;
133 if(!GetMaker(
"mGammaEventMaker"))
135 LOG_WARN <<
"StGammaEventMaker not in chain, no tree for you" << endm;
141 mGammaEvent = mGammaEventMaker->event();
145 if(!GetDataSet(
"MuDst"))
147 LOG_WARN <<
"No MuDst" << endm;
152 if(!StMuDst::numberOfPrimaryVertices() )
return kStOK;
155 if(mUseBemc) GetBarrel();
156 if(mUseEemc) GetEndcap();
165 void StGammaRawMaker::GetTracks()
169 if(!GetDataSet(
"MuDst"))
171 LOG_WARN <<
"No MuDst" << endm;
183 if(
track->momentum().perp() < mTrackCutoff)
continue;
188 mTracks.push_back(*mGammaTrack);
199 void StGammaRawMaker::GetBarrel()
207 if(!GetDataSet(
"MuDst"))
209 LOG_WARN <<
"No MuDst" << endm;
240 for(
int crate = 1; crate <= MAXCRATES; crate++)
242 StEmcCrateStatus crateStatus = btow->crateStatus(crate);
243 if(crateStatus == crateHeaderCorrupt)
274 for(
int m = 1; m <= 120; m++)
281 StEmcGeom *geom = StEmcGeom::instance(
"bemc");
284 StSPtrVecEmcRawHit& rawHits = module->hits();
287 for(UInt_t k = 0; k < rawHits.size(); k++)
292 int id = tempRawHit->softId(BTOW);
295 bool excludedTower =
false;
297 for(UInt_t i = 0; i < mExcludedBemcTowers.size(); ++i)
299 if(
id == mExcludedBemcTowers.at(i)) excludedTower =
true;
304 LOG_DEBUG <<
" Excluding tower "<<
id << endm;
309 int ADC = tempRawHit->adc();
310 double energy = tempRawHit->energy();
312 geom->getXYZ(
id, x, y, z);
313 TVector3 position(x, y, z);
314 position -= TVector3(
StMuDst::event()->primaryVertexPosition().xyz());
315 bEta = position.Eta();
316 bPhi = position.Phi();
318 mTables->
getPedestal(BTOW,
id, CAP, pedestal, rms);
321 if(ADC < 7)
continue;
322 double pADC = ADC - pedestal;
323 if(pADC < 5.0 * rms)
continue;
326 energy = (1.0 + mBemcGainShift) * energy;
329 if ( energy / TMath::CosH(bEta) < mTowerCutoff )
continue;
334 btower->layer = kBEmcTower;
335 btower->stat = status;
336 btower->fail = (int)(status != 1);
337 btower->energy = energy;
341 LOG_DEBUG <<
" TowerID =" <<
id <<
", Status = " << status <<
", Energy = " << energy << endm;
342 LOG_DEBUG <<
" Eta = " << bEta <<
", Phi = " << bPhi <<
", ADC = " << ADC << endm;
345 mTowers.push_back(*btower);
346 mBarrelEmcTower[btower->id] = btower;
356 StEmcGeom* geom = StEmcGeom::instance(
"bprs");
358 StSPtrVecEmcRawHit& rawHits=module->hits();
361 for(UInt_t k = 0; k < rawHits.size(); k++)
366 int ADC = tempRawHit->adc();
367 double energy = tempRawHit->energy();
368 id = tempRawHit->softId(BPRS);
370 geom->getXYZ(
id,x,y,z);
371 TVector3 position(x, y, z);
372 position -= TVector3(
StMuDst::event()->primaryVertexPosition().xyz());
373 bEta = position.Eta();
374 bPhi = position.Phi();
375 mTables->
getStatus(BPRS,
id, bprs_status);
376 mTables->
getPedestal(BPRS,
id, CAP, pedestal, rms);
379 if(ADC < 15)
continue;
380 double pADC = ADC-pedestal;
381 if ( pADC < 5.0 * rms )
continue;
386 mPreShower->layer = kBEmcPres;
387 mPreShower->stat = bprs_status;
388 mPreShower->fail = (int)(bprs_status != 1);
389 mPreShower->energy = energy;
390 mPreShower->eta = bEta;
391 mPreShower->phi = bPhi;
393 LOG_DEBUG <<
" BPRS ID =" <<
id <<
", Status = " << bprs_status <<
", Energy = " << energy << endm;
394 LOG_DEBUG <<
" Eta = " << bEta <<
", Phi = " << bPhi <<
", ADC = " << ADC << endm;
396 mPreshower1.push_back(*mPreShower);
408 StEmcGeom* geom = StEmcGeom::instance(
"bsmde");
410 for (UInt_t i = 1; i <= smde->numberOfModules(); i++)
417 StSPtrVecEmcRawHit& hits = module->hits();
419 for (
size_t k = 0; k < hits.size(); ++k)
423 Int_t smde_id = hit->softId(BSMDE);
425 Int_t smde_status = 0;
426 mTables->
getStatus(BSMDE,smde_id,smde_status);
427 mTables->
getPedestal(BSMDE, smde_id, hit->calibrationType(), pedestal, rms);
430 geom->getEta(smde_id,eta);
433 if( (hit->adc() - pedestal) < 3.0 * rms)
continue;
437 bstrip->index = smde_id;
438 bstrip->sector = hit->module();
439 bstrip->plane = kBEmcSmdEta;
440 bstrip->stat = smde_status;
441 bstrip->fail = (int)(smde_status != 1);
442 bstrip->energy = hit->energy();
443 bstrip->adc = hit->adc();
444 bstrip->position = 230.705 * sinh(eta);
446 double offset = fabs(eta) < 0.5 ? 0.73 : 0.94;
448 bstrip->right = bstrip->position + offset;
449 bstrip->left = bstrip->position - offset;
451 LOG_DEBUG <<
" eStrip ID =" << smde_id <<
", Status = " << smde_status << endm;
452 LOG_DEBUG <<
" Energy = " << hit->energy() <<
", Module = " << hit->module() << endm;
454 mStrips.push_back(*bstrip);
455 mBarrelSmdEtaStrip[bstrip->index] = bstrip;
469 StEmcGeom* geom = StEmcGeom::instance(
"bsmdp");
471 for (UInt_t i = 1; i <= smdp->numberOfModules(); i++)
478 StSPtrVecEmcRawHit& hits = module->hits();
479 for (
size_t k = 0; k < hits.size(); ++k)
483 Int_t smdp_id = hit->softId(BSMDP);
485 Int_t smdp_status = 0;
486 mTables->
getStatus(BSMDP,smdp_id,smdp_status);
487 mTables->
getPedestal(BSMDP, smdp_id, hit->calibrationType(), pedestal, rms);
490 geom->getPhi(smdp_id,phi);
493 if( (hit->adc() - pedestal) < 3.0 * rms)
continue;
497 bstrip->index = smdp_id;
498 bstrip->sector = hit->module();
499 bstrip->plane = kBEmcSmdPhi;
500 bstrip->stat = smdp_status;
501 bstrip->fail = (int)(smdp_status != 1);
502 bstrip->energy = hit->energy();
503 bstrip->adc = hit->adc();
504 bstrip->position = phi;
506 double offset = 0.00293;
508 bstrip->left = phi - offset;
509 bstrip->right = phi + offset;
511 LOG_DEBUG <<
" pStrip ID =" << smdp_id <<
", Status = " << smdp_status << endm;
512 LOG_DEBUG <<
" Energy = " << hit->energy() <<
", Module = " << hit->module() << endm;
514 mStrips.push_back(*bstrip);
515 mBarrelSmdPhiStrip[bstrip->index] = bstrip;
533 void StGammaRawMaker::GetEndcap()
539 if(!primaryVertex)
return;
541 Float_t zvertex = primaryVertex->position().z();
547 LOG_WARN <<
"GetEndcap() - StEEmcA2EMaker not found!" << endm;
552 const Float_t depths[] =
560 Int_t enumerations[] = {kEEmcTower, kEEmcPre1, kEEmcPre2, kEEmcPost };
563 for(Int_t layer = 0; layer < 4; layer++)
572 UInt_t sec = (UInt_t)tower.
sector();
574 UInt_t eta = (UInt_t)tower.
etabin();
577 Float_t depth = depths[layer];
578 center.SetMag( depth/center.CosTheta() );
580 Float_t R = center.Perp();
581 Float_t Z = depth-zvertex;
582 Float_t TanTheta = R/Z;
583 Float_t Phi = center.Phi();
586 Momentum.SetMagThetaPhi(tower.
energy(), TMath::ATan(TanTheta), Phi );
590 if(layer == 0) etower = mGammaEvent->
newTower();
595 etower->id = tower.
index();
596 etower->layer = enumerations[ tower.
layer() ];
597 etower->stat = tower.
stat();
598 etower->fail = tower.
fail();
599 etower->energy = tower.
energy();
600 etower->eta = Momentum.Eta();
601 etower->phi = Momentum.Phi();
603 if( Accept(*etower) )
609 if( etower->energy / TMath::CosH(etower->eta) > mTowerCutoff )
611 mTowers.push_back(*etower);
612 mEEtowers[ etower->id ][ etower->layer ] = etower;
618 mPreshower1.push_back(*etower);
619 mEEtowers[ etower->id ][ etower->layer ] = etower;
623 mPreshower2.push_back(*etower);
624 mEEtowers[ etower->id ][ etower->layer ] = etower;
628 mPostshower.push_back(*etower);
629 mEEtowers[ etower->id ][ etower->layer ] = etower;
640 Int_t smdenum[] = {kEEmcSmdu, kEEmcSmdv};
642 for (Int_t sector = 0; sector < 12; sector++)
645 for(Int_t plane = 0; plane < 2; plane++)
654 gstrip->index = strip.
index();
655 gstrip->sector = strip.
sector();
656 gstrip->plane = smdenum[ strip.
plane() ];
657 gstrip->stat = strip.
stat();
658 gstrip->fail = strip.
fail();
659 gstrip->energy = strip.
energy();
661 mEEstrips[ gstrip->sector ][ gstrip->plane ][ gstrip->index ]= gstrip;
663 if( Accept(*gstrip) )
665 mStrips.push_back(*gstrip);
690 if(tower.fail)
return false;
702 if(strip.fail)
return false;
714 return (track->
flag() > 0 &&
715 track->
flag() / 100 != 7 &&
716 track->
flag() / 100 != 8 &&
717 track->
flag() / 100 != 9 &&
736 StGammaTower *StGammaRawMaker::tower(Int_t
id, Int_t layer)
739 if (layer >= kEEmcTower && layer <= kEEmcPost)
741 return mEEtowers[id][layer];
743 else if(layer == kBEmcTower)
745 return mBarrelEmcTower[id];
747 else if(layer == kBEmcPres)
749 return mBarrelEmcPreshower[id];
768 StGammaStrip *StGammaRawMaker::strip(Int_t sec, Int_t plane, Int_t index)
771 if(plane == kEEmcSmdu || plane == kEEmcSmdv)
773 return mEEstrips[sec][plane][index];
775 else if(plane == kBEmcSmdEta)
777 return mBarrelSmdEtaStrip[index];
779 else if(plane == kBEmcSmdPhi)
781 return mBarrelSmdPhiStrip[index];
796 mTables->
getStatus(BSMDE, strip->index, smdStatus);
798 strip->stat = smdStatus;
799 strip->fail = (int)(smdStatus != 1);
801 mStrips.push_back(*strip);
802 mBarrelSmdEtaStrip[strip->index] = strip;
816 mTables->
getStatus(BSMDP, strip->index, smdStatus);
818 strip->stat = smdStatus;
819 strip->fail = (int)(smdStatus != 1);
821 mStrips.push_back(*strip);
822 mBarrelSmdPhiStrip[strip->index] = strip;
TVector3 getTowerCenter(const UInt_t sec, const UInt_t sub, const UInt_t etabin) const
static StMuPrimaryVertex * primaryVertex()
return pointer to current primary vertex
EEmc ADC –> energy maker.
void stat(unsigned s)
Set a status bit for this element.
StGammaTower * newPostshower()
Add a new preshower2 element.
StEEmcTower & hittower(Int_t hit, Int_t layer)
void layer(Int_t l)
Sets the layer, [0-4]=[TPQR].
void Clear(Option_t *opts="")
User defined functions.
StEEmcStrip & hitstrip(Int_t sec, Int_t pl, Int_t hit)
StGammaTower * newPreshower1()
Add a new tower.
Int_t numberOfHitStrips(Int_t sector, Int_t plane) const
StGammaTower * newTower()
Add a new track.
void getPedestal(Int_t det, Int_t softId, Int_t cap, Float_t &ped, Float_t &rms) const
Return pedestal mean and rms.
Int_t numberOfHitTowers(Int_t layer) const
Int_t etabin() const
Returns the etabin of this tower, pre- or postshower element.
void fail(unsigned f)
Set a fail bit for this element.
StGammaTower * newPreshower2()
Add a new preshower1 (bprs) element.
StGammaTrack * newTrack(StMuTrack *mutr=0)
Returns tower pT in eta range.
UShort_t nHitsFit() const
Return total number of hits used in fit.
Int_t subsector() const
Returns subsector of this tower, pre- or postshower element.
short flag() const
Returns flag, (see StEvent manual for type information)
void loadTables(StMaker *anyMaker)
load tables.
void index(Int_t i)
Sets the index for this SMD strip, 0..287.
Base class for representing tower, preshower and postshower elements.
StGammaStrip * newStrip()
Add a new postshower element.
static TObjArray * primaryTracks()
returns pointer to a list of tracks belonging to the selected primary vertex
Int_t sector() const
Returns sector of this tower, pre- or postshower element.
UShort_t nHitsPoss() const
Return number of possible hits on track.
static StMuEvent * event()
returns pointer to current StMuEvent (class holding the event wise information, e.g. event number, run number)
static StEmcCollection * emcCollection()
returns pointer to current StEmcCollection
void getStatus(Int_t det, Int_t softId, Int_t &status, const char *option="") const
Return status.
void plane(Int_t p)
Sets the plane for this SMD strip, 0=U, 1=V.
void sector(Int_t s)
Sets the sector for this SMD strip.
void energy(Float_t e)
Set the energy (adc-ped+0.5)/gain for this element.
Base class for describing an endcap SMD strip.