9#include <trg/ecl/modules/MCMatcherTRGECL/MCMatcherTRGECLModule.h>
12#include <framework/gearbox/Unit.h>
15#include <ecl/dataobjects/ECLHit.h>
18#include "trg/ecl/dataobjects/TRGECLHit.h"
19#include "trg/ecl/dataobjects/TRGECLDigi0.h"
22#include <mdst/dataobjects/MCParticle.h>
23#include <framework/datastore/RelationArray.h>
56 for (
int ii = 0; ii < 100; ii++) {
68 for (
int icont = 0 ; icont < 3 ; icont ++) {
145 eclPrimaryMap.
clear();
148 for (
int iPart = 0; iPart < nParticles ; ++iPart) {
149 if (mcParticles[iPart]->getMother() == NULL) {
162 && adhoc_StableInGeneratorFlag) {
163 if (mcParticles[iPart]->getArrayIndex() == -1) {
164 eclPrimaryMap.insert(pair<int, int>(iPart, iPart));
166 eclPrimaryMap.insert(pair<int, int>(mcParticles[iPart]->getArrayIndex(), mcParticles[iPart]->getArrayIndex()));
169 if (mcParticles[iPart]->getMother() == NULL)
continue;
170 if (eclPrimaryMap.find(mcParticles[iPart]->getMother()->getArrayIndex()) != eclPrimaryMap.end()) {
171 eclPrimaryMap.insert(
172 pair<int, int>(mcParticles[iPart]->getArrayIndex(), eclPrimaryMap[mcParticles[iPart]->getMother()->getArrayIndex()]));
183 RelationArray trgeclDigi0ToMCPart(trgeclDigi0Array, mcParticles);
188 const int NofTCDigiHit = trgeclDigi0Array.
getEntries();
192 for (
int ii = 0; ii < NofTCDigiHit; ii++) {
199 int itimeindex = (int)(
TCRawTiming[ihit] / 100 + 40);
202 for (
int hit = 0; hit < nHits_hit; hit++) {
204 ECLHit* aECLHit = eclHitArray[hit];;
207 if (hitE < 0.1) {
continue;}
208 int hitCellId = aECLHit->
getCellId() - 1;
210 int timeindex = (int)((aECLHit ->getTimeAve()) / 100 + 40);
213 if (hitTCId !=
TCId[ihit]) {
continue;}
214 if (itimeindex != timeindex) {
continue;}
224 XtalId[ihit][0] = hitCellId ;
234 XtalId[ihit][1] = hitCellId ;
243 XtalId[ihit][2] = hitCellId ;
249 for (
int index = 0; index < eclHitRel.
getEntries(); index++) {
250 int PrimaryIndex = -1;
252 map<int, int>::iterator iter = eclPrimaryMap.find(eclHitRel[index].getFromIndex());
254 if (iter != eclPrimaryMap.end()) {
255 PrimaryIndex = iter->second;
257 int eclhitRelSize = eclHitRel[index].getToIndices().size();
258 for (
int pri_hit = 0; pri_hit < eclhitRelSize ; pri_hit++) {
259 int ieclHitRel = eclHitRel[index].getToIndex(pri_hit);
260 if (
ieclhit[ihit][0] == ieclHitRel) {
264 if (
ieclhit[ihit][1] == ieclHitRel) {
268 if (
ieclhit[ihit][2] == ieclHitRel) {
295 mclist = mcParticles[
TCPrimaryIndex[ihit][0]]->getMother()-> getIndex();
297 if (mclist != 1 &&
mother[ihit][0] != 0 && (mcParticles[
TCPrimaryIndex[ihit][0]]->getMother()->getMother())) {
299 mclist = mcParticles[
TCPrimaryIndex[ihit][0]]->getMother()->getMother()-> getIndex();
302 if (mclist != 1 &&
gmother[ihit][0] != 0 && (mcParticles[
TCPrimaryIndex[ihit][0]]->getMother()->getMother()->getMother())) {
303 ggmother[ihit][0] = mcParticles[
TCPrimaryIndex[ihit][0]]->getMother()->getMother()->getMother() ->getPDG();
304 mclist = mcParticles[
TCPrimaryIndex[ihit][0]]->getMother()->getMother()->getMother()-> getIndex();
308 if (mclist != 1 &&
ggmother[ihit][0] != 0) {
309 if (mcParticles[
TCPrimaryIndex[ihit][0]]->getMother()->getMother()->getMother()->getMother()) {
310 gggmother[ihit][0] = mcParticles[
TCPrimaryIndex[ihit][0]]->getMother()->getMother()->getMother()->getMother() ->getPDG();
326 mclist = mcParticles[
TCPrimaryIndex[ihit][1]]->getMother()-> getIndex();
328 if (mclist != 1 &&
mother[ihit][1] != 0 && (mcParticles[
TCPrimaryIndex[ihit][1]]->getMother()->getMother())) {
330 mclist = mcParticles[
TCPrimaryIndex[ihit][1]]->getMother()->getMother()-> getIndex();
333 if (mclist != 1 &&
gmother[ihit][1] != 0 && (mcParticles[
TCPrimaryIndex[ihit][1]]->getMother()->getMother()->getMother())) {
334 ggmother[ihit][1] = mcParticles[
TCPrimaryIndex[ihit][1]]->getMother()->getMother()->getMother() ->getPDG();
335 mclist = mcParticles[
TCPrimaryIndex[ihit][1]]->getMother()->getMother()->getMother()-> getIndex();
337 if (mclist != 1 &&
ggmother[ihit][1] != 0) {
338 if (mcParticles[
TCPrimaryIndex[ihit][1]]->getMother()->getMother()->getMother()->getMother()) {
339 gggmother[ihit][1] = mcParticles[
TCPrimaryIndex[ihit][1]]->getMother()->getMother()->getMother()->getMother() ->getPDG();
353 mclist = mcParticles[
TCPrimaryIndex[ihit][2]]->getMother()-> getIndex();
355 if (mclist != 1 &&
mother[ihit][2] != 0 && (mcParticles[
TCPrimaryIndex[ihit][2]]->getMother()->getMother())) {
357 mclist = mcParticles[
TCPrimaryIndex[ihit][2]]->getMother()->getMother()-> getIndex();
360 if (mclist != 1 &&
gmother[ihit][2] != 0 && (mcParticles[
TCPrimaryIndex[ihit][2]]->getMother()->getMother()->getMother())) {
361 ggmother[ihit][2] = mcParticles[
TCPrimaryIndex[ihit][2]]->getMother()->getMother()->getMother() ->getPDG();
362 mclist = mcParticles[
TCPrimaryIndex[ihit][2]]->getMother()->getMother()->getMother()-> getIndex();
364 if (mclist != 1 &&
ggmother[ihit][2] != 0) {
365 if (mcParticles[
TCPrimaryIndex[ihit][2]]->getMother()->getMother()->getMother()->getMother()) {
366 gggmother[ihit][2] = mcParticles[
TCPrimaryIndex[ihit][2]]->getMother()->getMother()->getMother()->getMother() ->getPDG();
377 for (
int ii = 0; ii < ihit; ++ii) {
383 TCDigiArray[m_hitNum]->setEventId(
m_nEvent);
384 TCDigiArray[m_hitNum]->setTCId(
TCId[ii]);
386 TCDigiArray[m_hitNum]->setRawEnergy(
TCRawEnergy[ii]);
387 TCDigiArray[m_hitNum]->setRawTiming(
TCRawTiming[ii]);
388 TCDigiArray[m_hitNum]->setTrackId(
trackId[ii]);
389 TCDigiArray[m_hitNum]->setCellId(
XtalId[ii]);
391 TCDigiArray[m_hitNum]->setPDG(
pid[ii]);
392 TCDigiArray[m_hitNum]->setMother(
mother[ii]);
393 TCDigiArray[m_hitNum]->setGMother(
gmother[ii]);
394 TCDigiArray[m_hitNum]->setGGMother(
ggmother[ii]);
395 TCDigiArray[m_hitNum]->setGGGMother(
gggmother[ii]);
397 TCDigiArray[m_hitNum]->setPX(
px[ii]);
398 TCDigiArray[m_hitNum]->setPY(
py[ii]);
399 TCDigiArray[m_hitNum]->setPZ(
pz[ii]);
400 TCDigiArray[m_hitNum]->setMCEnergy(
MCEnergy[ii]);
401 TCDigiArray[m_hitNum]->setContribution(
contribution[ii]);
410 RelationArray trgeclHitToMCPart(trgeclHitArray, mcParticles);
411 const int NofTCHit = trgeclHitArray.
getEntries();
414 for (
int ii = 0; ii < NofTCHit; ++ii) {
416 TRGECLHit* aTRGECLHit = trgeclHitArray[ii];
420 int itimeindex = (int)(
TCHitTiming[ii] / 100 + 40);
422 for (
int index = 0; index < trgeclDigi0ToMCPart.
getEntries(); ++index) {
424 int idigitimeindex = (int)(
TCRawTiming[idigi] / 100 + 40);
426 if (itimeindex != idigitimeindex) {
continue;}
487 for (
int ii = 0; ii < trgeclHitArray.
getEntries(); ++ii) {
491 TCHitArray[m_hitNum]->setEventId(
m_nEvent);
492 TCHitArray[m_hitNum]-> setTCId(
TCIdHit[ii]);
493 TCHitArray[m_hitNum]->setCellId(
XtalIdHit[ii]);
494 TCHitArray[m_hitNum]->setEnergyDep(
TCHitEnergy[ii]);
495 TCHitArray[m_hitNum]-> setTimeAve(
TCHitTiming[ii]);
496 TCHitArray[m_hitNum]-> setTrackId(
trackIdHit[ii]);
497 TCHitArray[m_hitNum]-> setPDG(
pidHit[ii]);
498 TCHitArray[m_hitNum]->setMother(
motherHit[ii]);
499 TCHitArray[m_hitNum]->setGMother(
gmotherHit[ii]);
500 TCHitArray[m_hitNum]->setGGMother(
ggmotherHit[ii]);
502 TCHitArray[m_hitNum]->setPX(
pxHit[ii]);
503 TCHitArray[m_hitNum]->setPY(
pyHit[ii]);
504 TCHitArray[m_hitNum]->setPZ(
pzHit[ii]);
505 TCHitArray[m_hitNum]->setMCEnergy(
MCEnergyHit[ii]);
Class to store simulated hits which equate to average of ECLSImHit on crystals input for digitization...
int getCellId() const
Get Cell ID.
double getEnergyDep() const
Get deposit energy.
double pz[100][3]
Momentum Z of particle.
double BKGContribution[100]
Backgroun Contribution in a TC
int background_tag[100][3]
Beam background tag.
int gmotherHit[100][3]
Grand mother ID
double contributionHit[100][3]
particles contribution
int pidHit[100][3]
Particle ID.
virtual void initialize() override
Initialize variables, print info, and start CPU clock.
StoreArray< TRGECLDigi0MC > m_trgECLDigi0MC
output for TRGECLDigi0MC
double pxHit[100][3]
Momentum X of particle.
int ggmother[100][3]
Grand Grand Mother ID
double MCEnergy[100][3]
Raw Energy of particle
virtual void event() override
Actual digitization of all hits in the ECL.
double SignalContribution[100]
Signal Contibution in a TC.
double contribution[100][3]
particles contribution
double py[100][3]
Momentum Y of particle.
MCMatcherTRGECLModule()
Constructor.
virtual void endRun() override
Nothing so far.
double BKGContributionHit[100]
Backgroun Contribution in a TC
int gggmotherHit[100][3]
Grand Grand Grand Mother ID
virtual void terminate() override
Stopping of CPU clock.
int gmother[100][3]
Grand mother ID
double TCHitEnergy[100]
TC Hit energy
int TCPrimaryIndexHit[100][3]
Primary Index in TC hit
double pzHit[100][3]
Momentum X of particle.
double MCEnergyHit[100][3]
Raw Energy of particle
double SignalContributionHit[100]
Signal Contibution in a TC.
int pid[100][3]
Particle ID.
int ggmotherHit[100][3]
Grand Grand Mother ID
int TCPrimaryIndex[100][3]
Primary Index in TC hit
int mother[100][3]
Mother ID
int XtalId[100][3]
XtalId in TC
StoreArray< TRGECLHitMC > m_trgECLHitMC
output for TRGECLHitMC
int m_nEvent
Event number.
virtual void beginRun() override
Nothing so far.
int gggmother[100][3]
Grand Grand Grand Mother ID
double maxEnergy[100][3]
Energy of maximum contribtion particle
double TCRawEnergy[100]
TC raw energy.
int motherHit[100][3]
Mother ID
virtual ~MCMatcherTRGECLModule()
Destructor.
int trackId[100][3]
Track Id.
int background_tagHit[100][3]
Beam background tag.
int ieclhit[100][3]
eclhit id
double pyHit[100][3]
Momentum X of particle.
double TCHitTiming[100]
TC Hit Timking
std::map< int, int > PrimaryTrackMap
define a map for Primary Track
int XtalIdHit[100][3]
XtalId in TC
double TCRawTiming[100]
TC raw timing.
int trackIdHit[100][3]
Track Id.
double px[100][3]
Momentum X of particle.
TrgEclMapping * _TCMap
object of TC Mapping
@ c_PrimaryParticle
bit 0: Particle is primary particle.
@ c_StableInGenerator
bit 1: Particle is stable, i.e., not decaying in the generator.
void setDescription(const std::string &description)
Sets the description of the module.
void setPropertyFlags(unsigned int propertyFlags)
Sets the flags for the module properties.
@ c_ParallelProcessingCertified
This module can be run in parallel processing mode safely (All I/O must be done through the data stor...
Low-level class to create/modify relations between StoreArrays.
int getEntries() const
Get the number of elements.
void add(index_type from, index_type to, weight_type weight=1.0)
Add a new element to the relation.
virtual unsigned short getBackgroundTag() const
Get background tag.
Accessor to arrays stored in the data store.
T * appendNew()
Construct a new T object at the end of the array.
int getEntries() const
Get the number of objects in the array.
void clear() override
Delete all entries in this array.
bool registerRelationTo(const StoreArray< TO > &toArray, DataStore::EDurability durability=DataStore::c_Event, DataStore::EStoreFlags storeFlags=DataStore::c_WriteOut, const std::string &namedRelation="") const
Register a relation to the given StoreArray.
Raw TC result nefor digitizing.
double getRawTiming() const
Get raw TC timing.
int getTCId() const
Get TC id.
double getRawEnergy() const
Get Energy and Timing Get raw TC energy.
double getTimeAve() const
The method to get hit average time.
int getTCId() const
The method to get TC id.
int getTCIdFromXtalId(int)
get [TC ID] from [Xtal ID]
static const double us
[microsecond]
static const double GeV
Standard of [energy, momentum, mass].
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Abstract base class for different kinds of events.