Commit f7f9553e authored by Atlas-Software Librarian's avatar Atlas-Software Librarian Committed by Graeme Stewart
Browse files

'CMakeLists.txt' (InDetTruthTools-00-01-01)

parent 9a59215f
################################################################################
# Package: InDetTruthTools
################################################################################
# Declare the package name:
atlas_subdir( InDetTruthTools )
# Declare the package's dependencies:
atlas_depends_on_subdirs( PUBLIC
Control/AthenaBaseComps
InnerDetector/InDetTruth/InDetTruthInterfaces
Tracking/TrkEvent/TrkEventPrimitives
Tracking/TrkTruthTracks/TrkTruthTrackInterfaces
PRIVATE
DetectorDescription/AtlasDetDescr
GaudiKernel
InnerDetector/InDetDetDescr/InDetIdentifier
InnerDetector/InDetDetDescr/InDetReadoutGeometry
InnerDetector/InDetRawEvent/InDetSimData
InnerDetector/InDetRecEvent/InDetPrepRawData
Tracking/TrkEvent/TrkPrepRawData
Tracking/TrkEvent/TrkTruthData )
# External dependencies:
find_package( HepMC )
# Component(s) in the package:
atlas_add_component( InDetTruthTools
src/*.cxx
src/components/*.cxx
INCLUDE_DIRS ${HEPMC_INCLUDE_DIRS}
LINK_LIBRARIES ${HEPMC_LIBRARIES} AthenaBaseComps TrkEventPrimitives AtlasDetDescr GaudiKernel InDetIdentifier InDetReadoutGeometry InDetSimData InDetPrepRawData TrkPrepRawData TrkTruthData )
# Install files from the package:
atlas_install_headers( InDetTruthTools )
atlas_install_joboptions( share/*.py )
/*
Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
*/
// A tool to add a PRD to PRD_MultiTruthCollection.
#ifndef PRD_MULTITRUTHBUILDER_H
#define PRD_MULTITRUTHBUILDER_H
#include "AthenaBaseComps/AthAlgTool.h"
#include "InDetTruthInterfaces/IPRD_MultiTruthBuilder.h"
class PixelID;
namespace InDet {
class PRD_MultiTruthBuilder : virtual public IPRD_MultiTruthBuilder,
public AthAlgTool
{
public:
PRD_MultiTruthBuilder(const std::string& type, const std::string& name, const IInterface* parent);
virtual void addPrepRawDatum(PRD_MultiTruthCollection *prdTruth,
const Trk::PrepRawData* prd,
const InDetSimDataCollection* simDataMap,
bool pixels
);
virtual StatusCode initialize();
private:
const PixelID *m_idHelperPixel;
};
}
#endif/*PRD_MULTITRUTHBUILDER_H*/
/*
Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
*/
#ifndef INDETTRUTHTOOLS_PRD_TRUTHTRAJECTORYMANIPULATORID_H
#define INDETTRUTHTOOLS_PRD_TRUTHTRAJECTORYMANIPULATORID_H
///////////////////////////////////////////////////////////////////
// PRD_TruthTrajectoryManipulatorID.h, (c) ATLAS Detector software
///////////////////////////////////////////////////////////////////
#include "AthenaBaseComps/AthAlgTool.h"
#include "TrkTruthTrackInterfaces/IPRD_TruthTrajectoryManipulator.h"
class AtlasDetectorID;
namespace InDet {
/** @class PRD_TruthTrajectoryManipulatorID
--- description what this too does ! ----
@author : roland Jansky < email addy >
*/
class PRD_TruthTrajectoryManipulatorID : public AthAlgTool, virtual public Trk::IPRD_TruthTrajectoryManipulator {
public:
//** Constructor with parameters */
PRD_TruthTrajectoryManipulatorID( const std::string& t, const std::string& n, const IInterface* p );
/** Tool Destructor */
~PRD_TruthTrajectoryManipulatorID() {}
// Athena algtool's Hooks
StatusCode initialize();
StatusCode finalize();
/** Interface method from IPRD_TruthTrajectoryManipulator */
virtual bool manipulateTruthTrajectory( Trk::PRD_TruthTrajectory &) const;
private:
/**ID pixel helper*/
const AtlasDetectorID* m_atlasId;
};
}
#endif
/*
Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
*/
#ifndef INDETTRUTHTOOLS_PRD_TRUTHTRAJECTORYSELECTORID_H
#define INDETTRUTHTOOLS_PRD_TRUTHTRAJECTORYSELECTORID_H
///////////////////////////////////////////////////////////////////
// PRD_TruthTrajectorySelectorID.h, (c) ATLAS Detector software
///////////////////////////////////////////////////////////////////
#include "AthenaBaseComps/AthAlgTool.h"
#include "TrkTruthTrackInterfaces/IPRD_TruthTrajectorySelector.h"
//#include "EventPrimitives/EventPrimitives.h"
//#include "TrkEventPrimitives/LocalPosition.h"
//#include "TrkEventPrimitives/GlobalPosition.h"
#include "TrkEventPrimitives/LocalParameters.h"
//#include "CLHEP/Geometry/Transform3D.h"
class ThreePointCircle {
public :
/** The constructor from three points */
ThreePointCircle( const std::vector<Amg::Vector3D>& vecPoints);
/** Destructor */
~ThreePointCircle();
/** center */
const Amg::Vector2D& center() const;
/** Translation */
const Amg::Translation3D* frameTranslation() const;
/* d0, z0, eta, phi, pt */
double d0() const;
double z0() const;
double phi0() const;
double eta() const;
double pT() const;
double radius() const;
private :
void constructCircle(const Amg::Vector3D&, const Amg::Vector3D&, const Amg::Vector3D&);
// the reference point
Amg::Translation3D* m_translation;
// the parameters
double m_d0, m_z0, m_phi0, m_eta, m_pt;
double m_radius;
Amg::Vector2D m_center;
};
inline double ThreePointCircle::d0() const { return m_d0; }
inline double ThreePointCircle::z0() const { return m_z0; }
inline double ThreePointCircle::phi0() const { return m_phi0; }
inline double ThreePointCircle::eta() const { return m_eta; }
inline double ThreePointCircle::pT() const { return m_pt; }
inline double ThreePointCircle::radius() const { return m_radius; }
inline const Amg::Vector2D& ThreePointCircle::center() const { return m_center; }
inline const Amg::Translation3D* ThreePointCircle::frameTranslation() const { return m_translation; }
class AtlasDetectorID;
namespace InDet {
/** @class PRD_TruthTrajectorySelectorID
simple selector for truth trajectories
@author : Roland.Wolfgang.Jansky -at- cern.ch
*/
class PRD_TruthTrajectorySelectorID : public AthAlgTool, virtual public Trk::IPRD_TruthTrajectorySelector {
public:
//** Constructor with parameters */
PRD_TruthTrajectorySelectorID( const std::string& t, const std::string& n, const IInterface* p );
/** Tool Destructor */
~PRD_TruthTrajectorySelectorID() {}
// Athena algtool's Hooks
StatusCode initialize();
StatusCode finalize();
/** Interface method from IPRD_TruthTrajectorySelector */
virtual bool pass( const Trk::PRD_TruthTrajectory & ) const;
private:
/**ID pixel helper*/
const AtlasDetectorID* m_atlasId;
};
}
#endif
/*
Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
*/
#ifndef INDETTRUTHTOOLS_PRD_TRUTHTRAJECTORYSORTERID_H
#define INDETTRUTHTOOLS_PRD_TRUTHTRAJECTORYSORTERID_H
#include "AthenaBaseComps/AthAlgTool.h"
#include "TrkTruthTrackInterfaces/IPRD_TruthTrajectoryManipulator.h"
namespace InDet {
class PRD_TruthTrajectorySorterID : public AthAlgTool, virtual public Trk::IPRD_TruthTrajectoryManipulator {
public:
//** Constructor with parameters */
PRD_TruthTrajectorySorterID( const std::string& t, const std::string& n, const IInterface* p );
/** Tool Destructor */
~PRD_TruthTrajectorySorterID() {}
// Athena algtool's Hooks
StatusCode initialize();
StatusCode finalize();
virtual bool manipulateTruthTrajectory( Trk::PRD_TruthTrajectory &) const;
};
}
#endif
package InDetTruthTools
author Andrei Gaponenko <AGaponenko@lbl.gov>
public
use AtlasPolicy AtlasPolicy-*
use AthenaBaseComps AthenaBaseComps-* Control
use InDetTruthInterfaces InDetTruthInterfaces-* InnerDetector/InDetTruth
use TrkTruthTrackInterfaces TrkTruthTrackInterfaces-* Tracking/TrkTruthTracks
use TrkEventPrimitives TrkEventPrimitives-* Tracking/TrkEvent
private
use GaudiInterface GaudiInterface-* External
use TrkTruthData TrkTruthData-* Tracking/TrkEvent
use TrkPrepRawData TrkPrepRawData-* Tracking/TrkEvent
use InDetSimData InDetSimData-* InnerDetector/InDetRawEvent
use InDetPrepRawData InDetPrepRawData-* InnerDetector/InDetRecEvent
use InDetReadoutGeometry InDetReadoutGeometry-* InnerDetector/InDetDetDescr
use InDetIdentifier InDetIdentifier-* InnerDetector/InDetDetDescr
use AtlasDetDescr AtlasDetDescr-* DetectorDescription
use AtlasHepMC AtlasHepMC-* External
library InDetTruthTools *.cxx -s=components *.cxx
apply_pattern component_library
apply_pattern declare_joboptions files="*.py"
if not 'InDetTruthTools' in theApp.Dlls:
theApp.Dlls += ["InDetTruthTools"]
/*
Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
*/
#include "InDetTruthTools/PRD_MultiTruthBuilder.h"
#include "TrkTruthData/PRD_MultiTruthCollection.h"
#include "TrkPrepRawData/PrepRawData.h"
#include "InDetSimData/InDetSimDataCollection.h"
#include "InDetIdentifier/PixelID.h"
#include "InDetReadoutGeometry/SiDetectorElement.h"
#include "InDetReadoutGeometry/SiReadoutCellId.h"
#include <stdexcept>
#include <ext/functional>
namespace InDet {
#if 0 /*no-op block for better indentation in Emacs*/
}
#endif
//================================================================
PRD_MultiTruthBuilder::PRD_MultiTruthBuilder(const std::string& type, const std::string& name, const IInterface* parent)
: AthAlgTool(type,name,parent)
, m_idHelperPixel(0)
{
declareInterface<IPRD_MultiTruthBuilder>(this);
}
StatusCode PRD_MultiTruthBuilder::initialize() {
ATH_MSG_DEBUG("PRD truth builder initialization");
StatusCode sc = detStore()->retrieve(m_idHelperPixel, "PixelID");
if (!sc.isSuccess()) {
ATH_MSG_FATAL("Could not get Pixel ID helper");
return sc;
}
return sc;
}
//================================================================
void PRD_MultiTruthBuilder::addPrepRawDatum(PRD_MultiTruthCollection *prdTruth,
const Trk::PrepRawData* prd,
const InDetSimDataCollection* simDataMap,
bool pixels
)
{
if(!prdTruth || !prd || !simDataMap) {
return;
}
ATH_MSG_VERBOSE("addPrepRawDatum(): new PRD "<<prd<<", id="<<prd->identify()<<", number of RDOs: " << prd->rdoList().size());
bool gotSDO = false;
bool gotValidParticle = false;
// loop over RDOs
std::vector<Identifier>::const_iterator nextRDO = prd->rdoList().begin();
std::vector<Identifier>::const_iterator lastRDO = prd->rdoList().end();
for (; nextRDO!=lastRDO; ++nextRDO) {
InDetSimDataCollection::const_iterator iter(simDataMap->find(*nextRDO));
if ( pixels && iter == simDataMap->end() ) {
ATH_MSG_VERBOSE ( "Pixel=" << m_idHelperPixel->eta_index(*nextRDO) << "," <<m_idHelperPixel->phi_index(*nextRDO) << " does not match any SDO");
}
if(iter != simDataMap->end() ) {
gotSDO = true;
// Got an SDO. Try to associate the PRD to MC particles we have info about.
const InDetSimData& sdo = iter->second;
const std::vector< InDetSimData::Deposit >& deposits = sdo.getdeposits();
std::vector< InDetSimData::Deposit >::const_iterator nextdeposit = deposits.begin();
std::vector< InDetSimData::Deposit >::const_iterator lastdeposit = deposits.end();
for( ; nextdeposit!=lastdeposit; ++nextdeposit) {
const HepMcParticleLink& particleLink = nextdeposit->first;
ATH_MSG_VERBOSE("addPrepRawDatum(): Barcode " << particleLink.barcode());
if (particleLink.isValid()) {
gotValidParticle = true;
// Associate the particle to the PRD. But don't add duplicates.
// Note: it may be more efficient to filter out duplicates among particles for the current PRD, then check-and-add the reduced set to the large multimap.
// But may be not for the typically small RDO/PRD ratio.
typedef PRD_MultiTruthCollection::iterator truthiter;
std::pair<truthiter, truthiter> r = prdTruth->equal_range(prd->identify());
// FIXME: Is it OK to use the gcc (SGI) extensions of the STL?
if(r.second == std::find_if(r.first, r.second,
__gnu_cxx::compose1(std::bind2nd(std::equal_to<HepMcParticleLink>(), particleLink),
__gnu_cxx::select2nd<PRD_MultiTruthCollection::value_type>()
)
)
)
{
prdTruth->insert(std::make_pair(prd->identify(), particleLink));
}
}
}
}
}
if(gotSDO && !gotValidParticle) {
// Looked at all the deposits from all the SDOs, but did not find any valid particle link.
//prdTruth->insert(std::make_pair(prd, particleLinkUnknown));
ATH_MSG_VERBOSE("addPrepRawDatum(): got SDO but no particles");
}
// No SDO found after all, for all the RDOs? Then (and only then) this PRD is noise.
//if(!gotSDO) {
//prdTruth->insert(std::make_pair(prd, particleLinkNoise) );
//}
}
//================================================================
} // namespace InDet
/*
Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
*/
///////////////////////////////////////////////////////////////
// PRD_TruthTrajectoryManipulatorID.cxx, (c) ATLAS Detector software
///////////////////////////////////////////////////////////////////
#include "InDetTruthTools/PRD_TruthTrajectoryManipulatorID.h"
#include "InDetPrepRawData/SCT_Cluster.h"
#include "InDetPrepRawData/TRT_DriftCircle.h"
#include "InDetPrepRawData/PixelCluster.h"
// DetectorDescription
#include "AtlasDetDescr/AtlasDetectorID.h"
// HepMC
#include "HepMC/GenParticle.h"
InDet::PRD_TruthTrajectoryManipulatorID::PRD_TruthTrajectoryManipulatorID(const std::string& t, const std::string& n, const IInterface* p) :
AthAlgTool(t,n,p)
{
declareInterface<Trk::IPRD_TruthTrajectoryManipulator>(this);
}
StatusCode InDet::PRD_TruthTrajectoryManipulatorID::initialize() {
ATH_MSG_VERBOSE("Initializing ...");
StatusCode sc = detStore()->retrieve(m_atlasId, "AtlasID");
if (sc.isFailure()) {
msg(MSG::ERROR) << "Could not get AtlasID helper !" << endreq;
return StatusCode::FAILURE;
}
return StatusCode::SUCCESS;
}
StatusCode InDet::PRD_TruthTrajectoryManipulatorID::finalize() {
ATH_MSG_VERBOSE("Finalizing ...");
return StatusCode::SUCCESS;
}
bool InDet::PRD_TruthTrajectoryManipulatorID::manipulateTruthTrajectory( Trk::PRD_TruthTrajectory &prdvec) const {
if( (*prdvec.genParticle).barcode() < 100000){
srand( static_cast< unsigned int >( time( 0 ) ) );
const int pdg_id = (*prdvec.genParticle).pdg_id();
const double prob_pix = pdg_id == 2212 ? 4. : 0;
const double prob_sct = 4.375;
std::vector<const Trk::PrepRawData* >::iterator prdIter = prdvec.prds.begin();
std::vector<const Trk::PrepRawData* >::iterator prdIterE = prdvec.prds.end();
while( prdIter != prdIterE ){
if( m_atlasId->is_pixel((*prdIter)->identify()) ){
if( prob_pix > 0 ? rand()%100 <= prob_pix : false ){
prdIter = prdvec.prds.erase(prdIter);
prdIterE = prdvec.prds.end();
}
else ++ prdIter;
}
else if( m_atlasId->is_sct((*prdIter)->identify()) ){
if( rand()%100 <= prob_sct ){
prdIter = prdvec.prds.erase(prdIter);
prdIterE = prdvec.prds.end();
}
else ++ prdIter;
}
else ++ prdIter;
}
}
return true;
}
/*
Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
*/
///////////////////////////////////////////////////////////////
// PRD_TruthTrajectorySelectorID.cxx, (c) ATLAS Detector software
///////////////////////////////////////////////////////////////////
#include "InDetTruthTools/PRD_TruthTrajectorySelectorID.h"
#include "InDetPrepRawData/SCT_Cluster.h"
#include "InDetPrepRawData/TRT_DriftCircle.h"
#include "InDetPrepRawData/PixelCluster.h"
// DetectorDescription
#include "AtlasDetDescr/AtlasDetectorID.h"
// HepMC
#include "HepMC/GenParticle.h"
//for ThreePointCircle
//#include "TrkEventPrimitives/GlobalDirection.h"
InDet::PRD_TruthTrajectorySelectorID::PRD_TruthTrajectorySelectorID(const std::string& t, const std::string& n, const IInterface* p) :
AthAlgTool(t,n,p)
{
declareInterface<Trk::IPRD_TruthTrajectorySelector>(this);
}
StatusCode InDet::PRD_TruthTrajectorySelectorID::initialize() {
ATH_MSG_VERBOSE("Initializing ...");
StatusCode sc = detStore()->retrieve(m_atlasId, "AtlasID");
if (sc.isFailure()) {
msg(MSG::ERROR) << "Could not get AtlasID helper !" << endreq;
return StatusCode::FAILURE;
}
return StatusCode::SUCCESS;
}
StatusCode InDet::PRD_TruthTrajectorySelectorID::finalize() {
ATH_MSG_VERBOSE("Finalizing ...");
return StatusCode::SUCCESS;
}
ThreePointCircle::ThreePointCircle(const std::vector<Amg::Vector3D>& posVec) :
m_translation(0), m_d0(0.), m_z0(0.), m_phi0(0.), m_eta(0.), m_pt(0.)
{
if (posVec.size() <3)
std::cout << "[WARNING] not enough points" << std::endl;
else
constructCircle(posVec[0],posVec[1],posVec[2]);
}
ThreePointCircle::~ThreePointCircle()
{
delete m_translation;
}
void ThreePointCircle::constructCircle( const Amg::Vector3D& p1,
const Amg::Vector3D& p2,
const Amg::Vector3D& p3) {
// get the (x/y) coordinates of all points
double translationX = m_translation ? -(m_translation->translation()).x() : 0.;
double translationY = m_translation ? -(m_translation->translation()).y() : 0.;
double bx = p1.x()+translationX;
double by = p1.y()+translationY;
double cx = p2.x()+translationX;
double cy = p2.y()+translationY;
double dx = p3.x()+translationX;
double dy = p3.y()+translationY;
double temp = cx*cx+cy*cy;
// prepare the matrix for the linear equation
double bc = (bx*bx + by*by - temp)/2.0;
double cd = (temp - dx*dx - dy*dy)/2.0;
double det = (bx-cx)*(cy-dy)-(cx-dx)*(by-cy);
//if (fabs(det) < 1.0e-6) {
//std::cout << "[WARNING] Determinant is close to 0 " << std::endl;
//}
// invert the determinant
det = 1./det;
// x / y coordinate of the center
double cRx= (bc*(cy-dy)-cd*(by-cy))*det;
double cRy = ((bx-cx)*cd-(cx-dx)*bc)*det;
// get the center, radius and momentum
m_center = Amg::Vector2D(cRx,cRy);
m_radius = sqrt((cRx-bx)*(cRx-bx)+(cRy-by)*(cRy-by));
m_pt = m_radius*(0.3*2.);
// transverse parameters
m_d0 = sqrt(cRx*cRx+cRy*cRy)-m_radius;
// longitudinal parameters
double theta = 1./3*(p1.theta() + p2.theta() + p3.theta());
double