Linear Collider Forum



Home » Analysis and Reconstruction » Reconstruction » MC association
MC association [message #475] Fri, 19 May 2006 02:23 Go to next message
fabio
Messages: 3
Registered: April 2006
Hello,

I´m trying to get a MC association in order to associate reco photons to Monte Carlo ones.
My idea is to start from reco gamma, trying to get the corresponding reco clusters and so the reco calorimeter hits.
At this point it is possible to utilize the function "getRawHit()" of "CalorimeterHit" to get the Sim calorimeter hits and so all infos about the corresponding MC particle that caused these hits.

In particolar such a process is simply:
************************************************************ **
// loop over all RC gamma of the event
for(int i=0; i< num_RC_gamma ; i++){ ReconstructedParticleImpl* rp =
dynamic_cast<ReconstructedParticleImpl*>( RCcol->getElementAt( i ) ) ;
RecoClusters = rp->getClusters();
for( std::vector<Cluster*>::const_iterator cl =
RecoClusters.begin(); cl != RecoClusters.end(); cl++){
RecoCaloHits = (*cl)->getCalorimeterHits();

// loop over all Reco Cluster Hits
for( std::vector<CalorimeterHit*>::const_iterator ch = RecoCaloHits.begin(); ch != RecoCaloHits.end(); ch++){
LCObject* CSimhit = (*ch)->getRawHit();
SimCalorimeterHit* CRawHit = dynamic_cast<SimCalorimeterHit*> (CSimhit) ;

} // end for const_iterator ch
} // end for const_iterator cl
} // end for num_RC_gamma

************************************************************ ****
this code has been compliled without any problems... but when I try to execute it (qqbar collections) I have a segmentation violation..
In particular, using some printouts, CRawHit returns a pointer to 0...

1) Is it possible to arrange a MC association in this way? Has the calorimeter hit an associated RawCalorimeter hit?

2) how to associate CalorimeterHits to MCParticles?

Thanks in advance,
Fabio.

[Updated on: Fri, 19 May 2006 02:25]

Re: MC association [message #476 is a reply to message #475] Fri, 19 May 2006 06:16 Go to previous messageGo to next message
lima
Messages: 47
Registered: May 2004
Location: DeKalb, IL, USA
Hi Fabio,

At this point it is unlikely that you are using DigiSim upstream of your analysis code. Without DigiSim, there are no RawCalHits. If this is the case, all you need to do is to cast the CalorimeterHit* to a SimCalorimeterHit* and take it from there:

// loop over all Reco Cluster Hits
for( std::vector<CalorimeterHit*>::const_iterator ch = RecoCaloHits.begin(); ch != RecoCaloHits.end(); ch++){
// LCObject* CSimhit = (*ch)->getRawHit();
SimCalorimeterHit* simHit = dynamic_cast<SimCalorimeterHit*>(*ch) ;
int nMCcont = simHit->getNMCParticles();
//... and so on...
}

Please let us know if this does not work.

Good luck,
Guilherme
Re: MC association [message #477 is a reply to message #476] Fri, 19 May 2006 07:07 Go to previous messageGo to next message
gaede
Messages: 233
Registered: January 2004
Location: DESY, Hamburg
Hi Fabio,


I don't know which program has created the LCIO file that you are analyzing - however if it is writing 'proper' LCIO your code should fail where you try to cast a raw hit to a SimCalorimeterHit:

 
LCObject* CSimhit = (*ch)->getRawHit();
SimCalorimeterHit* CRawHit = dynamic_cast<SimCalorimeterHit*>(CSimhit) ;

This is only by convention - of course technically one could create CalorimeterHit objects where the pointer to a raw hit is
in facta SimCalorimeterHit (clearly not a raw hit...).

The basic idea is that there is no direct link between MC Truth objects (MCParticle, SimCalorimeterHit, SimTrackerHit) and the higher level objects, such as Cluster, CalorimeterHit etc. but that rather LCRelation is used between the MC Truth hits (SimCalorimeterHit ) and the data hits ( CalorimeterHit ).

The following is a code snippet from the cluster cheater in MarlinReco (where the digitizer module has created the LCRelation collection beforehand):

 relcol = evt->getCollection("RelationCaloHit");

 LCRelationNavigator navigate(relcol); 

for (unsigned int j=0; j < col->getNumberOfElements(); ++j) 
{
  CalorimeterHit * calhit = dynamic_cast<CalorimeterHit*>(col->getElementAt(j));
  
  const LCObjectVec& objectVec = navigate.getRelatedToObjects(calhit);
  
  if (objectVec.size() > 0) {
    
    SimCalorimeterHit * simhit=dynamic_cast<SimCalorimeterHit*>( objectVec[0] );


    // use simhit ...

  }
}


I hope this helps.

Cheers, Frank.


PS: Guilherme, I don't think your suggestion can work since dynamic_cast<SimCaloHit*> will allways be 0 for anything that is not a SimCaloHit such as CaloHit...


Re: MC association [message #478 is a reply to message #477] Fri, 19 May 2006 09:09 Go to previous messageGo to next message
lima
Messages: 47
Registered: May 2004
Location: DeKalb, IL, USA
Frank,

You are right, the SimCalorimeterHit object is not a CalorimeterHit in the LCIO library. In org.lcsim, SimCalorimeterHit extends CalorimeterHit indeed, and I had this relationship in mind when answering Fabio's question.

Sad Sorry for the inconvenience,
Guilherme
Re: MC association [message #481 is a reply to message #477] Mon, 22 May 2006 07:06 Go to previous messageGo to next message
fabio
Messages: 3
Registered: April 2006
Hi Frank,
your suggestion has been really useful..
Thanks to LCRelation I can access to SimCalorimeterHit infos.
I have now a new question:
in the corresponding class of SimCalorimeterHit there are several memb functions related to the possibility to get the number of MC contributions to the hit and the PDG code of the shower particle that caused each contribution.
The problem is that the number of MC contributions to the hit can be sometimes really large..
Among all these contributions how is it possible (if it is possible..) to establish the only MC particle responsible of the hit?


Thanks again,
Fabio.

Re: MC association [message #482 is a reply to message #481] Mon, 22 May 2006 08:44 Go to previous messageGo to next message
gaede
Messages: 233
Registered: January 2004
Location: DESY, Hamburg
Quote:

Among all these contributions how is it possible (if it is possible..) to establish the only MC particle responsible of the hit?


There can (and often will) be more than one particle contributing to one cell, i.e. one SimCalorimeterHit. If the collection parameter LCIO.CHBIT_PDG==0 then every MCContribution corresponds to one particle entering the calorimeter. To define the one that you want to use in your relation depends on your algorithm - you might simply take the particle with the largest contribution or you could keep the weights for your MC truth assignment.
If LCIO.CHBIT_PDG==1 then all simulator steps' contributions to the hit are stored, i.e. every shower particle is stored with it's PDG and energy - however the link getParticleCont() will still be to the MCParticle that entered the calorimeter - in this case you have to sum up all contributions for every particle to find the one that dominates...
The following code snippet should work in either case:

	  SimCalorimeterHit* sh =
	    dynamic_cast<SimCalorimeterHit*>( calVec->getElementAt(j) ) ; 

	  typedef std::map<MCParticle*,double>  MCPMap ; 
	  
	  MCPMap mcpMap ;
	  
	  for( int ii=0 ;ii< sh->getNMCContributions() ; ++ii){
	    mcpMap[  sh->getParticleCont(ii) ] += sh->getEnergyCont(ii) ;   
	  }
	  
	  double eMax(0.) ; 
	  MCParticle* mcp ;
	  
	  for( MCPMap::iterator it = mcpMap.begin() ;  it != mcpMap.end() ; ++it ){
	    
	    if( it->second > eMax ) {
	      mcp = it->first ;
	      eMax = it->second ;
	    }
	  } 
	  std::cout << " largest contribution " << eMax << " GeV from particle " 
		    << mcp << std::endl ;



Probably code like this should go into some sort of utility method in a future LCIO release ...

Frank.
Re: MC association [message #492 is a reply to message #482] Mon, 29 May 2006 04:51 Go to previous message
fabio
Messages: 3
Registered: April 2006
Hi frank,

your suggestion works really fine for a MC association between reco particles (with associated reco cluster) and simulated ones..
Immagine now to arrange a similar MC association for reco particles with tracks (for example electrons).
Using the same approach it is possible to follow a similar chain to get associated sim info:

reco electrons -> reco trak associated -> reco track hits -> sim track hits -> MC particle

Again, in principle, it can be possible that for the same reco track more than one MC particle is associated.
how to establish the only one responsible of the track hit?

Is it possible to study again a contribution in energy/momentum?
Or, we should simply analyse the frequency of appearance of a particular MC particle among all related to differents sim hits?

Is somewhere already an implemented code that can work about it?

Thanks again,
Fabio

[Updated on: Mon, 29 May 2006 05:09]

Previous Topic:Problem with Marlin 0.9.4
Next Topic:org.lcsim Tracking Examples
Goto Forum:
  

[ PDF ]

Current Time: Thu Oct 17 06:25:41 Pacific Daylight Time 2019
.:: Contact :: Home ::.

Powered by: FUDforum 3.0.1.
Copyright ©2001-2010 FUDforum Bulletin Board Software