diff --git a/InnerDetector/InDetAlignTools/InDetAlignGenTools/InDetAlignGenTools/AlignSiModuleList.h b/InnerDetector/InDetAlignTools/InDetAlignGenTools/InDetAlignGenTools/AlignSiModuleList.h new file mode 100755 index 0000000000000000000000000000000000000000..3e66745e8fbfd4d882f18aee6f983ad3a881f904 --- /dev/null +++ b/InnerDetector/InDetAlignTools/InDetAlignGenTools/InDetAlignGenTools/AlignSiModuleList.h @@ -0,0 +1,30 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +#ifndef INDETALIGNGENALGS_ALIGNSIMODULELIST_H +#define INDETALIGNGENALGS_ALIGNSIMODULELIST_H +// AlignSiModuleList.h +// Richard Hawkings, started 15/4/04 +// Encapsulation of vector of AlignSiModules and modlist +// for InDetAlignGenAlgs allowing it to be stored in TDS + +#ifndef CLIDSVC_CLASSDEF_H +#include "CLIDSvc/CLASS_DEF.h" +#endif +#include <vector> +#include <map> +#include "InDetAlignTrkInfo/AlignSiModule.h" +#include "Identifier/Identifier.h" + +class AlignSiModuleList { + public: + std::vector<AlignSiModule> vec; + typedef std::less<Identifier> lessp; + typedef std::map<Identifier,int,lessp> ModuleIDMap; + ModuleIDMap idmap; +}; + +CLASS_DEF( AlignSiModuleList , 230164606 , 1 ) + +#endif // INDETALIGNGENALGS_ALIGNSIMODULELIST_H diff --git a/InnerDetector/InDetAlignTools/InDetAlignGenTools/InDetAlignGenTools/IInDetAlignDBTool.h b/InnerDetector/InDetAlignTools/InDetAlignGenTools/InDetAlignGenTools/IInDetAlignDBTool.h new file mode 100755 index 0000000000000000000000000000000000000000..3f55c797a5541eabf7a6c6444631403d5aad8252 --- /dev/null +++ b/InnerDetector/InDetAlignTools/InDetAlignGenTools/InDetAlignGenTools/IInDetAlignDBTool.h @@ -0,0 +1,70 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +#ifndef INDETALIGNTOOLS_ALIGNDBTOOL_IH +#define INDETALIGNTOOLS_ALIGNDBTOOL_IH +// IInDetAlignDBTool.h +// an AlgTool to manage the inner detector alignment database classes +// abstract interface class, actual implementation and doc +// is in InDetAlignDBTool.h +// Richard Hawkings, started 11/4/05 + +#include "GaudiKernel/IAlgTool.h" +#include "EventPrimitives/EventPrimitives.h" +#include "GeoPrimitives/GeoPrimitives.h" + + +class Identifier; + +static const InterfaceID +IID_INDETALIGN_IInDetAlignDBTool("IInDetAlignDBTool",1,0); + +class IInDetAlignDBTool: virtual public IAlgTool { + public: + static const InterfaceID& interfaceID(); + + virtual void createDB() const =0; + virtual void dispGroup(const int, const int, const int, const int, const int, + const float, const float, const float, + const int, const int, const int) const =0; + + virtual void writeFile(const bool, const std::string) const =0; + virtual void readTextFile(const std::string) const =0; + virtual void readNtuple(const std::string) const =0; + + virtual bool idToDetSet(const Identifier, + int&,int&,int&,int&,int&,int&) const =0; + virtual std::string dirkey(const Identifier&, const int) const =0; + virtual std::string dirkey(const int,const int,const int, const int) const =0; + virtual std::string dirkey(const int,const int,const int, const int, const int) const =0; // new function + virtual std::string DBMkey(const int,const int,const int, const int) const =0; // new function + + virtual bool setTrans(const Identifier&, const int, const Amg::Transform3D& ) + const =0; + virtual bool setTrans(const Identifier& ident, const int level, + const Amg::Vector3D & translate, double alpha, double beta, double gamma) const = 0; + virtual bool tweakTrans(const Identifier&, const int, const Amg::Transform3D&) + const =0; + virtual bool tweakTrans(const Identifier& ident, const int level, + const Amg::Vector3D& translate, double alpha, + double beta, double gamma) const = 0; + virtual Identifier getL1L2fromL3Identifier( const Identifier& ident + , const int& level + ) const=0 ; + virtual Amg::Transform3D getTransL123( const Identifier& ident ) const=0 ; + virtual Amg::Transform3D getTrans(const Identifier&, const int) const=0; + virtual StatusCode outputObjs() const=0; + virtual void fillDB(const std::string, const unsigned int,const unsigned int, + const unsigned int, const unsigned int) const=0; + virtual void printDB(const int) const=0; + virtual void sortTrans() const=0; + virtual void extractAlphaBetaGamma(const Amg::Transform3D & trans, + double& alpha, double& beta, double &gamma) const=0; +}; + +inline const InterfaceID& IInDetAlignDBTool::interfaceID() +{ return IID_INDETALIGN_IInDetAlignDBTool; } + + +#endif // INDETALIGNTOOLS_ALIGNDBTOOL_IH diff --git a/InnerDetector/InDetAlignTools/InDetAlignGenTools/InDetAlignGenTools/IInDetAlignFillSiCluster.h b/InnerDetector/InDetAlignTools/InDetAlignGenTools/InDetAlignGenTools/IInDetAlignFillSiCluster.h new file mode 100644 index 0000000000000000000000000000000000000000..0439277c84d6dacf8fcbb4d144066d6c7c6498c0 --- /dev/null +++ b/InnerDetector/InDetAlignTools/InDetAlignGenTools/InDetAlignGenTools/IInDetAlignFillSiCluster.h @@ -0,0 +1,28 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +#ifndef INDETALIGNTOOLS_INDETALIGNFILLSICLUSTER_IH +#define INDETALIGNTOOLS_INDETALIGNFILLSICLUSTER_IH +// IInDetAlignFillSiCluster.h +// is in InDetAlignFillSiCluster.h +// Carlos Escobar, started 08/03/2008 + +#include "GaudiKernel/IAlgTool.h" + +static const InterfaceID + IID_INDETALIGN_IInDetAlignFillSiCluster("IInDetAlignFillSiCluster",1,0); + +class IInDetAlignFillSiCluster: virtual public IAlgTool { + public: + static const InterfaceID& interfaceID(); + + virtual StatusCode FillSiCluster() = 0; + +}; + +inline const InterfaceID& IInDetAlignFillSiCluster::interfaceID() +{ return IID_INDETALIGN_IInDetAlignFillSiCluster; } + + +#endif // INDETALIGNTOOLS_INDETALIGNFILLSICLUSTER_IH diff --git a/InnerDetector/InDetAlignTools/InDetAlignGenTools/InDetAlignGenTools/IInDetAlignFillTrack.h b/InnerDetector/InDetAlignTools/InDetAlignGenTools/InDetAlignGenTools/IInDetAlignFillTrack.h new file mode 100644 index 0000000000000000000000000000000000000000..1e0cf75dd2988b114135d0139bcd06618c3e9024 --- /dev/null +++ b/InnerDetector/InDetAlignTools/InDetAlignGenTools/InDetAlignGenTools/IInDetAlignFillTrack.h @@ -0,0 +1,34 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +#ifndef INDETALIGNTOOLS_INDETALIGNFILLTRACK_IH +#define INDETALIGNTOOLS_INDETALIGNFILLTRACK_IH +// IInDetAlignFillTrack.h +// is in InDetAlignFillTrack.h +// Carlos Escobar, started 27/12/2007 + +#include "GaudiKernel/IAlgTool.h" + +static const InterfaceID + IID_INDETALIGN_IInDetAlignFillTrack("IInDetAlignFillTrack",1,0); + +class IInDetAlignFillTrack: virtual public IAlgTool { + public: + static const InterfaceID& interfaceID(); + + virtual StatusCode FillTrack() = 0; + + virtual int GetTrks() const = 0; + virtual int GetTrkHits() const = 0; + virtual int GetTrkPixHits() const = 0; + virtual int GetTrkSCTHits() const = 0; + virtual int GetTrkTRTHits() const = 0; + +}; + +inline const InterfaceID& IInDetAlignFillTrack::interfaceID() +{ return IID_INDETALIGN_IInDetAlignFillTrack; } + + +#endif // INDETALIGNTOOLS_INDETALIGNFILLTRACK_IH diff --git a/InnerDetector/InDetAlignTools/InDetAlignGenTools/InDetAlignGenTools/IInDetAlignHitQualSelTool.h b/InnerDetector/InDetAlignTools/InDetAlignGenTools/InDetAlignGenTools/IInDetAlignHitQualSelTool.h new file mode 100644 index 0000000000000000000000000000000000000000..0705fdca91d825bd5c57f6b8f612ef19db7544e5 --- /dev/null +++ b/InnerDetector/InDetAlignTools/InDetAlignGenTools/InDetAlignGenTools/IInDetAlignHitQualSelTool.h @@ -0,0 +1,43 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +#ifndef INDETALIGNTOOLS_IHITQUALSELTOOL_IH +#define INDETALIGNTOOLS_IHITQUALSELTOOL_IH + +#include "GaudiKernel/IAlgTool.h" + +/** Must declare this, with name of interface*/ +static const InterfaceID IID_INDETALIGN_IInDetAlignHitQualSelTool( "IInDetAlignHitQualSelTool", 1, 0 ) ; + +/** @class InDetAlignHitQualSelTool + @brief The InDetAlignHitQualSelTool is to select good quality hits for alignment to + build residuals with possible cuts on outlier hits, hits which are too + large clusters, hits with large incidence angles, etc. + @author Oleg Brandt <http://consult.cern.ch/xwho> + */ + // Many thanks to Carlos Escobar and Sebastian Fleischmann! + +namespace Trk { + class TrackStateOnSurface ; + class RIO_OnTrack ; +} + +class IInDetAlignHitQualSelTool : virtual public IAlgTool { + public: + /** for ToolHandle functionality */ + static const InterfaceID& interfaceID() ; + + /** main method: from a TrackStateOnSurface select a good hit cutting on outlier hits, hits with + too many pixels/strips, hits with large incidence angles */ + virtual const Trk::RIO_OnTrack* getGoodHit( const Trk::TrackStateOnSurface* tsos ) const = 0 ; + /** from a TrackStateOnSurface select a good hole in track cutting on + large incidence angles only */ + virtual bool getGoodHole( const Trk::TrackStateOnSurface* tsos ) const = 0 ; +} ; + +inline const InterfaceID& IInDetAlignHitQualSelTool::interfaceID() { + return IID_INDETALIGN_IInDetAlignHitQualSelTool ; +} + +#endif // INDETALIGNTOOLS_IHITQUALSELTOOL_IH diff --git a/InnerDetector/InDetAlignTools/InDetAlignGenTools/InDetAlignGenTools/IInDetAlignOverlapTool.h b/InnerDetector/InDetAlignTools/InDetAlignGenTools/InDetAlignGenTools/IInDetAlignOverlapTool.h new file mode 100644 index 0000000000000000000000000000000000000000..34cbcf694f23630a9c17abbc3b2327a1a9ded144 --- /dev/null +++ b/InnerDetector/InDetAlignTools/InDetAlignGenTools/InDetAlignGenTools/IInDetAlignOverlapTool.h @@ -0,0 +1,44 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +#ifndef INDETALIGNTOOLS_IOVERLAPTOOL_IH +#define INDETALIGNTOOLS_IOVERLAPTOOL_IH + +#include "GaudiKernel/IAlgTool.h" +#include "InDetAlignGenTools/IInDetAlignFillTrack.h" +//#include "DataModel/DataVector.h" +#include <vector> +/** @class InDetAlignOverlapTool + @brief The InDetAlignOverlapTool is to select good quality tracks for alignment to + build residuals with possible cuts on overlap hits, overlaps in SCT and PIXELs (in Phi and Eta). +*/ + +class AlignTrk; +class AlignSiHit; + +static const InterfaceID +IID_INDETALIGN_IInDetAlignOverlapTool( "IInDetAlignOverlapTool", 1, 0 ) ; + + +class IInDetAlignOverlapTool : virtual public IAlgTool { + public: + static const InterfaceID& interfaceID() ; + + /** main method: from a TrackStateOnSurface select an overlap hit */ + virtual int getNumberOverlapPIX( const AlignTrk& )const = 0 ; + virtual int getNumberOverlapSCT( const AlignTrk& )const = 0 ; + // virtual const DataVector<AlignSiHit>* getOverlapHit(const AlignTrk& )const = 0 ; + //virtual std::vector<AlignSiHit> getOverlapHit(const AlignTrk& ) =0;//const = 0 ; + virtual std::vector<AlignSiHit> getOverlapHit(const AlignTrk& ) = 0;//const = 0 ; + + +} ; + +inline const InterfaceID& IInDetAlignOverlapTool::interfaceID() { + return IID_INDETALIGN_IInDetAlignOverlapTool ; +} + +#endif // INDETALIGNTOOLS_IOVERLAPTOOL_IH + + diff --git a/InnerDetector/InDetAlignTools/InDetAlignGenTools/InDetAlignGenTools/IInDetAlignTrackSelTool.h b/InnerDetector/InDetAlignTools/InDetAlignGenTools/InDetAlignGenTools/IInDetAlignTrackSelTool.h new file mode 100755 index 0000000000000000000000000000000000000000..bc1cb48be843d796acf353a782d44bdb4cdd36dc --- /dev/null +++ b/InnerDetector/InDetAlignTools/InDetAlignGenTools/InDetAlignGenTools/IInDetAlignTrackSelTool.h @@ -0,0 +1,35 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +// IInDetAlignTrackSelTool.h + +// Sergio Gonzalez Sevilla, started 04/7/05 +// Miguel Olivo Gomez, extended 07/6/06 + +// AlgTool to select high quality tracks for the inner detector +// (Pixel+SCT) alignment algorithms. +// This AlgTool provides a track selection based on the following cut variables: +// Momentum, pt, number of shared hits, number of holes and chi2 probability. +// Returns 0 in case a track is not accepted, otherwise 1 + +#ifndef INDETALIGNTOOLS_ALIGNTRACKSELTOOL_IH +#define INDETALIGNTOOLS_ALIGNTRACKSELTOOL_IH + +#include "GaudiKernel/IAlgTool.h" +#include "TrkTrack/Track.h" + +static const InterfaceID + IID_INDETALIGN_IInDetAlignTrackSelTool("IInDetAlignTrackSelTool",1,0); + +class IInDetAlignTrackSelTool: virtual public IAlgTool { + public: + static const InterfaceID& interfaceID(); + + virtual int getStatus(const Trk::Track&) const = 0; +}; + +inline const InterfaceID& IInDetAlignTrackSelTool::interfaceID() +{ return IID_INDETALIGN_IInDetAlignTrackSelTool; } + +#endif // INDETALIGNTOOLS_ALIGNTRACKSELTOOL_IH diff --git a/InnerDetector/InDetAlignTools/InDetAlignGenTools/InDetAlignGenTools/InDetAlignDBTool.h b/InnerDetector/InDetAlignTools/InDetAlignGenTools/InDetAlignGenTools/InDetAlignDBTool.h new file mode 100755 index 0000000000000000000000000000000000000000..09a552aaeb031054eaf9b661578ec22469ea60fc --- /dev/null +++ b/InnerDetector/InDetAlignTools/InDetAlignGenTools/InDetAlignGenTools/InDetAlignDBTool.h @@ -0,0 +1,191 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +#ifndef INDETALIGNGENTOOLS_ALIGNDBTOOL_H +#define INDETALIGNGENTOOLS_ALIGNDBTOOL_H +// InDetAlignDBTool.h +// an AlgTool to manage the inner detector alignment database classes +// Richard Hawkings, started 8/4/05 +// +// This AlgTool provides utilities to manage the ID alignment conditions +// data in the transient store, and to read and write it in a variety of +// formats +// Tool methods provided: +// createDB - create a 'null' set of AlignableTransforms, with no shifts +// dispGroup - apply a pattern to displace a group of modules randomly +// or systematically +// printDB - printout details of all the module transforms +// writeFile - export database contents on a text file or ntuple +// readTextFile - import database from text file +// readNtuple - import database from ntuple +// +// The following utility methods are also exposed +// idToDetSet - convert an SCT/pixel Identifier to set of integers +// specifying detector,barrel/ec,layer,ring, sector (and side for SCT) +// dirkey - return the AlignableTransform TDS key for a given module +// (by Identifier or set of integers) +// setTrans - set the transformation for a particular module, identified by +// a given Identifier +// tweakTrans - tweak (add to existing transform) the transformation for +// a particular module, identified by a given Identifier + +#include <vector> +#include <string> +#include "AthenaBaseComps/AthAlgTool.h" +#include "GaudiKernel/ServiceHandle.h" +#include "InDetAlignGenTools/IInDetAlignDBTool.h" + +class PixelID; +class SCT_ID; +class AlignableTransform; + +namespace InDetDD { + class PixelDetectorManager; + class SCT_DetectorManager; +} + +class InDetAlignDBTool: virtual public IInDetAlignDBTool, public AthAlgTool { + public: + InDetAlignDBTool(const std::string& type, const std::string& name, + const IInterface* parent); + virtual ~InDetAlignDBTool(); + + virtual StatusCode initialize(); + virtual StatusCode finalize(); + + // tool methods + // create null database in TDS for subsequent manipulation + + virtual void createDB() const; + + // displace a group of modules + // modules to effect are specified by dettype (1=pixel, 2=SCT, -1 both), + // bec (barrel=0, 2 for both endcaps, -1 for all), layer (0 to n-1, -1=all) + // ring (-1=all) and sector (-1=all) + // rphidisp,rdisp and zdisp specify the offsets to be used + // syst specifies what will be done, 1=systematic displacement of affected + // modules, 2=random (rphidisp etc are RMS), 3,4=systematic and random + // displacements but interpreting r/phi/zdisp as x/y/z, 5 chooses the + // same randomised numbers for all modules; + // level specifies the level oill look for the right balance of CERN staff supervisors and ATLAS + //USERS supervisors; f transformations affected + // skip=n optionally skips n random numbers to enable different patterns + virtual void dispGroup(const int dettype, const int bec, const int layer, + const int ring, const int sector, + const float rphidisp, const float rdisp, const float zdisp, + const int syst, const int level, const int skip) const; + + // write database contents to flat text file or ntuple + // for flat file, file gives filename, for ntuple, file is e.g. + // /NTUPLES/FILE1 + virtual void writeFile(const bool ntuple, const std::string file) const; + + // read back database from text file + virtual void readTextFile(const std::string file) const; + + // read back database from text file + void readOldTextFile(const std::string file) const; + + // read back database from ntuple + virtual void readNtuple(const std::string file) const; + + // convert an Identifier to a set of integers specifying detector, barrel/ec + // (0 for barrel, -1 for endcap C, +1 for endcap A), layer, ring (eta), + // sector (phi), and side (0/1, pixel always 0) + // same as the fields of the Identifier, except bec is +-1 not +-2 for endcap + virtual bool idToDetSet(const Identifier ident, + int& det,int& bec,int& layer,int& ring,int& sector,int& side) const; + + // return the AlignableTransform name where the transform for the given + // module can be found, 3 levels corresponding to alignment hierarchy + virtual std::string dirkey(const Identifier&, const int) const; + virtual std::string dirkey(const int,const int,const int, const int) const; + virtual std::string dirkey(const int,const int,const int, const int, const int) const; + virtual std::string DBMkey(const int,const int,const int, const int) const; + + // set a particular transform specified by Identifier + virtual bool setTrans(const Identifier& ident, const int level, + const Amg::Transform3D& trans) const; + + // As above but takes translation and rotations alpha, beta, gamma (rotations around x,y,z axes) as input. + // Calculates transform as HepGeom::Translate3D(translate) * HepGeom::RotateX3D(alpha) * HepGeom::RotateY3D(beta) * HepGeom::RotateZ3D(gamma) + virtual bool setTrans(const Identifier& ident, const int level, + const Amg::Vector3D& translate, double alpha, double beta, double gamma) const; + + + // tweak a particular transform specified by Identifier (i.e. add to + // already existing transformation) + virtual bool tweakTrans(const Identifier& ident, const int level, + const Amg::Transform3D& trans) const; + + // As above but takes translation and rotations alpha, beta, gamma (rotations around x,y,z axes) as input. + // Calculates transform as HepGeom::Translate3D(translate) * HepGeom::RotateX3D(alpha) * HepGeom::RotateY3D(beta) * HepGeom::RotateZ3D(gamma) + virtual bool tweakTrans(const Identifier& ident, const int level, + const Amg::Vector3D& translate, double alpha, double beta, double gamma) const; + + /** convert L3 module identifier to L1 or L2 */ + virtual Identifier getL1L2fromL3Identifier( const Identifier& ident + , const int& level + ) const ; + + /** get cumulative L1, L2, L3 trafo for (L3-) module. Result is in local frame. */ + virtual Amg::Transform3D getTransL123( const Identifier& ident ) const ; + + /** return value of particular transform specified by identifier and level + calculates L1 and L2 identifiers automatically by getL1L2fromL3Identifier + if L3 identifier passed. L1, L2 are in global, L3 in local frame. */ + virtual Amg::Transform3D getTrans( const Identifier& ident + , const int level + ) const ; + + // write the transforms to outputstream + virtual StatusCode outputObjs() const; + + // write the transform IOVs to IOVDB + virtual void fillDB(const std::string tag, + const unsigned int run1, const unsigned int event1, + const unsigned int run2, const unsigned int event2) const; + + // print the transforms, level=1 lists sizes, level=2 lists all transforms + virtual void printDB(const int level) const; + + // sort all the AlignableTransform objects so searches/inserts work properly + virtual void sortTrans() const; + + //calculate three rotations around locX,locY,locY = alpha,beta,gamma out of an HepGeom::Transform3D + void extractAlphaBetaGamma(const Amg::Transform3D & trans, + double& alpha, double& beta, double &gamma) const; + + private: + ServiceHandle < IToolSvc > p_toolsvc; + + const PixelID* m_pixid; + const SCT_ID* m_sctid; + const InDetDD::PixelDetectorManager* m_pixman; + const InDetDD::SCT_DetectorManager* m_sctman; + + std::vector<std::string> m_alignobjs; + std::vector<int> m_alignchans; + + bool par_newdb; // create database using new (collection) format + bool par_scttwoside; // create structures with separated SCT module sides + int par_fake; // set to 1 to fake full ATLAS geom, 2 to fake CTB geom + std::string par_condstream; // name of conditions stream + /** name of the root folder for constants, which can be set via + the <key> syntax. Default: /Indet/Align. */ + std::string par_dbroot; + /** the base part of the key for loading AlignableTransform objects + from the Transient Data Store. Default: /Indet/Align */ + std::string par_dbkey; + bool par_oldTextFile; // Input text file using old format + + AlignableTransform* getTransPtr(const std::string key) const; + const AlignableTransform* cgetTransPtr(const std::string key) const; + void fakeGeom(const int nbpix, const int necpix, + const int nbsct, const int necsct); + + +}; + +#endif // INDETALIGNGENTOOLS_ALIGNDBTOOL_H diff --git a/InnerDetector/InDetAlignTools/InDetAlignGenTools/InDetAlignGenTools/InDetAlignFillSiCluster.h b/InnerDetector/InDetAlignTools/InDetAlignGenTools/InDetAlignGenTools/InDetAlignFillSiCluster.h new file mode 100644 index 0000000000000000000000000000000000000000..1bf0b4fdbff3fa8bffd990b6efe677e60aae9c00 --- /dev/null +++ b/InnerDetector/InDetAlignTools/InDetAlignGenTools/InDetAlignGenTools/InDetAlignFillSiCluster.h @@ -0,0 +1,100 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +#ifndef INDETALIGNGENTOOLS_INDETALIGNFILLSICLUSTER_H +#define INDETALIGNGENTOOLS_INDETALIGNFILLSICLUSTER_H +// ================================================ +// InDetAlignFillSiCluster +// ================================================ +// +// InDetAlignFillSiCluster.h +// Headerfile for InDetAlignFillSiCluster +// +// Carlos Escobar, started 08/03/2008 +// +// AthAlgTool to create store Silicon Cluster information in a ntuple + +#include <string> + +#include "AthenaBaseComps/AthAlgTool.h" + +#include "InDetPrepRawData/PixelClusterContainer.h" +#include "InDetPrepRawData/SCT_ClusterContainer.h" + +#include "InDetAlignGenTools/IInDetAlignFillSiCluster.h" + +// Forward declaration +class PixelID; +class SCT_ID; + +class InDetAlignFillSiCluster: virtual public IInDetAlignFillSiCluster, public AthAlgTool { + + typedef InDet::PixelClusterContainer PixelClusterContainer; + typedef InDet::SCT_ClusterContainer SCT_ClusterContainer; + + public: + InDetAlignFillSiCluster(const std::string& type, + const std::string& name, + const IInterface* parent); + virtual ~InDetAlignFillSiCluster(); + + virtual StatusCode initialize(); + virtual StatusCode finalize(); + + virtual StatusCode FillSiCluster(); + + private: + INTupleSvc* m_ntupleSvc; + + const PixelID* m_pixelid; //!< Pixel ID helper + const SCT_ID* m_sctID; //!< SCT ID helper + + std::string m_Pixel_SiClustersName; + std::string m_Sct_SiClustersName; + + const InDet::PixelClusterContainer* m_Pixel_clcontainer; + const InDet::SCT_ClusterContainer* m_Sct_clcontainer; + + // methods + void bookNtuple(); + StatusCode RetrieveSCTSiClusters(); + StatusCode RetrievePixelSiClusters(); + void FillSCTSiNtuple(); + void FillPixelSiNtuple(); + + // variables + std::string m_ntupleName; //!< ntuple name + + NTuple::Item<long> m_pixel_nclusters; //!< number of Pixel Clusters + NTuple::Item<long> m_sct_nclusters; //!< number of SCT Clusters + + // ---------------------------------------------------------------------- + // Pixel Clusters + NTuple::Array<float> m_pixel_clx; //!< Pixel Cluster X + NTuple::Array<float> m_pixel_cly; //!< Pixel Cluster Y + NTuple::Array<float> m_pixel_clz; //!< Pixel Cluster Z + NTuple::Array<long> m_pixel_groupsize; //!< Pixel Cluster Group Size + NTuple::Array<long> m_pixel_layer; //!< Pixel Cluster layer + NTuple::Array<long> m_pixel_phi; //!< Pixel Cluster phi + NTuple::Array<float> m_pixel_LocX; //!< Pixel Cluster Local X + NTuple::Array<float> m_pixel_LocY; //!< Pixel Cluster Local Y + NTuple::Matrix<long> m_pixel_clrow; // list of pixel row + NTuple::Matrix<long> m_pixel_clcol; // list of pixel col + // ---------------------------------------------------------------------- + + // ---------------------------------------------------------------------- + // SCT Clusters + NTuple::Array<float> m_sct_clx; //!< SCT Cluster X + NTuple::Array<float> m_sct_cly; //!< SCT Cluster Y + NTuple::Array<float> m_sct_clz; //!< SCT Cluster Z + NTuple::Array<long> m_sct_groupsize; //!< SCT Cluster Group Size + NTuple::Array<long> m_sct_layer; //!< SCT Cluster layer + NTuple::Array<long> m_sct_eta; //!< SCT Cluster eta + NTuple::Array<long> m_sct_phi; //!< SCT Cluster phi + NTuple::Array<long> m_sct_side; //!< SCT Cluster side + // ---------------------------------------------------------------------- + +}; + +#endif // INDETALIGNGENTOOLS_INDETALIGNFILLSICLUSTER_H diff --git a/InnerDetector/InDetAlignTools/InDetAlignGenTools/InDetAlignGenTools/InDetAlignFillTrack.h b/InnerDetector/InDetAlignTools/InDetAlignGenTools/InDetAlignGenTools/InDetAlignFillTrack.h new file mode 100644 index 0000000000000000000000000000000000000000..036a9f2b3799d86608feb57e0d9de8e7041f1afc --- /dev/null +++ b/InnerDetector/InDetAlignTools/InDetAlignGenTools/InDetAlignGenTools/InDetAlignFillTrack.h @@ -0,0 +1,245 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +#ifndef INDETALIGNGENTOOLS_INDETALIGNFILLTRACK_H +#define INDETALIGNGENTOOLS_INDETALIGNFILLTRACK_H +// ================================================ +// InDetAlignFillTrack +// ================================================ +// +// InDetAlignFillTrack.h +// Headerfile for InDetAlignFillTrack +// +// Carlos Escobar, started 27/12/2007 +// +// AthAlgTool to create store Trk::Track and MCTrack information in a ntuple + +#include <string> + +#include "HepPDT/ParticleDataTable.hh" + +#include "TrkTrack/TrackCollection.h" + +#include "AthenaBaseComps/AthAlgTool.h" +#include "GaudiKernel/NTuple.h" +#include "GaudiKernel/ToolHandle.h" + +#include "InDetAlignGenTools/IInDetAlignFillTrack.h" +#include "TrkExInterfaces/IExtrapolator.h" + + +namespace Trk { + class Track; + class ITruthToTrack; //!< Produces perigee track parameters from generated parameters +} + +namespace HepMC { + class ParticleDataTable; +} + +class INTupleSvc; +class IExtrapolator; + +class InDetAlignFillTrack: virtual public IInDetAlignFillTrack, public AthAlgTool { + public: + InDetAlignFillTrack(const std::string& type, const std::string& name, + const IInterface* parent); + virtual ~InDetAlignFillTrack(); + + virtual StatusCode initialize(); + virtual StatusCode finalize(); + + virtual StatusCode FillTrack(); + + virtual int GetTrks() const { return m_totaltrks; } + virtual int GetTrkHits() const { return m_totalhits; } + virtual int GetTrkPixHits() const { return m_totalPixhits; } + virtual int GetTrkSCTHits() const { return m_totalSCThits; } + virtual int GetTrkTRTHits() const { return m_totalTRThits; } + + private: + INTupleSvc* m_ntupleSvc; + HepPDT::ParticleDataTable *m_mctable; + + // member functions + void bookNtuple(); + void bookUpNtuple(); + void bookLowNtuple(); + void bookMatchingNtuple(); + int dumpTrackCol(const TrackCollection*); + int dumpTrackCol(const TrackCollection*, const std::string); + void dumpTrack(int, const Trk::Track*, const std::string); + void dumpPerigee(const Trk::TrackParameters*, int); + StatusCode dumpMatching(const TrackCollection*,const TrackCollection*); + + // variables + bool m_doMatching; //!< switch on/off the matching information + bool m_doTruth; //!< switch on/off the truth information + + std::string m_inputCol; + std::string m_inputUpCol; + std::string m_inputLowCol; + std::string m_TruthTrkCol; + std::string m_ntupleName; //!< ntuple name + + int m_totaltrks; + int m_totalhits; + int m_totalPixhits; + int m_totalSCThits; + int m_totalTRThits; + int m_totalUptrks; + int m_totalUphits; + int m_totalUpPixhits; + int m_totalUpSCThits; + int m_totalUpTRThits; + int m_totalLowtrks; + int m_totalLowhits; + int m_totalLowPixhits; + int m_totalLowSCThits; + int m_totalLowTRThits; + int m_events; + + float m_matchedRcut; + float m_mindR; + + ToolHandle < Trk::ITruthToTrack > m_truthToTrack; //!< Pointer to TruthToTrack + ToolHandle < Trk::IExtrapolator > m_extrapolator; //!< Pointer to IExtrapolator + ToolHandle < Trk::ITrackSummaryTool > m_trackSumTool; //!< Pointer to Trk::ITrackSummaryTool + + NTuple::Item<long> nt_ntracks; //!< number of tracks + NTuple::Item<long> nt_nUptracks; //!< number of Up tracks + NTuple::Item<long> nt_nLowtracks; //!< number of Low tracks + NTuple::Item<long> nt_nmctracks; //!< number of mc tracks + NTuple::Item<long> nt_matchingTrk; //!< matching tracks + // ---------------------------------------------------------------------- + // Trk::Track parameters + NTuple::Array<float> nt_Trk_d0; //!< d0 parameter + NTuple::Array<float> nt_Trk_z0; //!< z0 parameter + NTuple::Array<float> nt_Trk_phi0; //!< phi0 parameter + NTuple::Array<float> nt_Trk_theta0; //!< theta0 parameter + NTuple::Array<float> nt_Trk_qoverp; //!< q/p parameter + NTuple::Array<float> nt_Trk_pt; //!< pt parameter + // ---------------------------------------------------------------------- + + // ---------------------------------------------------------------------- + NTuple::Array<long> nt_Trk_nHits; //!< number of hits + NTuple::Array<long> nt_Trk_nhitspix; //!< number of Pixel hits + NTuple::Array<long> nt_Trk_nhitssct; //!< number of SCT hits + NTuple::Array<long> nt_Trk_nhitstrt; //!< number of TRT hits + + NTuple::Array<long> nt_Trk_nsharedPixels; //!< number of Pixel shared hits + NTuple::Array<long> nt_Trk_nsharedSCT; //!< number of SCT shared hits + NTuple::Array<long> nt_Trk_nshared; //!< number of shared hits + + NTuple::Array<long> nt_Trk_nholesPixels; //!< number of Pixel holes + NTuple::Array<long> nt_Trk_nholesSCT; //!< number of SCT holes + NTuple::Array<long> nt_Trk_nholes; //!< number of holes + + NTuple::Array<float> nt_Trk_chi2; //!< number of chi2 + NTuple::Array<long> nt_Trk_ndof; //!< number of ndof + NTuple::Array<float> nt_Trk_chi2Prob; //!< number of chi2 probability + // ---------------------------------------------------------------------- + + // ---------------------------------------------------------------------- + // Trk::Track parameters + NTuple::Array<float> nt_Trk_d0_Up; //!< d0 parameter (Up track) + NTuple::Array<float> nt_Trk_z0_Up; //!< z0 parameter (Up track) + NTuple::Array<float> nt_Trk_phi0_Up; //!< phi0 parameter (Up track) + NTuple::Array<float> nt_Trk_theta0_Up; //!< theta0 parameter (Up track) + NTuple::Array<float> nt_Trk_qoverp_Up; //!< q/p parameter (Up track) + NTuple::Array<float> nt_Trk_pt_Up; //!< pt parameter (Up track) + // ---------------------------------------------------------------------- + + // ---------------------------------------------------------------------- + NTuple::Array<long> nt_Trk_nHits_Up; //!< number of hits (Up track) + NTuple::Array<long> nt_Trk_nhitspix_Up; //!< number of Pixel hits (Up track) + NTuple::Array<long> nt_Trk_nhitssct_Up; //!< number of SCT hits (Up track) + NTuple::Array<long> nt_Trk_nhitstrt_Up; //!< number of TRT hits (Up track) + + NTuple::Array<long> nt_Trk_nsharedPixels_Up; //!< number of Pixel shared hits (Up track) + NTuple::Array<long> nt_Trk_nsharedSCT_Up; //!< number of SCT shared hits (Up track) + NTuple::Array<long> nt_Trk_nshared_Up; //!< number of shared hits (Up track) + + NTuple::Array<long> nt_Trk_nholesPixels_Up; //!< number of Pixel holes (Up track) + NTuple::Array<long> nt_Trk_nholesSCT_Up; //!< number of SCT holes (Up track) + NTuple::Array<long> nt_Trk_nholes_Up; //!< number of holes (Up track) + + NTuple::Array<float> nt_Trk_chi2_Up; //!< number of chi2 (Up track) + NTuple::Array<long> nt_Trk_ndof_Up; //!< number of ndof (Up track) + NTuple::Array<float> nt_Trk_chi2Prob_Up; //!< number of chi2 probability (Up track) + // ---------------------------------------------------------------------- + + // ---------------------------------------------------------------------- + // Trk::Track parameters + NTuple::Array<float> nt_Trk_d0_Low; //!< d0 parameter (Low track) + NTuple::Array<float> nt_Trk_z0_Low; //!< z0 parameter (Low track) + NTuple::Array<float> nt_Trk_phi0_Low; //!< phi0 parameter (Low track) + NTuple::Array<float> nt_Trk_theta0_Low; //!< theta0 parameter (Low track) + NTuple::Array<float> nt_Trk_qoverp_Low; //!< q/p parameter (Low track) + NTuple::Array<float> nt_Trk_pt_Low; //!< pt parameter (Low track) + // ---------------------------------------------------------------------- + + // ---------------------------------------------------------------------- + NTuple::Array<long> nt_Trk_nHits_Low; //!< number of hits (Low track) + NTuple::Array<long> nt_Trk_nhitspix_Low; //!< number of Pixel hits (Low track) + NTuple::Array<long> nt_Trk_nhitssct_Low; //!< number of SCT hits (Low track) + NTuple::Array<long> nt_Trk_nhitstrt_Low; //!< number of TRT hits (Low track) + + NTuple::Array<long> nt_Trk_nsharedPixels_Low; //!< number of Pixel shared hits (Low track) + NTuple::Array<long> nt_Trk_nsharedSCT_Low; //!< number of SCT shared hits (Low track) + NTuple::Array<long> nt_Trk_nshared_Low; //!< number of shared hits (Low track) + + NTuple::Array<long> nt_Trk_nholesPixels_Low; //!< number of Pixel holes (Low track) + NTuple::Array<long> nt_Trk_nholesSCT_Low; //!< number of SCT holes (Low track) + NTuple::Array<long> nt_Trk_nholes_Low; //!< number of holes (Low track) + + NTuple::Array<float> nt_Trk_chi2_Low; //!< number of chi2 (Low track) + NTuple::Array<long> nt_Trk_ndof_Low; //!< number of ndof (Low track) + NTuple::Array<float> nt_Trk_chi2Prob_Low; //!< number of chi2 probability (Low track) + // ---------------------------------------------------------------------- + + // ---------------------------------------------------------------------- + NTuple::Array<float> nt_mc_trkistruth; //!< Has the Track an associated truth track? + + // generated particle parameters + NTuple::Array<float> nt_mc_Trk_genParticlePt; //!< generated pt + NTuple::Array<float> nt_mc_Trk_genParticleEta; //!< generated eta + NTuple::Array<float> nt_mc_Trk_genParticlePhi; //!< generated phi + + // MonteCarlo Track parameters + NTuple::Array<float> nt_mc_Trk_d0; //!< MonteCarlo d0 parameter + NTuple::Array<float> nt_mc_Trk_z0; //!< MonteCarlo z0 parameter + NTuple::Array<float> nt_mc_Trk_phi0; //!< MonteCarlo phi0 parameter + NTuple::Array<float> nt_mc_Trk_theta0; //!< MonteCarlo theta0 parameter + NTuple::Array<float> nt_mc_Trk_eta; //!< MonteCarlo eta parameter + NTuple::Array<float> nt_mc_Trk_qoverp; //!< MonteCarlo q/p parameter + NTuple::Array<float> nt_mc_Trk_qoverpt; //!< MonteCarlo q/pt parameter + + NTuple::Array<float> nt_mc_Trk_pt; //!< MonteCarlo pt parameter + NTuple::Array<float> nt_mc_Trk_charge; //!< MonteCarlo charge parameter + NTuple::Array<float> nt_mc_Trk_prob; //!< MonteCarlo prob parameter + NTuple::Array<float> nt_mc_Trk_pdg; //!< MonteCarlo pdg parameter + + NTuple::Array<float> nt_mc_Trk_vtxX; //!< MonteCarlo Vertex.X parameter + NTuple::Array<float> nt_mc_Trk_vtxY; //!< MonteCarlo Vertex.Y parameter + NTuple::Array<float> nt_mc_Trk_vtxZ; //!< MonteCarlo Vertex.Z parameter + // ---------------------------------------------------------------------- + + + // ---------------------------------------------------------------------- + // Matching Trk::Track_Up and Trk::Track_Low + NTuple::Array<float> nt_Trk_delta_d0; //!< d0 parameter + NTuple::Array<float> nt_Trk_delta_z0; //!< z0 parameter + NTuple::Array<float> nt_Trk_delta_phi0; //!< phi0 parameter + //** + NTuple::Array<float> nt_Trk_delta_theta0; //!< theta parameter + NTuple::Array<float> nt_Trk_delta_eta; //!< eta parameter + NTuple::Array<float> nt_Trk_delta_qoverpt; //!< q/pt parameter + NTuple::Array<float> nt_Trk_delta_pt; //!< pt parameter + NTuple::Array<float> nt_Trk_delta_charge; //!< charge parameter + // ---------------------------------------------------------------------- + +}; + +#endif // INDETALIGNGENTOOLS_INDETALIGNFILLTRACK_H diff --git a/InnerDetector/InDetAlignTools/InDetAlignGenTools/InDetAlignGenTools/InDetAlignHitQualSelTool.h b/InnerDetector/InDetAlignTools/InDetAlignGenTools/InDetAlignGenTools/InDetAlignHitQualSelTool.h new file mode 100644 index 0000000000000000000000000000000000000000..687b43f90f3cde22d4023f762f5940db34f0da2b --- /dev/null +++ b/InnerDetector/InDetAlignTools/InDetAlignGenTools/InDetAlignGenTools/InDetAlignHitQualSelTool.h @@ -0,0 +1,87 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +#ifndef INDETALIGNTOOLS_HITQUALSELTOOL_H +#define INDETALIGNTOOLS_HITQUALSELTOOL_H + +#include "AthenaBaseComps/AthAlgTool.h" +// #include "GaudiKernel/ToolHandle.h" +#include "InDetAlignGenTools/IInDetAlignHitQualSelTool.h" +// #include "AthenaBaseComps/AthMsgStreamMacros.h" + + +/** @class InDetAlignHitQualSelTool + @brief The InDetAlignHitQualSelTool is to select good quality hits to + build residuals for alignment. The following cuts are possible: + outlier hits, hits which have too large clusters, + hits with large incidence angles, ganged pixels, edge channels. + @author Oleg Brandt <http://consult.cern.ch/xwho> + */ + // Many thanks to Carlos Escobar and Sebastian Fleischmann! + +namespace Trk { + class TrackStateOnSurface ; + class RIO_OnTrack ; + class PrepRawData ; +} +namespace InDetDD { + class PixelDetectorManager ; + class SCT_DetectorManager ; +} +class AtlasDetectorID ; +class Identifier ; + +class InDetAlignHitQualSelTool : virtual public IInDetAlignHitQualSelTool, public AthAlgTool { + public: + InDetAlignHitQualSelTool( const std::string&, const std::string&, const IInterface* ) ; + virtual ~InDetAlignHitQualSelTool() ; + + virtual StatusCode initialize() ; + virtual StatusCode finalize () ; + + /** main method: from a TrackStateOnSurface select a good hit cutting on outlier + hits, hits with too many pixels/strips, hits with large incidence angles */ + const Trk::RIO_OnTrack* getGoodHit( const Trk::TrackStateOnSurface* tsos ) const ; + /** from a TrackStateOnSurface select a good hole in track cutting on + large incidence angles only */ + bool getGoodHole( const Trk::TrackStateOnSurface* tsos ) const ; + private: + // methods + /** check, whether cluster contains a ganged pixel */ + bool isGangedPixel( const Trk::PrepRawData* prd ) const ; + /** check, whether cluster size within limits of m_maxClusterSize */ + bool isGoodClusterSize( const std::vector<Identifier>& idVec ) const ; + /** check, whether the strip/pixel is an edge channel */ + bool isEdgeChannel( const std::vector<Identifier>& idVec ) const ; + /** calculate track incidence angle in local x-z frame */ + double incidAngle( const Trk::TrackParameters* trkPar + , const InDetDD::SiDetectorElement* detEle + ) const ; + /** check whether track incidence angle within limits of m_maxIncidAngle */ + bool isGoodAngle( const Trk::TrackParameters* trkPar + , const InDetDD::SiDetectorElement* detEle + ) const ; + // steering parameters + /** reject hits labeled as outliers by the track fitter */ + bool m_rejectOutliers ; + /* maximum size of the cluster forming the hit. typically, clusters with more than + 5 pixels/hits are due to brem and are of less quality for alignment. */ + int m_maxClusterSize ; + /** reject clusters containing edge channels */ + bool m_rejectEdgeChannels ; + /** reject clusters containing ganged pixels */ + bool m_rejectGangedPixels ; + /** maximum incidence angle of a track (to which the hit belongs) on the Si-module. + It is caculated in the local xz frame to ensure that the only the angle component + perpendicular to the strips is considered. It is defined w/r/t the local z-axis. */ + float m_maxIncidAngle ; + // "infrastructure" members + + const InDetDD::PixelDetectorManager* m_PIXManager ; //!< to get pixel phi and eta identifiers + const InDetDD::SCT_DetectorManager* m_SCTManager ; //!< to get strip numbers + const PixelID* m_pixelid ; //!< Pixel id helper + const SCT_ID* m_sctID ; //!< Pixel id helper +} ; + +#endif // INDETALIGNTOOLS_HITQUALSELTOOL_H diff --git a/InnerDetector/InDetAlignTools/InDetAlignGenTools/InDetAlignGenTools/InDetAlignOverlapTool.h b/InnerDetector/InDetAlignTools/InDetAlignGenTools/InDetAlignGenTools/InDetAlignOverlapTool.h new file mode 100644 index 0000000000000000000000000000000000000000..c9b10f13c09db2d70f92ee7047b2aa259d984b51 --- /dev/null +++ b/InnerDetector/InDetAlignTools/InDetAlignGenTools/InDetAlignGenTools/InDetAlignOverlapTool.h @@ -0,0 +1,63 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +#ifndef INDETALIGNTOOLS_OVERLAPTOOL_H +#define INDETALIGNTOOLS_OVERLAPTOOL_H + +#include "AthenaBaseComps/AthAlgTool.h" +#include "InDetAlignGenTools/IInDetAlignOverlapTool.h" + +#include "InDetAlignGenTools/IInDetAlignFillTrack.h" +#include "InDetAlignGenTools/AlignSiModuleList.h" +#include <vector> +//#include "DataModel/DataVector.h" + +/** @class InDetAlignOverlapTool + @brief The InDetAlignOverlapTool is to select good quality tracks for alignment to + build residuals with possible cuts on overlap hits, overlaps in SCT and PIXELs (in Phi and Eta). +*/ + +class AlignTrk; +class AlignSiHit; +class AlignSiModuleList; + + + +class InDetAlignOverlapTool : virtual public IInDetAlignOverlapTool, public AthAlgTool { + public: + //standard constructor + InDetAlignOverlapTool( const std::string&, const std::string&, const IInterface* ) ; + virtual ~InDetAlignOverlapTool() ; + + virtual StatusCode initialize() ; + virtual StatusCode finalize () ; + + /**main method to get an overlap hit **/ + + int getNumberOverlapPIX( const AlignTrk& ) const ; + int getNumberOverlapSCT( const AlignTrk& ) const ; + // bool getOverlapHit( const AlignSiHit&, const AlignTrk& ) const ; + //AlignSiHit getOverlapHit( const AlignTrk& ) const ; + //const DataVector<AlignSiHit>* getOverlapHit( const AlignTrk& ) const ; + //std::vector<AlignSiHit> getOverlapHit( const AlignTrk& ) ;//const ; + std::vector<AlignSiHit> getOverlapHit( const AlignTrk& );// const ; + + + private: + //ServiceHandle < IToolSvc > p_toolsvc; + //StatusCode FindOverlap(const AlignTrk& ); + //AlignSiHit* Overlaphits; + //DataVector<AlignSiHit>* m_Overlaphits; + + std::vector<AlignSiHit> m_Overlaphits; + + int m_minPixelOverlapPerTrack; + int m_minSCTOverlapPerTrack; + + const AlignSiModuleList* p_modlist; + + +} ; + +#endif // INDETALIGNTOOLS_OVERLAPTOOL_H diff --git a/InnerDetector/InDetAlignTools/InDetAlignGenTools/InDetAlignGenTools/InDetAlignTrackSelTool.h b/InnerDetector/InDetAlignTools/InDetAlignGenTools/InDetAlignGenTools/InDetAlignTrackSelTool.h new file mode 100755 index 0000000000000000000000000000000000000000..c335920cc3bd62575f905f67b15761d87088cd48 --- /dev/null +++ b/InnerDetector/InDetAlignTools/InDetAlignGenTools/InDetAlignGenTools/InDetAlignTrackSelTool.h @@ -0,0 +1,64 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +// InDetAlignTrackSelTool.h + +// Sergio Gonzalez Sevilla, started 04/7/05 +// Miguel Olivo Gomez, extended 07/6/06 + +// AlgTool to select high quality tracks for the inner detector +// (Pixel+SCT) alignment algorithms. +// This AlgTool provides a track selection based on the following cut variables: +// Momentum, pt, number of shared hits, number of holes and chi2 probability. +// Returns 0 in case a track is not accepted, otherwise 1 + +/** @class InDetAlignTrackSelTool + @brief The InDetAlignTrackSelTool serves to select high quality tracks for the inner detector + (Pixel+SCT) alignment algorithms. This AlgTool provides a track selection based on the + following cut variables: + Momentum, pt, number of shared hits, number of holes and chi2 probability. + Returns 0 in case a track is not accepted, otherwise 1 + @author Sergio Gonzalez Sevilla, Miguel Olivo Gomez <http://consult.cern.ch/xwho> + */ + +#ifndef INDETALIGNGENTOOLS_ALIGNTRACKSELTOOL_H +#define INDETALIGNGENTOOLS_ALIGNTRACKSELTOOL_H + +#include "AthenaBaseComps/AthAlgTool.h" +#include "GaudiKernel/ToolHandle.h" +#include "TrkToolInterfaces/ITrackSummaryTool.h" +#include "InDetAlignGenTools/IInDetAlignTrackSelTool.h" + +namespace Trk { + class Track; + class ITrackSummaryTool; +} + +class InDetAlignTrackSelTool: virtual public IInDetAlignTrackSelTool, public AthAlgTool { + public: + InDetAlignTrackSelTool(const std::string& type, const std::string& name, + const IInterface* parent); + virtual ~InDetAlignTrackSelTool(); + virtual StatusCode initialize(); + virtual StatusCode finalize(); + + virtual int getStatus(const Trk::Track&) const; + + private: + double Momentum(const Trk::Track&) const; + double Pt(const Trk::Track&) const; + int nShared(const Trk::Track&) const; + int nHoles(const Trk::Track&) const; + double chi2Prob(const Trk::Track&) const; + + ToolHandle < Trk::ITrackSummaryTool > m_trackSumTool; //!< Pointer to Trk::ITrackSummaryTool + + double m_minMomentum; //!< Minimum value of the momentum of the track + double m_minPt; //!< Minimum value of the transverse momentum of the track + int m_maxShared; //!< Maximum number of shared hits of the track + int m_maxHoles; //!< Maximum number of holes of the track + double m_minChi2Prob; //!< Minimum value of the chi2 Probality of the track +}; + +#endif // INDETALIGNTOOLS_ALIGNTRACKSELTOOL_H diff --git a/InnerDetector/InDetAlignTools/InDetAlignGenTools/InDetAlignGenTools/RefitSiOnlyTool.h b/InnerDetector/InDetAlignTools/InDetAlignGenTools/InDetAlignGenTools/RefitSiOnlyTool.h new file mode 100755 index 0000000000000000000000000000000000000000..d0081c4504d1a58564e35a70ac54bd1c7db2d474 --- /dev/null +++ b/InnerDetector/InDetAlignTools/InDetAlignGenTools/InDetAlignGenTools/RefitSiOnlyTool.h @@ -0,0 +1,82 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +#ifndef INDETALIGNGENTOOLS_REFITSIONLYTOOL_H +#define INDETALIGNGENTOOLS_REFITSIONLYTOOL_H + +////////////////////////////////////////////////////////////////////////// +// ================================================ +// RefitSiOnly +// ================================================ +// +// RefitSiOnlyTool.h +// Header file for RefitSiOnlyTool +// +// Namespace SiRobustAlign +// +// AlgTool to create refitted tracks with only silicon hits present + + +#include "AthenaBaseComps/AthAlgTool.h" +#include "GaudiKernel/ToolHandle.h" +#include "InDetAlignToolInterfaces/ICollectionProcessor.h" + +#include "TrkEventPrimitives/ParticleHypothesis.h" +#include "TrkEventUtils/TrkParametersComparisonFunction.h" + +namespace Trk { + class Track; + class ITrackFitter; +} + +class AtlasDetectorID; + +namespace InDetAlignment { + class RefitSiOnlyTool : virtual public InDetAlignment::ICollectionProcessor, public AthAlgTool { + public: + RefitSiOnlyTool (const std::string& type, const std::string& name, const IInterface* parent); //constructor + virtual ~RefitSiOnlyTool(); //destructor + + virtual StatusCode initialize(); ///< Athena Tool methods + virtual StatusCode finalize(); ///< Athena Tool methods + + /** does necessary configuration to the tool */ + //virtual void configure(); + virtual StatusCode configure(); + /** main processing of track collection - + refits a track using only filtered hits (e.g. only silicon hits) */ + virtual const DataVector<Trk::Track>* processTrackCollection(const DataVector<Trk::Track>* trks); + const AtlasDetectorID * m_idHelper; + + void dumpTrackCol(DataVector<Trk::Track>* tracks); + + private: + ///refits a track not using non-silicon hits + Trk::Track* refit(const Trk::Track* tr) const; + + ///determines from flags, if this hit should be included in the fit or not + bool filterHit(const Identifier& id) const; + + + ToolHandle<Trk::ITrackFitter> m_ITrkFitter; + + Trk::TrkParametersComparisonFunction* m_comPar; + + bool m_doCosmic; //!< Flag for cosmics + bool m_OutlierRemoval; //!< flag for outlier removal + + bool m_stripSCTHits; ///< flag to not consider SCT hits in the refit + bool m_stripPixelHits; ///< flag to not consider Pixel hits in the refit + bool m_stripTRTHits; ///< flag to not consider TRT hits in the refit + + int m_ParticleNumber; + Trk::ParticleHypothesis m_ParticleHypothesis; + int m_trkindex; + + + + }; +} //namespace + +#endif // INDETALIGNGENTOOLS_REFITSIONLYTOOL_H diff --git a/InnerDetector/InDetAlignTools/InDetAlignGenTools/cmt/requirements b/InnerDetector/InDetAlignTools/InDetAlignGenTools/cmt/requirements new file mode 100755 index 0000000000000000000000000000000000000000..7824b52874e03acc12e99c1e33346fdff7a61b6f --- /dev/null +++ b/InnerDetector/InDetAlignTools/InDetAlignGenTools/cmt/requirements @@ -0,0 +1,55 @@ +package InDetAlignGenTools + +author Richard Hawkings <Richard.Hawkings@cern.ch> +author Sergio Gonzalez Sevilla <segonzal@ific.uv.es> +author Carlos Escobar <cescobar@ific.uv.es> + +use AtlasPolicy AtlasPolicy-* +use GaudiInterface GaudiInterface-* External +use AthenaBaseComps AthenaBaseComps-* Control +use HepPDT v* LCG_Interfaces +use InDetPrepRawData InDetPrepRawData-* InnerDetector/InDetRecEvent +use InDetAlignToolInterfaces InDetAlignToolInterfaces-* InnerDetector/InDetAlignTools +use TrkToolInterfaces TrkToolInterfaces-* Tracking/TrkTools +use TrkTrack TrkTrack-* Tracking/TrkEvent +use TrkEventPrimitives TrkEventPrimitives-* Tracking/TrkEvent +use TrkEventUtils TrkEventUtils-* Tracking/TrkEvent +use TrkExInterfaces TrkExInterfaces-* Tracking/TrkExtrapolation +use InDetAlignTrkInfo InDetAlignTrkInfo-* InnerDetector/InDetAlignEvent +use CLIDSvc CLIDSvc-* Control +#use InDetAlignGenAlgs InDetAlignGenAlgs-* InnerDetector/InDetAlignAlgs +use Identifier Identifier-* DetectorDescription +use EventPrimitives EventPrimitives-* Event +use GeoPrimitives GeoPrimitives-* DetectorDescription + +private +use AtlasHepMC AtlasHepMC-* External +use GaudiInterface GaudiInterface-* External +use AthenaKernel AthenaKernel-* Control +use DataModel DataModel-* Control +use RegistrationServices RegistrationServices-* Database +use AtlasDetDescr AtlasDetDescr-* DetectorDescription +use DetDescrConditions DetDescrConditions-* DetectorDescription/DetDescrCond +#use InDetPriVxFinderTool InDetPriVxFinderTool-* InnerDetector/InDetRecTools +#use InDetRecToolInterfaces InDetRecToolInterfaces-* InnerDetector/InDetRecTools +use InDetIdentifier InDetIdentifier-* InnerDetector/InDetDetDescr +use InDetReadoutGeometry InDetReadoutGeometry-* InnerDetector/InDetDetDescr +use TrackRecord TrackRecord-* Simulation/G4Sim +use TrkTruthData TrkTruthData-* Tracking/TrkEvent +#use VxVertex VxVertex-* Tracking/TrkEvent +#use TrkVertexOnTrack TrkVertexOnTrack-* Tracking/TrkEvent +use TrkTrackSummary TrkTrackSummary-* Tracking/TrkEvent +use TrkMeasurementBase TrkMeasurementBase-* Tracking/TrkEvent +use TrkParameters TrkParameters-* Tracking/TrkEvent +use TrkPrepRawData TrkPrepRawData-* Tracking/TrkEvent +use TrkRIO_OnTrack TrkRIO_OnTrack-* Tracking/TrkEvent +#use TrkTruthToTrack TrkTruthToTrack-* Tracking/TrkTools +#use TrkVertexFitterInterfaces TrkVertexFitterInterfaces-* Tracking/TrkVertexFitter +use TrkFitterInterfaces TrkFitterInterfaces-* Tracking/TrkFitter + +use AtlasCLHEP AtlasCLHEP-* External + + +public +library InDetAlignGenTools *.cxx components/*.cxx +apply_pattern component_library diff --git a/InnerDetector/InDetAlignTools/InDetAlignGenTools/doc/mainpage.h b/InnerDetector/InDetAlignTools/InDetAlignGenTools/doc/mainpage.h new file mode 100755 index 0000000000000000000000000000000000000000..235e459fbfb32a016a0a279971ab3b633491ad39 --- /dev/null +++ b/InnerDetector/InDetAlignTools/InDetAlignGenTools/doc/mainpage.h @@ -0,0 +1,734 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +/** +\mainpage The InDetAlignGenTools package + +\section introductionInDetAlignGenTools Introduction + +This package provides a place to put AlgTools related to ID alignment. +The available AlgTools (so far InDetAlignDBTool and InDetAlignTrackSelTool) +are discussed below. + +\section InDetAlignDBTool + +The InDetAlignDBTool AlgTool provides tools (methods) to manipulate the +set of AlignableTransform objects which make up the ID alignment database. +For more details on the object format, see documentation for the InDetAlignAlgs +package. Note that the format changed at release 10.5.0 when COOL was +introduced and the writing of seperate transforms for each side of SCT +modules was suppressed. + +On initialisation, the InDetAlignDBTool builds a list of AlignableTransform +names appropriate to the current geometry (e.g. CTB or full ATLAS, with or +without middle pixel layers). It provides the following methods (see the +header file for detailed parameter list). + + - createDB: Create a null database in memory, with all modules at default + positions. + + - dispGroup: displace a group of modules (specified by detector and layer), + either randomly or systematically. + + - printDB: Print the contents of the database to the screen (transform IDs + and displacements). + + - writeFile: Write the contents of the database to a text file or ntuple. + + - readTextFile: Read back the database from a text file, and fill the + corresponding structures in the TDS + + - readNtuple: Read back the database from an ntuple (9002), and fill the + corresondping structures in the TDS. + + - dirkey: return the leafname of the AlignableTransform corresponding to + a particular SCT or pixel identifier (e.g. SCTB3 or PIXEA1). + + - setTrans: set the transformation for a particular module to a given value. + + - tweakTrans: modify (i.e. add to existing transform) the transformation + for a particular module, allowing incremental changes to the alignment. + + - get Trans: return the value of the transform for a particular module and + level. + + - outputObjs: Stream out the alignment database objects from the TDS + to a output stream, using the AthenaOutputStreamTool. + + - fillDB: Register the alignment database objects to the conditionsDB, using + the FillIOVTool provided by the database group. Note that the objects + must already have been streamed out using outputObjs. + +When reading/writing to an ntuple, the ntuple is numbered 9002. The format +is similar to the old ntuple 9000 alignment ntuples, but the transformation +'level' has been added, allowing the global transforms for the whole +detector or detector layers to be written and read back too. Unlike ntuple +9000, the r-phi and stereo sides of SCT modules are not included as separate +modules. + +Note that when reading back an alignment database through the tool (e.g. from +a text file or ntuple, the database will be +loaded into memory, but there is currently no method to tell the ID detector +descrption to 'pick up' the new alignment constants (no IOVSvc +call backs available). This can only be done when the database is stored +as objects in a POOL file. The readback can then be triggered either +by loading the POOL file directly using the CondProxyProviderSvc, or by +loading it via entries in the conditions database. See the joboptions +ReadDBS.py and ReadPool.py in InDetAlignAlgs for examples of how to do this. + +The tool has the following properties: + + - NewDB: Flag to use new COOL DB format (default: TRUE) + - FakeDB: Flag allowing geometry information to be faked if GeoModel + information is not available. This allows limited database manipulation to + be performed in an Athena job without GeoModel (database I/O but not + creation/modification of databases). FakeDB may be set to 1 (full ATLAS + geometry with 3 pixel layers) or 2 (combined testbeam geometry). + - CondStream: Name of conditions output stream on which outputObjs will + write conditions data (default: "Condstream1") + - DBRoot: Name of the conditions database directory and Storegate container + where the alignment data is stored (default: "/Indet/Align"). + - properties related to the DispCSC method (see below). + +The algorithm InDetAlignWrt in InDetAlignAlgs provides some examples of how +to use the methods in this AlgTool. + +The DBRoot parameter can be changed to allow the tool to manipulate a different +instance of the alignment constants in the TDS and conditions database. +This can be used for example to store the SCT/pixel survey data +(also consisting of AlignableTransforms for every module) in the conditions +database in the folder /Indet/SiSurvey. For more details of this, see +the documentation of the InDetAlignWrt algorithm in InDetAlignGenAlgs. + +\section DispCSC + +The DispCSC method inside InDetAlignDBTool allows to generate the misalignments +of the Pixel and SCT detectors for the CSC exercise. The misalignments are generated at +three different levels, where, depending on which, the reference frame for the transformation +to be applied may change. +These levels are quoted below: +- LEVEL 3: module level (local frame) +- LEVEL 2: layers/disks level (global frame) +- LEVEL 1: subdetector level (global frame) + +For more details you may take a look on the wiki topic SiliconMisaCSC and related links inside. + +This DispCSC method is triggered while setting the property: +<pre> +InDetAlignWrt.DispCSC +</pre> +of the InDetAlignWrt algorithm in InDetAlignGenAlgs to the desired level of misalignment. +A combination of different levels is valid. For example, if we want to generate misalignments +for all silicon modules (level 3) and silicon subdetectors (level 1), without any +misalignment for silicon layers nor disks (level 2), we should set: +<pre> +InDetAlignWrt.DispCSC=31 +</pre> +For the CSC sets, misalignments at all different levels need to be generated. Thus, it +remains mandatory to set: +\verbatim +InDetAlignWrt.DispCSC=321 +\endverbatim + +The properties the user can set are almost identical for all the three levels. +The misalignment constants can be generated in four different ways: +- from fixed values +- from a flat distribution +- from a gaussian distribution +- from a truncated gaussian distribution + +The properties +<pre> +DistTypePixelBaModules +DistTypePixelEcModules +DistTypeSCTBaModules +DistTypeSCTEcModules +DistTypePixelLayers +DistTypeSCTLayers +DistTypePixelDisks +DistTypeSCTDisks +DistTypeID +</pre> +allow to select between the different ways explained above and have an effect for +Pixel barrel modules (DistTypePixelBaModules), Pixel endcap modules (DistTypePixelEcModules), +SCT barrel modules (DistTypeSCTBaModules), SCT endcap modules (DistTypeSCTEcModules), +Pixel barrel layers (DistTypePixelLayers), SCT barrel layers (DistTypeSCTLayers), +Pixel endcap disks (DistTypePixelDisks), SCT endcap disks (DistTypeSCTDisks) and +ID subdetectors (DistTypeID). + +The only allowed values for these properties is an integer number in the range [0;3] +where each number means: +<pre> +0: fixed numbers +1: flat distribution +2: gaussian distribution +3: truncated gaussian distribution +</pre> + +If other value different from the above is set, StatusCode::FAILURE is returned. + +If the numbers are chosen to be generated from a distribution, then additional properties may be set, +as the mean and the sigma of the given distribution. The distribution characteristics can be specified +for all 6 variables considered for the transformations leading up to the desired misalignment: 3 translations and +3 rotations (defined by three angles: alpha, beta, gamma). + +As an example for the pixel endcap modules, if the following properties are set as: +<pre> +DistTypePixelEcModules(2) + +Pix_Mod_Ec_Mean_x = 10 +Pix_Mod_Ec_Mean_y = 10 +Pix_Mod_Ec_Mean_z = 30 +Pix_Mod_Ec_Mean_a = 0 +Pix_Mod_Ec_Mean_b = 0 +Pix_Mod_Ec_Mean_g = 0 + +Pix_Mod_Ec_Sigma_x = 0.020 +Pix_Mod_Ec_Sigma_y = 0.020 +Pix_Mod_Ec_Sigma_z = 0.060 +Pix_Mod_Ec_Sigma_a = 0.001 +Pix_Mod_Ec_Sigma_b = 0.001 +Pix_Mod_Ec_Sigma_g = 0.001 +</pre> +the in-plane translations will be generated according to a gaussian distribution centered in 10 microns and +with sigma 20 microns, the out-of plane translation according to a gaussian centered in 30 microns and +with sigma 60 microns, and the three rotations according to gaussians centered in zero and with sigma equal +to 1 mrad. + +If a truncated gaussian is chose, additional variables may be set: +<pre> +TruncGaussCut_X +TruncGaussCut_Y +TruncGaussCut_Z +TruncGaussCut_A +TruncGaussCut_B +TruncGaussCut_G +</pre> +which specify, for each of the different movements, +the limit value times the sigma out-of which the generated numbers are rejected. The default values are: +\verbatim +InDetDBTool.TruncGaussCut_X(2.) +InDetDBTool.TruncGaussCut_Y(2.) +InDetDBTool.TruncGaussCut_Z(2.) +InDetDBTool.TruncGaussCut_A(2.) +InDetDBTool.TruncGaussCut_B(2.) +InDetDBTool.TruncGaussCut_G(2.) +\endverbatim +which means that the generated numbers that do not remain in the range [-2*sigma; 2*sigma] are rejected. + +For completeness, all available properties including their default values +for all three levels are shown below. + +All displacements are given in mm and rotations in mrad. + +Note the different suffix for level 3 misalignments (lowercase) with respect level 2 and 1 +(uppercase). This is intended for differentiating the local to the global frame. + +\subsection subsection1 LEVEL 3 + -# Pixel barrel modules: + \verbatim + InDetDBTool.DistTypePixelBaModules(1) - distribution type, Pixel barrel modules (in this case: flat) + + InDetDBTool.Pix_Mod_Ba_Mean_x(0) - mean value of distribution for x translation, Pixel barrel modules + InDetDBTool.Pix_Mod_Ba_Mean_y(0) - mean value of distribution for y translation, Pixel barrel modules + InDetDBTool.Pix_Mod_Ba_Mean_z(0) - mean value of distribution for z translation, Pixel barrel modules + InDetDBTool.Pix_Mod_Ba_Mean_a(0) - mean value of distribution for alpha rotation, Pixel barrel modules + InDetDBTool.Pix_Mod_Ba_Mean_b(0) - mean value of distribution for beta rotation, Pixel barrel modules + InDetDBTool.Pix_Mod_Ba_Mean_g(0) - mean value of distribution for gamma rotation, Pixel barrel modules + + InDetDBTool.Pix_Mod_Ba_Sigma_x(0.030) - sigma value of distribution for x translation, Pixel barrel modules + InDetDBTool.Pix_Mod_Ba_Sigma_y(0.030) - sigma value of distribution for y translation, Pixel barrel modules + InDetDBTool.Pix_Mod_Ba_Sigma_z(0.050) - sigma value of distribution for z translation, Pixel barrel modules + InDetDBTool.Pix_Mod_Ba_Sigma_a(0.001) - sigma value of distribution for alpha rotation, Pixel barrel modules + InDetDBTool.Pix_Mod_Ba_Sigma_b(0.001) - sigma value of distribution for beta rotation, Pixel barrel modules + InDetDBTool.Pix_Mod_Ba_Sigma_g(0.001) - sigma value of distribution for gamma rotation, Pixel barrel modules + + InDetDBTool.Pix_Mod_Ba_Fix_x(0) - fixed value for x translation, Pixel barrel modules + InDetDBTool.Pix_Mod_Ba_Fix_y(0) - fixed value for y translation, Pixel barrel modules + InDetDBTool.Pix_Mod_Ba_Fix_z(0) - fixed value for z translation, Pixel barrel modules + InDetDBTool.Pix_Mod_Ba_Fix_a(0) - fixed value for alpha rotation, Pixel barrel modules + InDetDBTool.Pix_Mod_Ba_Fix_b(0) - fixed value for beta rotation, Pixel barrel modules + InDetDBTool.Pix_Mod_Ba_Fix_g(0) - fixed value for gamma rotation, Pixel barrel modules + \endverbatim + -# Pixel endcap modules: + \verbatim + InDetDBTool.DistTypePixelEcModules(1) - distribution type, Pixel endcap modules (in this case: flat) + + InDetDBTool.Pix_Mod_Ec_Mean_x(0) - mean value of distribution for x translation, Pixel endcap modules + InDetDBTool.Pix_Mod_Ec_Mean_y(0) - mean value of distribution for y translation, Pixel endcap modules + InDetDBTool.Pix_Mod_Ec_Mean_z(0) - mean value of distribution for z translation, Pixel endcap modules + InDetDBTool.Pix_Mod_Ec_Mean_a(0) - mean value of distribution for alpha rotation, Pixel endcap modules + InDetDBTool.Pix_Mod_Ec_Mean_b(0) - mean value of distribution for beta rotation, Pixel endcap modules + InDetDBTool.Pix_Mod_Ec_Mean_g(0) - mean value of distribution for gamma rotation, Pixel endcap modules + + InDetDBTool.Pix_Mod_Ec_Sigma_x(0.030) - sigma value of distribution for x translation, Pixel endcap modules + InDetDBTool.Pix_Mod_Ec_Sigma_y(0.030) - sigma value of distribution for y translation, Pixel endcap modules + InDetDBTool.Pix_Mod_Ec_Sigma_z(0.050) - sigma value of distribution for z translation, Pixel endcap modules + InDetDBTool.Pix_Mod_Ec_Sigma_a(0.001) - sigma value of distribution for alpha rotation, Pixel endcap modules + InDetDBTool.Pix_Mod_Ec_Sigma_b(0.001) - sigma value of distribution for beta rotation, Pixel endcap modules + InDetDBTool.Pix_Mod_Ec_Sigma_g(0.001) - sigma value of distribution for gamma rotation, Pixel endcap modules + + InDetDBTool.Pix_Mod_Ec_Fix_x(0) - fixed value for x translation, Pixel endcap modules + InDetDBTool.Pix_Mod_Ec_Fix_y(0) - fixed value for y translation, Pixel endcap modules + InDetDBTool.Pix_Mod_Ec_Fix_z(0) - fixed value for z translation, Pixel endcap modules + InDetDBTool.Pix_Mod_Ec_Fix_a(0) - fixed value for alpha rotation, Pixel endcap modules + InDetDBTool.Pix_Mod_Ec_Fix_b(0) - fixed value for beta rotation, Pixel endcap modules + InDetDBTool.Pix_Mod_Ec_Fix_g(0) - fixed value for gamma rotation, Pixel endcap modules + \endverbatim + -# SCT barrel modules: + \verbatim + InDetDBTool.DistTypeSCTBaModules(1) - distribution type, SCT barrel modules (in this case: flat) + + InDetDBTool.SCT_Mod_Ba_Mean_x(0) - mean value of distribution for x translation, SCT barrel modules + InDetDBTool.SCT_Mod_Ba_Mean_y(0) - mean value of distribution for y translation, SCT barrel modules + InDetDBTool.SCT_Mod_Ba_Mean_z(0) - mean value of distribution for z translation, SCT barrel modules + InDetDBTool.SCT_Mod_Ba_Mean_a(0) - mean value of distribution for alpha rotation, SCT barrel modules + InDetDBTool.SCT_Mod_Ba_Mean_b(0) - mean value of distribution for beta rotation, SCT barrel modules + InDetDBTool.SCT_Mod_Ba_Mean_g(0) - mean value of distribution for gamma rotation, SCT barrel modules + + InDetDBTool.SCT_Mod_Ba_Sigma_x(0.150) - sigma value of distribution for x translation, SCT barrel modules + InDetDBTool.SCT_Mod_Ba_Sigma_y(0.150) - sigma value of distribution for y translation, SCT barrel modules + InDetDBTool.SCT_Mod_Ba_Sigma_z(0.150) - sigma value of distribution for z translation, SCT barrel modules + InDetDBTool.SCT_Mod_Ba_Sigma_a(0.001) - sigma value of distribution for alpha rotation, SCT barrel modules + InDetDBTool.SCT_Mod_Ba_Sigma_b(0.001) - sigma value of distribution for beta rotation, SCT barrel modules + InDetDBTool.SCT_Mod_Ba_Sigma_g(0.001) - sigma value of distribution for gamma rotation, SCT barrel modules + + InDetDBTool.SCT_Mod_Ba_Fix_x(0) - fixed value for x translation, SCT barrel modules + InDetDBTool.SCT_Mod_Ba_Fix_y(0) - fixed value for y translation, SCT barrel modules + InDetDBTool.SCT_Mod_Ba_Fix_z(0) - fixed value for z translation, SCT barrel modules + InDetDBTool.SCT_Mod_Ba_Fix_a(0) - fixed value for alpha rotation, SCT barrel modules + InDetDBTool.SCT_Mod_Ba_Fix_b(0) - fixed value for beta rotation, SCT barrel modules + InDetDBTool.SCT_Mod_Ba_Fix_g(0) - fixed value for gamma rotation, SCT barrel modules + \endverbatim + -# SCT endcap modules: + \verbatim + InDetDBTool.DistTypeSCTEcModules(1) - distribution type, SCT endcap modules (in this case: flat) + + InDetDBTool.SCT_Mod_Ec_Mean_x(0) - mean value of distribution for x translation, SCT endcap modules + InDetDBTool.SCT_Mod_Ec_Mean_y(0) - mean value of distribution for y translation, SCT endcap modules + InDetDBTool.SCT_Mod_Ec_Mean_z(0) - mean value of distribution for z translation, SCT endcap modules + InDetDBTool.SCT_Mod_Ec_Mean_a(0) - mean value of distribution for alpha rotation, SCT endcap modules + InDetDBTool.SCT_Mod_Ec_Mean_b(0) - mean value of distribution for beta rotation, SCT endcap modules + InDetDBTool.SCT_Mod_Ec_Mean_g(0) - mean value of distribution for gamma rotation, SCT endcap modules + + InDetDBTool.SCT_Mod_Ec_Sigma_x(0.100) - sigma value of distribution for x translation, SCT endcap modules + InDetDBTool.SCT_Mod_Ec_Sigma_y(0.100) - sigma value of distribution for y translation, SCT endcap modules + InDetDBTool.SCT_Mod_Ec_Sigma_z(0.150) - sigma value of distribution for z translation, SCT endcap modules + InDetDBTool.SCT_Mod_Ec_Sigma_a(0.001) - sigma value of distribution for alpha rotation, SCT endcap modules + InDetDBTool.SCT_Mod_Ec_Sigma_b(0.001) - sigma value of distribution for beta rotation, SCT endcap modules + InDetDBTool.SCT_Mod_Ec_Sigma_g(0.001) - sigma value of distribution for gamma rotation, SCT endcap modules + + InDetDBTool.SCT_Mod_Ec_Fix_x(0) - fixed value for x translation, SCT endcap modules + InDetDBTool.SCT_Mod_Ec_Fix_y(0) - fixed value for y translation, SCT endcap modules + InDetDBTool.SCT_Mod_Ec_Fix_z(0) - fixed value for z translation, SCT endcap modules + InDetDBTool.SCT_Mod_Ec_Fix_a(0) - fixed value for alpha rotation, SCT endcap modules + InDetDBTool.SCT_Mod_Ec_Fix_b(0) - fixed value for beta rotation, SCT endcap modules + InDetDBTool.SCT_Mod_Ec_Fix_g(0) - fixed value for gamma rotation, SCT endcap modules + \endverbatim + +\subsection subsection2 LEVEL 2 + -# Pixel Barrel layers + \verbatim + InDetDBTool.DistTypePixelLayers(0) - distribution type, Pixel barrel layers (in this case: fixed values) + + InDetDBTool.Pix_Layers_Mean_X(0) - mean value of distribution for X translation, Pixel barrel layers + InDetDBTool.Pix_Layers_Mean_Y(0) - mean value of distribution for Y translation, Pixel barrel layers + InDetDBTool.Pix_Layers_Mean_Z(0) - mean value of distribution for Z translation, Pixel barrel layers + InDetDBTool.Pix_Layers_Mean_A(0) - sigma value of distribution for Alpha rotation, Pixel barrel layers + InDetDBTool.Pix_Layers_Mean_B(0) - mean value of distribution for Beta rotation, Pixel barrel layers + InDetDBTool.Pix_Layers_Mean_G(0) - mean value of distribution for Gamma rotation, Pixel barrel layers + + InDetDBTool.Pix_Layers_Sigma_X(0) - sigma value of distribution for X translation, Pixel barrel layers + InDetDBTool.Pix_Layers_Sigma_Y(0) - sigma value of distribution for Y translation, Pixel barrel layers + InDetDBTool.Pix_Layers_Sigma_Z(0) - sigma value of distribution for Z translation, Pixel barrel layers + InDetDBTool.Pix_Layers_Sigma_A(0) - sigma value of distribution for Alpha rotation, Pixel barrel layers + InDetDBTool.Pix_Layers_Sigma_B(0) - sigma value of distribution for Beta rotation, Pixel barrel layers + InDetDBTool.Pix_Layers_Sigma_G(0) - sigma value of distribution for Gamma rotation, Pixel barrel layers + + InDetDBTool.Pix_Ba0_X(0.020) - fixed value for X translation, pixel barrel layer 0 + InDetDBTool.Pix_Ba0_Y(0.010) - fixed value for Y translation, pixel barrel layer 0 + InDetDBTool.Pix_Ba0_Z(0) - fixed value for Z translation, pixel barrel layer 0 + InDetDBTool.Pix_Ba0_A(0) - fixed value for Alpha rotation, pixel barrel layer 0 + InDetDBTool.Pix_Ba0_B(0) - fixed value for Beta rotation, pixel barrel layer 0 + InDetDBTool.Pix_Ba0_G(0.006) - fixed value for Gamma rotation, pixel barrel layer 0 + + InDetDBTool.Pix_Ba1_X(-0.030) - fixed value for X translation, pixel barrel layer 1 + InDetDBTool.Pix_Ba1_Y(0.030) - fixed value for Y translation, pixel barrel layer 1 + InDetDBTool.Pix_Ba1_Z(0) - fixed value for Z translation, pixel barrel layer 1 + InDetDBTool.Pix_Ba1_A(0) - fixed value for Alpha rotation, pixel barrel layer 1 + InDetDBTool.Pix_Ba1_B(0) - fixed value for Beta rotation, pixel barrel layer 1 + InDetDBTool.Pix_Ba1_G(0.005) - fixed value for Gamma rotation, pixel barrel layer 1 + + InDetDBTool.Pix_Ba2_X(-0.020) - fixed value for X translation, pixel barrel layer 2 + InDetDBTool.Pix_Ba2_Y(0.030) - fixed value for Y translation, pixel barrel layer 2 + InDetDBTool.Pix_Ba2_Z(0) - fixed value for Z translation, pixel barrel layer 2 + InDetDBTool.Pix_Ba2_A(0) - fixed value for Alpha rotation, pixel barrel layer 2 + InDetDBTool.Pix_Ba2_B(0) - fixed value for Beta rotation, pixel barrel layer 2 + InDetDBTool.Pix_Ba2_G(0.004) - fixed value for Gamma rotation, pixel barrel layer 2 + \endverbatim + + -# Pixel Endcap disks + \verbatim + InDetDBTool.DistTypePixelDisks(1) - distribution type, Pixel endcap disks (in this case: flat) + + InDetDBTool.Pix_Disks_Mean_X(0) - mean value of distribution for X translation, Pixel endcap disks + InDetDBTool.Pix_Disks_Mean_Y(0) - mean value of distribution for Y translation, Pixel endcap disks + InDetDBTool.Pix_Disks_Mean_Z(0) - mean value of distribution for Z translation, Pixel endcap disks + InDetDBTool.Pix_Disks_Mean_A(0) - mean value of distribution for Alpha rotation, Pixel endcap disks + InDetDBTool.Pix_Disks_Mean_B(0) - mean value of distribution for Beta rotation, Pixel endcap disks + InDetDBTool.Pix_Disks_Mean_G(0) - mean value of distribution for Gamma rotation, Pixel endcap disks + + InDetDBTool.Pix_Disks_Sigma_X(0.150) - sigma value of distribution for X translation, Pixel endcap disks + InDetDBTool.Pix_Disks_Sigma_Y(0.150) - sigma value of distribution for Y translation, Pixel endcap disks + InDetDBTool.Pix_Disks_Sigma_Z(0.200) - sigma value of distribution for Z translation, Pixel endcap disks + InDetDBTool.Pix_Disks_Sigma_A(0.001) - sigma value of distribution for Alpha rotation, Pixel endcap disks + InDetDBTool.Pix_Disks_Sigma_B(0.001) - sigma value of distribution for Beta rotation, Pixel endcap disks + InDetDBTool.Pix_Disks_Sigma_G(0.001) - sigma value of distribution for Gamma rotation, Pixel endcap disks + \endverbatim + Pixel Endcap-A disks + \verbatim + InDetDBTool.Pix_EcA_Disk1_X(0) - fixed value for X translation, Pixel-EcA disk 1 + InDetDBTool.Pix_EcA_Disk1_Y(0) - fixed value for X translation, Pixel-EcA disk 1 + InDetDBTool.Pix_EcA_Disk1_Z(0) - fixed value for X translation, Pixel-EcA disk 1 + InDetDBTool.Pix_EcA_Disk1_A(0) - fixed value for Alpha rotation, Pixel-EcA disk 1 + InDetDBTool.Pix_EcA_Disk1_B(0) - fixed value for Beta rotation, Pixel-EcA disk 1 + InDetDBTool.Pix_EcA_Disk1_G(0) - fixed value for Gamma rotation, Pixel-EcA disk 1 + + InDetDBTool.Pix_EcA_Disk2_X(0) - fixed value for X translation, Pixel-EcA disk 2 + InDetDBTool.Pix_EcA_Disk2_Y(0) - fixed value for X translation, Pixel-EcA disk 2 + InDetDBTool.Pix_EcA_Disk2_Z(0) - fixed value for X translation, Pixel-EcA disk 2 + InDetDBTool.Pix_EcA_Disk2_A(0) - fixed value for Alpha rotation, Pixel-EcA disk 2 + InDetDBTool.Pix_EcA_Disk2_B(0) - fixed value for Beta rotation, Pixel-EcA disk 2 + InDetDBTool.Pix_EcA_Disk2_G(0) - fixed value for Gamma rotation, Pixel-EcA disk 2 + + InDetDBTool.Pix_EcA_Disk3_X(0) - fixed value for X translation, Pixel-EcA disk 3 + InDetDBTool.Pix_EcA_Disk3_Y(0) - fixed value for X translation, Pixel-EcA disk 3 + InDetDBTool.Pix_EcA_Disk3_Z(0) - fixed value for X translation, Pixel-EcA disk 3 + InDetDBTool.Pix_EcA_Disk3_A(0) - fixed value for Alpha rotation, Pixel-EcA disk 3 + InDetDBTool.Pix_EcA_Disk3_B(0) - fixed value for Beta rotation, Pixel-EcA disk 3 + InDetDBTool.Pix_EcA_Disk3_G(0) - fixed value for Gamma rotation, Pixel-EcA disk 3 + \endverbatim + Pixel Endcap-C disks + \verbatim + InDetDBTool.Pix_EcC_Disk1_X(0) - fixed value for X translation, Pixel-EcC disk 1 + InDetDBTool.Pix_EcC_Disk1_Y(0) - fixed value for X translation, Pixel-EcC disk 1 + InDetDBTool.Pix_EcC_Disk1_Z(0) - fixed value for X translation, Pixel-EcC disk 1 + InDetDBTool.Pix_EcC_Disk1_A(0) - fixed value for Alpha rotation, Pixel-EcC disk 1 + InDetDBTool.Pix_EcC_Disk1_B(0) - fixed value for Beta rotation, Pixel-EcC disk 1 + InDetDBTool.Pix_EcC_Disk1_G(0) - fixed value for Gamma rotation, Pixel-EcC disk 1 + + InDetDBTool.Pix_EcC_Disk2_X(0) - fixed value for X translation, Pixel-EcC disk 2 + InDetDBTool.Pix_EcC_Disk2_Y(0) - fixed value for X translation, Pixel-EcC disk 2 + InDetDBTool.Pix_EcC_Disk2_Z(0) - fixed value for X translation, Pixel-EcC disk 2 + InDetDBTool.Pix_EcC_Disk2_A(0) - fixed value for Alpha rotation, Pixel-EcC disk 2 + InDetDBTool.Pix_EcC_Disk2_B(0) - fixed value for Beta rotation, Pixel-EcC disk 2 + InDetDBTool.Pix_EcC_Disk2_G(0) - fixed value for Gamma rotation, Pixel-EcC disk 2 + + InDetDBTool.Pix_EcC_Disk3_X(0) - fixed value for X translation, Pixel-EcC disk 3 + InDetDBTool.Pix_EcC_Disk3_Y(0) - fixed value for X translation, Pixel-EcC disk 3 + InDetDBTool.Pix_EcC_Disk3_Z(0) - fixed value for X translation, Pixel-EcC disk 3 + InDetDBTool.Pix_EcC_Disk3_A(0) - fixed value for Alpha rotation, Pixel-EcC disk 3 + InDetDBTool.Pix_EcC_Disk3_B(0) - fixed value for Beta rotation, Pixel-EcC disk 3 + InDetDBTool.Pix_EcC_Disk3_G(0) - fixed value for Gamma rotation, Pixel-EcC disk 3 + \endverbatim + -# SCT Barrel layers + \verbatim + InDetDBTool.DistTypeSCTLayers(0) - distribution type, SCT barrel layers (in this case: fixed values) + + InDetDBTool.SCT_Layers_Mean_X(0) - mean value of distribution for X translation, SCT barrel layers + InDetDBTool.SCT_Layers_Mean_Y(0) - mean value of distribution for Y translation, SCT barrel layers + InDetDBTool.SCT_Layers_Mean_Z(0) - mean value of distribution for Z translation, SCT barrel layers + InDetDBTool.SCT_Layers_Mean_A(0) - mean value of distribution for Alpha rotation, SCT barrel layers + InDetDBTool.SCT_Layers_Mean_B(0) - mean value of distribution for Beta rotation, SCT barrel layers + InDetDBTool.SCT_Layers_Mean_G(0) - mean value of distribution for Gamma rotation, SCT barrel layers + + InDetDBTool.SCT_Layers_Sigma_X(0) - sigma value of distribution for X translation, SCT barrel layers + InDetDBTool.SCT_Layers_Sigma_Y(0) - sigma value of distribution for Y translation, SCT barrel layers + InDetDBTool.SCT_Layers_Sigma_Z(0) - sigma value of distribution for Z translation, SCT barrel layers + InDetDBTool.SCT_Layers_Sigma_A(0) - sigma value of distribution for Alpha rotation, SCT barrel layers + InDetDBTool.SCT_Layers_Sigma_B(0) - sigma value of distribution for Beta rotation, SCT barrel layers + InDetDBTool.SCT_Layers_Sigma_G(0) - sigma value of distribution for Gamma rotation, SCT barrel layers + + InDetDBTool.SCT_Ba0_X(0) - fixed value for X translation, SCT barrel layer 0 + InDetDBTool.SCT_Ba0_Y(0) - fixed value for Y translation, SCT barrel layer 0 + InDetDBTool.SCT_Ba0_Z(0) - fixed value for Z translation, SCT barrel layer 0 + InDetDBTool.SCT_Ba0_A(0) - fixed value for Alpha rotation, SCT barrel layer 0 + InDetDBTool.SCT_Ba0_B(0) - fixed value for Gamma rotation, SCT barrel layer 0 + InDetDBTool.SCT_Ba0_G(-0.001) - fixed value for Gamma rotation, SCT barrel layer 0 + + InDetDBTool.SCT_Ba1_X(0.050) - fixed value for X translation, SCT barrel layer 1 + InDetDBTool.SCT_Ba1_Y(0.040) - fixed value for Y translation, SCT barrel layer 1 + InDetDBTool.SCT_Ba1_Z(0) - fixed value for Z translation, SCT barrel layer 1 + InDetDBTool.SCT_Ba1_A(0) - fixed value for Alpha rotation, SCT barrel layer 1 + InDetDBTool.SCT_Ba1_B(0) - fixed value for Gamma rotation, SCT barrel layer 1 + InDetDBTool.SCT_Ba1_G(0.009) - fixed value for Gamma rotation, SCT barrel layer 1 + + InDetDBTool.SCT_Ba2_X(0.070) - fixed value for X translation, SCT barrel layer 2 + InDetDBTool.SCT_Ba2_Y(0.080) - fixed value for Y translation, SCT barrel layer 2 + InDetDBTool.SCT_Ba2_Z(0) - fixed value for Z translation, SCT barrel layer 2 + InDetDBTool.SCT_Ba2_A(0) - fixed value for Alpha rotation, SCT barrel layer 2 + InDetDBTool.SCT_Ba2_B(0) - fixed value for Gamma rotation, SCT barrel layer 2 + InDetDBTool.SCT_Ba2_G(0.008) - fixed value for Gamma rotation, SCT barrel layer 2 + + InDetDBTool.SCT_Ba3_X(0.100) - fixed value for X translation, SCT barrel layer 3 + InDetDBTool.SCT_Ba3_Y(0.090) - fixed value for Y translation, SCT barrel layer 3 + InDetDBTool.SCT_Ba3_Z(0) - fixed value for Z translation, SCT barrel layer 3 + InDetDBTool.SCT_Ba3_A(0) - fixed value for Alpha rotation, SCT barrel layer 3 + InDetDBTool.SCT_Ba3_B(0) - fixed value for Gamma rotation, SCT barrel layer 3 + InDetDBTool.SCT_Ba3_G(0.007) - fixed value for Gamma rotation, SCT barrel layer 3 + \endverbatim + + -# SCT endcap disks + \verbatim + InDetDBTool.DistTypeSCTDisks(0) - distribution type, SCT endcap disks (in this case: fixed values) + + InDetDBTool.SCT_Disks_Mean_X(0) - mean value of distribution for X translation, SCT disks + InDetDBTool.SCT_Disks_Mean_Y(0) - mean value of distribution for Y translation, SCT disks + InDetDBTool.SCT_Disks_Mean_Z(0) - mean value of distribution for Z translation, SCT disks + InDetDBTool.SCT_Disks_Mean_A(0) - mean value of distribution for Alpha rotation, SCT disks + InDetDBTool.SCT_Disks_Mean_B(0) - mean value of distribution for Beta rotation, SCT disks + InDetDBTool.SCT_Disks_Mean_G(0) - mean value of distribution for Gamma rotation, SCT disks + + InDetDBTool.SCT_Disks_Sigma_X(0) - sigma value of distribution for X translation, SCT disks + InDetDBTool.SCT_Disks_Sigma_Y(0) - sigma value of distribution for Y translation, SCT disks + InDetDBTool.SCT_Disks_Sigma_Z(0) - sigma value of distribution for Z translation, SCT disks + InDetDBTool.SCT_Disks_Sigma_A(0) - sigma value of distribution for Alpha rotation, SCT disks + InDetDBTool.SCT_Disks_Sigma_B(0) - sigma value of distribution for Beta rotation, SCT disks + InDetDBTool.SCT_Disks_Sigma_G(0) - sigma value of distribution for Gamma rotation, SCT disks + \endverbatim + SCT Endcap-C disks + \verbatim + InDetDBTool.SCT_EcA_Disk1_X(0.050) - fixed value for X translation, SCT-EcA disk 1 + InDetDBTool.SCT_EcA_Disk1_Y(0.040) - fixed value for Y translation, SCT-EcA disk 1 + InDetDBTool.SCT_EcA_Disk1_Z(0) - fixed value for Z translation, SCT-EcA disk 1 + InDetDBTool.SCT_EcA_Disk1_A(0) - fixed value for Alpha rotation, SCT-EcA disk 1 + InDetDBTool.SCT_EcA_Disk1_B(0) - fixed value for Beta rotation, SCT-EcA disk 1 + InDetDBTool.SCT_EcA_Disk1_G(-0.001) - fixed value for Gamma rotation, SCT-EcA disk 1 + + InDetDBTool.SCT_EcA_Disk2_X(0.010) - fixed value for X translation, SCT-EcA disk 2 + InDetDBTool.SCT_EcA_Disk2_Y(-0.080) - fixed value for Y translation, SCT-EcA disk 2 + InDetDBTool.SCT_EcA_Disk2_Z(0) - fixed value for Z translation, SCT-EcA disk 2 + InDetDBTool.SCT_EcA_Disk2_A(0) - fixed value for Alpha rotation, SCT-EcA disk 2 + InDetDBTool.SCT_EcA_Disk2_B(0) - fixed value for Beta rotation, SCT-EcA disk 2 + InDetDBTool.SCT_EcA_Disk2_G(0) - fixed value for Gamma rotation, SCT-EcA disk 2 + + InDetDBTool.SCT_EcA_Disk3_X(-0.050) - fixed value for X translation, SCT-EcA disk 3 + InDetDBTool.SCT_EcA_Disk3_Y(0.020) - fixed value for Y translation, SCT-EcA disk 3 + InDetDBTool.SCT_EcA_Disk3_Z(0) - fixed value for Z translation, SCT-EcA disk 3 + InDetDBTool.SCT_EcA_Disk3_A(0) - fixed value for Alpha rotation, SCT-EcA disk 3 + InDetDBTool.SCT_EcA_Disk3_B(0) - fixed value for Beta rotation, SCT-EcA disk 3 + InDetDBTool.SCT_EcA_Disk3_G(0.001) - fixed value for Gamma rotation, SCT-EcA disk 3 + + InDetDBTool.SCT_EcA_Disk4_X(-0.080) - fixed value for X translation, SCT-EcA disk 4 + InDetDBTool.SCT_EcA_Disk4_Y(0.060) - fixed value for Y translation, SCT-EcA disk 4 + InDetDBTool.SCT_EcA_Disk4_Z(0) - fixed value for Z translation, SCT-EcA disk 4 + InDetDBTool.SCT_EcA_Disk4_A(0) - fixed value for Alpha rotation, SCT-EcA disk 4 + InDetDBTool.SCT_EcA_Disk4_B(0) - fixed value for Beta rotation, SCT-EcA disk 4 + InDetDBTool.SCT_EcA_Disk4_G(0.002) - fixed value for Gamma rotation, SCT-EcA disk 4 + + InDetDBTool.SCT_EcA_Disk5_X(0.050) - fixed value for X translation, SCT-EcA disk 5 + InDetDBTool.SCT_EcA_Disk5_Y(0.040) - fixed value for Y translation, SCT-EcA disk 5 + InDetDBTool.SCT_EcA_Disk5_Z(0) - fixed value for Z translation, SCT-EcA disk 5 + InDetDBTool.SCT_EcA_Disk5_A(0) - fixed value for Alpha rotation, SCT-EcA disk 5 + InDetDBTool.SCT_EcA_Disk5_B(0) - fixed value for Beta rotation, SCT-EcA disk 5 + InDetDBTool.SCT_EcA_Disk5_G(0.003) - fixed value for Gamma rotation, SCT-EcA disk 5 + + InDetDBTool.SCT_EcA_Disk6_X(-0.050) - fixed value for X translation, SCT-EcA disk 6 + InDetDBTool.SCT_EcA_Disk6_Y(0.030) - fixed value for Y translation, SCT-EcA disk 6 + InDetDBTool.SCT_EcA_Disk6_Z(0) - fixed value for Z translation, SCT-EcA disk 6 + InDetDBTool.SCT_EcA_Disk6_A(0) - fixed value for Alpha rotation, SCT-EcA disk 6 + InDetDBTool.SCT_EcA_Disk6_B(0) - fixed value for Beta rotation, SCT-EcA disk 6 + InDetDBTool.SCT_EcA_Disk6_G(0.004) - fixed value for Gamma rotation, SCT-EcA disk 6 + + InDetDBTool.SCT_EcA_Disk7_X(-0.030) - fixed value for X translation, SCT-EcA disk 7 + InDetDBTool.SCT_EcA_Disk7_Y(-0.020) - fixed value for Y translation, SCT-EcA disk 7 + InDetDBTool.SCT_EcA_Disk7_Z(0) - fixed value for Z translation, SCT-EcA disk 7 + InDetDBTool.SCT_EcA_Disk7_A(0) - fixed value for Alpha rotation, SCT-EcA disk 7 + InDetDBTool.SCT_EcA_Disk7_B(0) - fixed value for Beta rotation, SCT-EcA disk 7 + InDetDBTool.SCT_EcA_Disk7_G(0.005) - fixed value for Gamma rotation, SCT-EcA disk 7 + + InDetDBTool.SCT_EcA_Disk8_X(0.060) - fixed value for X translation, SCT-EcA disk 8 + InDetDBTool.SCT_EcA_Disk8_Y(0.030) - fixed value for Y translation, SCT-EcA disk 8 + InDetDBTool.SCT_EcA_Disk8_Z(0) - fixed value for Z translation, SCT-EcA disk 8 + InDetDBTool.SCT_EcA_Disk8_A(0) - fixed value for Alpha rotation, SCT-EcA disk 8 + InDetDBTool.SCT_EcA_Disk8_B(0) - fixed value for Beta rotation, SCT-EcA disk 8 + InDetDBTool.SCT_EcA_Disk8_G(0.006) - fixed value for Gamma rotation, SCT-EcA disk 8 + + InDetDBTool.SCT_EcA_Disk9_X(0.080) - fixed value for X translation, SCT-EcA disk 9 + InDetDBTool.SCT_EcA_Disk9_Y(-0.050) - fixed value for Y translation, SCT-EcA disk 9 + InDetDBTool.SCT_EcA_Disk9_Z(0) - fixed value for Z translation, SCT-EcA disk 9 + InDetDBTool.SCT_EcA_Disk9_A(0) - fixed value for Alpha rotation, SCT-EcA disk 9 + InDetDBTool.SCT_EcA_Disk9_B(0) - fixed value for Beta rotation, SCT-EcA disk 9 + InDetDBTool.SCT_EcA_Disk9_G(0.007) - fixed value for Gamma rotation, SCT-EcA disk 9 + \endverbatim + SCT Endcap-C disks + \verbatim + InDetDBTool.SCT_EcC_Disk1_X(-0.050) - fixed value for X translation, SCT-EcC disk 1 + InDetDBTool.SCT_EcC_Disk1_Y(0.050) - fixed value for Y translation, SCT-EcC disk 1 + InDetDBTool.SCT_EcC_Disk1_Z(0) - fixed value for Z translation, SCT-EcC disk 1 + InDetDBTool.SCT_EcC_Disk1_A(0) - fixed value for Alpha rotation, SCT-EcC disk 1 + InDetDBTool.SCT_EcC_Disk1_B(0) - fixed value for Beta rotation, SCT-EcC disk 1 + InDetDBTool.SCT_EcC_Disk1_G(0.008) - fixed value for Gamma rotation, SCT-EcC disk 1 + + InDetDBTool.SCT_EcC_Disk2_X(0) - fixed value for X translation, SCT-EcC disk 2 + InDetDBTool.SCT_EcC_Disk2_Y(0.080) - fixed value for Y translation, SCT-EcC disk 2 + InDetDBTool.SCT_EcC_Disk2_Z(0) - fixed value for Z translation, SCT-EcC disk 2 + InDetDBTool.SCT_EcC_Disk2_A(0) - fixed value for Alpha rotation, SCT-EcC disk 2 + InDetDBTool.SCT_EcC_Disk2_B(0) - fixed value for Beta rotation, SCT-EcC disk 2 + InDetDBTool.SCT_EcC_Disk2_G(0) - fixed value for Gamma rotation, SCT-EcC disk 2 + + InDetDBTool.SCT_EcC_Disk3_X(0.020) - fixed value for X translation, SCT-EcC disk 3 + InDetDBTool.SCT_EcC_Disk3_Y(0.010) - fixed value for Y translation, SCT-EcC disk 3 + InDetDBTool.SCT_EcC_Disk3_Z(0) - fixed value for Z translation, SCT-EcC disk 3 + InDetDBTool.SCT_EcC_Disk3_A(0) - fixed value for Alpha rotation, SCT-EcC disk 3 + InDetDBTool.SCT_EcC_Disk3_B(0) - fixed value for Beta rotation, SCT-EcC disk 3 + InDetDBTool.SCT_EcC_Disk3_G(0.001) - fixed value for Gamma rotation, SCT-EcC disk 3 + + InDetDBTool.SCT_EcC_Disk4_X(0.040) - fixed value for X translation, SCT-EcC disk 4 + InDetDBTool.SCT_EcC_Disk4_Y(-0.080) - fixed value for Y translation, SCT-EcC disk 4 + InDetDBTool.SCT_EcC_Disk4_Z(0) - fixed value for Z translation, SCT-EcC disk 4 + InDetDBTool.SCT_EcC_Disk4_A(0) - fixed value for Alpha rotation, SCT-EcC disk 4 + InDetDBTool.SCT_EcC_Disk4_B(0) - fixed value for Beta rotation, SCT-EcC disk 4 + InDetDBTool.SCT_EcC_Disk4_G(-0.008) - fixed value for Gamma rotation, SCT-EcC disk 4 + + InDetDBTool.SCT_EcC_Disk5_X(0) - fixed value for X translation, SCT-EcC disk 5 + InDetDBTool.SCT_EcC_Disk5_Y(0.030) - fixed value for Y translation, SCT-EcC disk 5 + InDetDBTool.SCT_EcC_Disk5_Z(0) - fixed value for Z translation, SCT-EcC disk 5 + InDetDBTool.SCT_EcC_Disk5_A(0) - fixed value for Alpha rotation, SCT-EcC disk 5 + InDetDBTool.SCT_EcC_Disk5_B(0) - fixed value for Beta rotation, SCT-EcC disk 5 + InDetDBTool.SCT_EcC_Disk5_G(0.003) - fixed value for Gamma rotation, SCT-EcC disk 5 + + InDetDBTool.SCT_EcC_Disk6_X(0.010) - fixed value for X translation, SCT-EcC disk 6 + InDetDBTool.SCT_EcC_Disk6_Y(0.030) - fixed value for Y translation, SCT-EcC disk 6 + InDetDBTool.SCT_EcC_Disk6_Z(0) - fixed value for Z translation, SCT-EcC disk 6 + InDetDBTool.SCT_EcC_Disk6_A(0) - fixed value for Alpha rotation, SCT-EcC disk 6 + InDetDBTool.SCT_EcC_Disk6_B(0) - fixed value for Beta rotation, SCT-EcC disk 6 + InDetDBTool.SCT_EcC_Disk6_G(-0.004) - fixed value for Gamma rotation, SCT-EcC disk 6 + + InDetDBTool.SCT_EcC_Disk7_X(0) - fixed value for X translation, SCT-EcC disk 7 + InDetDBTool.SCT_EcC_Disk7_Y(-0.060) - fixed value for Y translation, SCT-EcC disk 7 + InDetDBTool.SCT_EcC_Disk7_Z(0) - fixed value for Z translation, SCT-EcC disk 7 + InDetDBTool.SCT_EcC_Disk7_A(0) - fixed value for Alpha rotation, SCT-EcC disk 7 + InDetDBTool.SCT_EcC_Disk7_B(0) - fixed value for Beta rotation, SCT-EcC disk 7 + InDetDBTool.SCT_EcC_Disk7_G(0.004) - fixed value for Gamma rotation, SCT-EcC disk 7 + + InDetDBTool.SCT_EcC_Disk8_X(0.030) - fixed value for X translation, SCT-EcC disk 8 + InDetDBTool.SCT_EcC_Disk8_Y(0.030) - fixed value for Y translation, SCT-EcC disk 8 + InDetDBTool.SCT_EcC_Disk8_Z(0) - fixed value for Z translation, SCT-EcC disk 8 + InDetDBTool.SCT_EcC_Disk8_A(0) - fixed value for Alpha rotation, SCT-EcC disk 8 + InDetDBTool.SCT_EcC_Disk8_B(0) - fixed value for Beta rotation, SCT-EcC disk 8 + InDetDBTool.SCT_EcC_Disk8_G(0.006) - fixed value for Gamma rotation, SCT-EcC disk 8 + + InDetDBTool.SCT_EcC_Disk9_X(0.040) - fixed value for X translation, SCT-EcC disk 9 + InDetDBTool.SCT_EcC_Disk9_Y(0.050) - fixed value for Y translation, SCT-EcC disk 9 + InDetDBTool.SCT_EcC_Disk9_Z(0) - fixed value for Z translation, SCT-EcC disk 9 + InDetDBTool.SCT_EcC_Disk9_A(0) - fixed value for Alpha rotation, SCT-EcC disk 9 + InDetDBTool.SCT_EcC_Disk9_B(0) - fixed value for Beta rotation, SCT-EcC disk 9 + InDetDBTool.SCT_EcC_Disk9_G(-0.007) - fixed value for Gamma rotation, SCT-EcC disk 9 + \endverbatim + +\subsection subsection3 LEVEL 1 + +\verbatim + InDetDBTool.DistTypeID(0) - distribution type, ID subdetectors (in this case: fixed values) + + InDetDBTool.ID_Mean_X(0) - mean value of distribution for X translation, ID subdetectors + InDetDBTool.ID_Mean_Y(0) - mean value of distribution for Y translation, ID subdetectors + InDetDBTool.ID_Mean_Z(0) - mean value of distribution for Z translation, ID subdetectors + InDetDBTool.ID_Mean_A(0) - mean value of distribution for Alpha rotation, ID subdetectors + InDetDBTool.ID_Mean_B(0) - mean value of distribution for Beta rotation, ID subdetectors + InDetDBTool.ID_Mean_G(0) - mean value of distribution for Gamma rotation, ID subdetectors + + InDetDBTool.ID_Sigma_X(0) - sigma value of distribution for X translation, ID subdetectors + InDetDBTool.ID_Sigma_Y(0) - sigma value of distribution for Y translation, ID subdetectors + InDetDBTool.ID_Sigma_Z(0) - sigma value of distribution for Z translation, ID subdetectors + InDetDBTool.ID_Sigma_A(0) - sigma value of distribution for Alpha rotation, ID subdetectors + InDetDBTool.ID_Sigma_B(0) - sigma value of distribution for Beta rotation, ID subdetectors + InDetDBTool.ID_Sigma_G(0) - sigma value of distribution for Gamma rotation, ID subdetectors + + InDetDBTool.Pix_X(0.60) - fixed value for X translation, Pixel (whole detector: barrel+endcaps) + InDetDBTool.Pix_Y(1.05) - fixed value for Y translation, Pixel (whole detector: barrel+endcaps) + InDetDBTool.Pix_Z(1.15) - fixed value for Z translation, Pixel (whole detector: barrel+endcaps) + InDetDBTool.Pix_A(-0.10*0.001) - fixed value for Alpha rotation, Pixel (whole detector: barrel+endcaps) + InDetDBTool.Pix_B(0.25*0.001) - fixed value for Beta rotation, Pixel (whole detector: barrel+endcaps) + InDetDBTool.Pix_G(0.65*0.001) - fixed value for Gamma rotation, Pixel (whole detector: barrel+endcaps) + + InDetDBTool.SCT_Ba_X(0.70) - fixed value for X translation, SCT barrel + InDetDBTool.SCT_Ba_Y(1.20) - fixed value for Y translation, SCT barrel + InDetDBTool.SCT_Ba_Z(1.30) - fixed value for Z translation, SCT barrel + InDetDBTool.SCT_Ba_A(0.10*0.001) - fixed value for Alpha rotation, SCT barrel + InDetDBTool.SCT_Ba_B(0.05*0.001) - fixed value for Beta rotation, SCT barrel + InDetDBTool.SCT_Ba_G(0.80*0.001) - fixed value for Gamma rotation, SCT barrel + + InDetDBTool.SCT_EcA_X(2.10) - fixed value for X translation, SCT-EcA + InDetDBTool.SCT_EcA_Y(-0.80) - fixed value for Y translation, SCT-EcA + InDetDBTool.SCT_EcA_Z(1.80) - fixed value for Z translation, SCT-EcA + InDetDBTool.SCT_EcA_A(-0.25*0.001) - fixed value for Alpha rotation, SCT-EcA + InDetDBTool.SCT_EcA_B(0) - fixed value for Beta rotation, SCT-EcA + InDetDBTool.SCT_EcA_G(-0.50*0.001) - fixed value for Gamma rotation, SCT-EcA + + InDetDBTool.SCT_EcC_X(-1.90) - fixed value for X translation, SCT-EcA + InDetDBTool.SCT_EcC_Y(2.00) - fixed value for Y translation, SCT-EcA + InDetDBTool.SCT_EcC_Z(-3.10) - fixed value for Z translation, SCT-EcA + InDetDBTool.SCT_EcC_A(-0.10*0.001) - fixed value for Alpha rotation, SCT-EcA + InDetDBTool.SCT_EcC_B(0.05*0.001) - fixed value for Beta rotation, SCT-EcA + InDetDBTool.SCT_EcC_G(0.40*0.001) - fixed value for Gamma rotation, SCT-EcA +\endverbatim + +\section InDetAlignTrackSelTool + +This tool provides a standard track selection for 'alignment quality' tracks +(for the SCT and pixels), +allowing this to be standardised among the various alignment algorithms +being developed. It currently makes a selection based on the track +transverse momentum (p_T), number of 'holes' on the track in the SCT and +pixel detectors (i.e. layers crossed by the track but having no associated +hit), and number of 'shared' hits assigned to more than one track. The latter +quantity is calcuated internally by the tool, by analysing the set of tracks +it is passed in a preliminary fillMap call (since the information is not +calculated in the TrkTrack framework). + +The tool is used as follows: + + - create map by passing as argument to the tool a collection of tracks: +<pre> +// track collection +const DataVector<Trk::Track>* trks; + +//retrieve pre-existing collection from StoreGate +sc = m_StoreGate->retrieve(trks, m_trackCollName); + +// fill the map +m_alignTrackSelTool->fillMap(trks); +</pre> + + - loop over tracks and check if each track passes cuts or not: +<pre> +DataVector<Trk::Track>::const_iterator trkit; +for (trkit=trks->begin(); trkit!=trks->end();++trkit) { + int stat = m_alignTrackSelTool->getStatus(**trkit); + ... + } +</pre> +The method getStatus returns 1 if the track passes cuts, or 0 otherwise. +Note that the fillMap step is essential to correctly fill the information +on shared hits. + +An example of the use of this tool can be found in the InDetAlignNt +algorithm in InnerDetector/InDetAlignment/InDetAlignAlgs. + +@section InDetAlignTrackSelToolJobOptions InDetAlignTrackSelTool JobOptions + +@include InDetAlignTrackSelTool_jobOptions.py + +@section requirements requirements + +@include requirements + +@section uses Packages used + +@htmlinclude used_packages.html + +*/ + diff --git a/InnerDetector/InDetAlignTools/InDetAlignGenTools/share/ConfiguredInDetAlignTrackSelTool.py b/InnerDetector/InDetAlignTools/InDetAlignGenTools/share/ConfiguredInDetAlignTrackSelTool.py new file mode 100755 index 0000000000000000000000000000000000000000..4969d1530333cbda03b444ab2745132cfcefaafb --- /dev/null +++ b/InnerDetector/InDetAlignTools/InDetAlignGenTools/share/ConfiguredInDetAlignTrackSelTool.py @@ -0,0 +1,133 @@ +##################################################################### +# ConfiguredInDetAlignTrackSelTool +##################################################################### +# Python Setup Class InDetAlignGenTools/InDetAlignTrackSelTool +# +# Sergio Gonzalez-Sevilla +# Date: 21/12/2005 +# +##################################################################### + +class ConfiguredInDetAlignTrackSelTool : + def __init__ (self, + instname='InDetAlignTrackSelTool', + propagatorType='RungeKutta') : + + # Load the Dlls + if 'InDetAlignGenTools' not in theApp.Dlls: + theApp.Dlls += [ 'InDetAlignGenTools' ] + + # Set names + self.__instname__ = instname + self.__svcname__ = 'InDetAlignTrackSelTool' + + self.__toolname__ = 'ToolSvc.'+self.__instname__ + self.__thistool__ = Service( self.__toolname__ ) + + # Set options for tool + self.__thistool__.MinMomentum = 0 + self.__thistool__.MinPt = 2 + self.__thistool__.MaxShared = 0 + self.__thistool__.MaxHoles = 1 + self.__thistool__.MinChi2Prob = 0.2 + + ToolSvc = Service("ToolSvc") + + # Propagator + if propagatorType is "STEP": + from TrkExSTEP_Propagator.TrkExSTEP_PropagatorConf import Trk__STEP_Propagator as Propagator + else: + from TrkExRungeKuttaPropagator.TrkExRungeKuttaPropagatorConf import Trk__RungeKuttaPropagator as Propagator + InDetAlignTrackSelPropagator = Propagator(name = 'InDetAlignTrackSelPropagator') + ToolSvc += InDetAlignTrackSelPropagator + print InDetAlignTrackSelPropagator + + # Navigator + from TrkExTools.TrkExToolsConf import Trk__Navigator + InDetAlignTrackSelNavigator = Trk__Navigator(name = 'InDetAlignTrackSelNavigator') + ToolSvc += InDetAlignTrackSelNavigator + print InDetAlignTrackSelNavigator + + # Material updator + from TrkExTools.TrkExToolsConf import Trk__MaterialEffectsUpdator + InDetAlignTrackSelMaterialUpdator = Trk__MaterialEffectsUpdator(name = "InDetAlignTrackSelMaterialEffectsUpdator") + ToolSvc += InDetAlignTrackSelMaterialUpdator + print InDetAlignTrackSelMaterialUpdator + + # Extrapolator + from TrkExTools.TrkExToolsConf import Trk__Extrapolator + InDetAlignTrackSelExtrapolator = Trk__Extrapolator(name = 'InDetAlignTrackSelExtrapolator', + Propagators = [ InDetAlignTrackSelPropagator.getType() ], + PropagatorInstances = [ InDetAlignTrackSelPropagator.getName() ], + Navigator = InDetAlignTrackSelNavigator, + MaterialEffectsUpdator = InDetAlignTrackSelMaterialUpdator) + ToolSvc += InDetAlignTrackSelExtrapolator + print InDetAlignTrackSelExtrapolator + + # associator + from InDetAssociationTools.InDetAssociationToolsConf import InDet__InDetPRD_AssociationToolGangedPixels + InDetAlignTrackSelAssociationTool = InDet__InDetPRD_AssociationToolGangedPixels() + ToolSvc += InDetAlignTrackSelAssociationTool + + # summary + from InDetTrackSummaryHelperTool.InDetTrackSummaryHelperToolConf import InDet__InDetTrackSummaryHelperTool + InDetAlignTrackSelSummaryTool = InDet__InDetTrackSummaryHelperTool(name = "InDetAlignTrackSelSummaryTool", + Extrapolator = InDetAlignTrackSelExtrapolator, + AssoTool = InDetAlignTrackSelAssociationTool, + DoSharedHits = True) + ToolSvc += InDetAlignTrackSelSummaryTool + print InDetAlignTrackSelSummaryTool + + # Set message level + def msgStreamLevel(self, level): + self.__thistool__.OutputLevel = level + + # Set track cuts: min Momentum(GeV) + def MinMomentum(self,minMom=0): + self.__thistool__.MinMomentum=minMom + + # Set track cuts: min pt(GeV) + def MinPt(self,minpt=2): + self.__thistool__.MinPt=minpt + + # Set track cuts: max shared hits + def MaxShared(self,maxshared=0): + self.__thistool__.MaxShared=maxshared + + # Set track cuts: max holes + def MaxHoles(self,maxholes=1): + self.__thistool__.MaxHoles=maxholes + + # Set track cuts: min chi2Prob + def MinChi2Prob(self,minchi2prob=0.2): + self.__thistool__.MinChi2Prob=minchi2prob + + # Set track cuts: all together + def SetCuts(self,minMom=0,minpt=2,maxshared=0,maxholes=1,minchi2prob=0.2): + self.__thistool__.MinMomentum=minMom + self.__thistool__.MinPt=minpt + self.__thistool__.MaxShared=maxshared + self.__thistool__.MaxHoles=maxholes + self.__thistool__.MinChi2Prob=minchi2prob + + # Return method for the service name + def name(self): + return self.__svcname__ + + # Return method for the instance name + def instance(self): + return self.__instname__ + + # Output info + def printInfo(self) : + print '***** ConfiguredInDetAlignTrackSelTool ***************************************' + print '* - ToolName: ', self.__svcname__,' (Instance: ', self.__instname__, ')' + print '* --------------------------------------------------------------------------------' + print '* - Selected track cuts:' + print '* o Min Momentum(GeV) = ', self.__thistool__.MinMomentum + print '* o Min Pt(GeV) = ', self.__thistool__.MinPt + print '* o Max Shared = ', self.__thistool__.MaxShared + print '* o Max Holes = ', self.__thistool__.MaxHoles + print '* o Min Chi2Prob = ', self.__thistool__.MinChi2Prob + print '*******************************************************************************' + diff --git a/InnerDetector/InDetAlignTools/InDetAlignGenTools/src/InDetAlignDBTool.cxx b/InnerDetector/InDetAlignTools/InDetAlignGenTools/src/InDetAlignDBTool.cxx new file mode 100755 index 0000000000000000000000000000000000000000..c46f904ca0907282ffd6ac876be7bc2bc9072312 --- /dev/null +++ b/InnerDetector/InDetAlignTools/InDetAlignGenTools/src/InDetAlignDBTool.cxx @@ -0,0 +1,1356 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +// InDetAlignDBTool.cxx +// AlgTool to manage the SCT/pixel AlignableTransforms in the conditions store +// Richard Hawkings, started 8/4/05 + +#include <fstream> +#include <iostream> + +#include "GaudiKernel/NTuple.h" +#include "GaudiKernel/INTupleSvc.h" +#include "GaudiKernel/SmartDataPtr.h" + +#include "GaudiKernel/RndmGenerators.h" +#include "GaudiKernel/IRndmGenSvc.h" + +#include "GeoPrimitives/CLHEPtoEigenConverter.h" + +#include "AthenaKernel/IAthenaOutputStreamTool.h" +#include "InDetReadoutGeometry/SiDetectorManager.h" +#include "InDetReadoutGeometry/PixelDetectorManager.h" +#include "InDetReadoutGeometry/SCT_DetectorManager.h" +#include "InDetReadoutGeometry/SiDetectorElementCollection.h" +#include "InDetReadoutGeometry/SiDetectorElement.h" + +#include "InDetIdentifier/PixelID.h" +#include "InDetIdentifier/SCT_ID.h" + +#include "DetDescrConditions/AlignableTransformContainer.h" +#include "RegistrationServices/IIOVRegistrationSvc.h" + +#include "InDetAlignGenTools/InDetAlignDBTool.h" + +// alignment DBS ntuple 9002 definition +NTuple::Item<long> nt_dettype; +NTuple::Item<long> nt_bec; +NTuple::Item<long> nt_layer; +NTuple::Item<long> nt_ring; +NTuple::Item<long> nt_sector; +NTuple::Item<long> nt_side; +NTuple::Item<long> nt_level; +NTuple::Item<float> nt_xofs; +NTuple::Item<float> nt_yofs; +NTuple::Item<float> nt_zofs; +NTuple::Item<float> nt_phi; +NTuple::Item<float> nt_theta; +NTuple::Item<float> nt_psi; +NTuple::Item<float> nt_alpha; +NTuple::Item<float> nt_beta; +NTuple::Item<float> nt_gamma; + +// name of the folder for ID alignment information +#define IND_ALIGN "/Indet/Align" + +InDetAlignDBTool::InDetAlignDBTool(const std::string& type, + const std::string& name, const IInterface* parent) + : AthAlgTool(type,name,parent), + p_toolsvc("ToolSvc",name), + m_pixid(0), + m_sctid(0), + m_pixman(0), + m_sctman(0), + par_newdb(true), + par_scttwoside(false), + par_fake(0), + par_condstream("CondStream1"), + par_dbroot( IND_ALIGN ), + par_dbkey( IND_ALIGN ), + par_oldTextFile(false) +{ + declareInterface<IInDetAlignDBTool>(this); + declareProperty("IToolSvc", p_toolsvc); + declareProperty("NewDB", par_newdb); + declareProperty("SCTTwoSide", par_scttwoside); + declareProperty("FakeDB", par_fake); + declareProperty("CondStream", par_condstream); + declareProperty("DBRoot", par_dbroot,"name of the root folder for constants"); + declareProperty("DBKey", par_dbkey,"base part of the key for loading AlignableTransforms"); + declareProperty("OldTextFile", par_oldTextFile); + +} + +InDetAlignDBTool::~InDetAlignDBTool() +{} + +StatusCode InDetAlignDBTool::initialize() +{ + + ATH_MSG_DEBUG("InDetAlignDBTool initialize instance: " << name() ); + + // get storegate access to conditions store + if (detStore().retrieve().isFailure()){ + ATH_MSG_FATAL("Detector store not found"); + } + + if ( p_toolsvc.retrieve().isFailure() ) { + ATH_MSG_FATAL( "Failed to retrieve service " << p_toolsvc ); + return StatusCode::FAILURE; + } else + ATH_MSG_DEBUG( "Retrieved service " << p_toolsvc ); + + // attempt to get ID helpers from detector store + // (relying on GeoModel to put them) + m_alignobjs.clear(); + m_alignchans.clear(); + int ndet[2]; + ndet[0]=0; + ndet[1]=0; + if ((StatusCode::SUCCESS!=detStore()->retrieve(m_pixid)) || + (StatusCode::SUCCESS!=detStore()->retrieve(m_sctid))) { + // attempt failed - optionally fake up the geometry + if (par_fake==1) { + ATH_MSG_DEBUG("Initialising fake full ATLAS geometry"); + fakeGeom(3,3,4,9); + } else if (par_fake==2) { + ATH_MSG_DEBUG("Initialising fake CTB geometry"); + fakeGeom(3,0,4,0); + } else { + ATH_MSG_FATAL( "Problem retrieving ID helpers"); + return StatusCode::FAILURE; + } + } else { + // setup list of alignable transforms from geometry + int chan[3]; + for (int i=0;i<3;++i) chan[i]=100*i; + std::string man_name; + InDetDD::SiDetectorElementCollection::const_iterator iter,itermin,itermax; + for (int idet=1;idet<3;++idet) { + if (idet==1) { + if (StatusCode::SUCCESS!=detStore()->retrieve(m_pixman,"Pixel") || m_pixman==0) { + ATH_MSG_WARNING( "Could not find pixel manager "); + return StatusCode::FAILURE; + } else { + itermin=m_pixman->getDetectorElementBegin(); + itermax=m_pixman->getDetectorElementEnd(); + } + } else { + if (StatusCode::SUCCESS!=detStore()->retrieve(m_sctman,"SCT") || m_sctman==0) { + ATH_MSG_FATAL("Could not find SCT manager "); + return StatusCode::FAILURE; + } else { + itermin=m_sctman->getDetectorElementBegin(); + itermax=m_sctman->getDetectorElementEnd(); + } + } + for (iter=itermin;iter!=itermax;++iter) { + const InDetDD::SiDetectorElement* element=*iter; + if (element!=0) { + const Identifier ident=element->identify(); + int det,bec,layer,ring,sector,side; + if (idToDetSet(ident,det,bec,layer,ring,sector,side)) { + if (det==idet) { + std::string level[3]; + for (int i=0;i<3;++i) { + level[i]=dirkey(det,bec,layer,1+i,sector); + // add this to list if not seen already + std::vector<std::string>::const_iterator ix= + find(m_alignobjs.begin(),m_alignobjs.end(),level[i]); + if (ix==m_alignobjs.end()) { + m_alignobjs.push_back(level[i]); + m_alignchans.push_back(chan[i]++); + } + } + ++ndet[idet-1]; + } else { + ATH_MSG_ERROR("Detector of wrong type in list " << idet ); + } + } else { + ATH_MSG_ERROR("Si detector element type " << idet << " has no detset conversion" ); + } + } + } + } + } + if (msgLvl(MSG::DEBUG)) { + ATH_MSG_DEBUG( "Database root folder " << par_dbroot ); + ATH_MSG_DEBUG( "Geometry initialisation sees " << ndet[0] << + " pixel and " << ndet[1] << " SCT modules giving " << m_alignobjs.size() + << " alignment keys" ); + ATH_MSG_DEBUG("Keys/channels are:"); + + for (unsigned int i=0;i<m_alignobjs.size();++i) + ATH_MSG_DEBUG( " " << m_alignobjs[i] << " [" << m_alignchans[i] << "]" ); + + if (par_newdb) + ATH_MSG_DEBUG("Assuming new COOL alignment DB model based on AlignableTransformContainer"); + else + ATH_MSG_DEBUG("Assuming old (Lisbon) alignment DB model based on separate AlignableTransforms"); + } + + return StatusCode::SUCCESS; +} + +StatusCode InDetAlignDBTool::finalize() +{ + ATH_MSG_DEBUG( "InDetAlignDBTool finalize method called" ); + return StatusCode::SUCCESS; +} + +void InDetAlignDBTool::createDB() const +{ + + ATH_MSG_DEBUG("createDB method called"); + // check not running in fake mode (need real geometry here) + if (par_fake) { + ATH_MSG_FATAL("Cannot create new database when geometry is faked"); + } + AlignableTransform* pat; + AlignableTransformContainer* patc=0; + // loop over all SiDetectorElements (pixel and SCT) and fill corresponding + // AlignableTransform objects with default values + + // first create the empty AlignableTransform objects in TDS + if (par_newdb) { + // check object does not already exist + if (detStore()->contains<AlignableTransformContainer>(par_dbroot)) { + ATH_MSG_WARNING("createDB: AlignableTransformContainer already exists"); + return; + } + // put them in a collection /Indet/Align + ATH_MSG_DEBUG( "Setup database structures in AlignableTransformContainer"); + patc=new AlignableTransformContainer; + } + else { + ATH_MSG_DEBUG( "Setup separate AlignableTransform for each layer"); + } + + if (msgLvl(MSG::DEBUG)) { + if (par_scttwoside) ATH_MSG_DEBUG( "Produce separate transforms for each side of SCT modules" ); + else ATH_MSG_DEBUG( "Treat both sides of SCT module as single entity" ); + } + + for (unsigned int i=0;i<m_alignobjs.size();++i) { + pat=new AlignableTransform(m_alignobjs[i]); + if (par_newdb) { + // add to collection and set corresponding channel number + patc->push_back(pat); + patc->add(m_alignchans[i]); + } else { + // store directly in SG + // first check object not already there + if (detStore()->contains<AlignableTransform>(m_alignobjs[i])) { + ATH_MSG_FATAL( "create DB: AlignableTransform " << m_alignobjs[i] << " already exists" ); + return; + } + if (StatusCode::SUCCESS!=detStore()->record(pat,m_alignobjs[i])) + ATH_MSG_ERROR( "Could not record AlignableTransform "<< m_alignobjs[i] ); + } + } + if (par_newdb) { + // record collection in SG + if (StatusCode::SUCCESS!=detStore()->record(patc,par_dbroot)) + ATH_MSG_ERROR("Could not record AlignableTransformContainer"); + ATH_MSG_DEBUG( "Collection has size " << patc->size() ); + } + + // now loop over all detector modules and add null level 3 transforms + std::vector<std::string> level2; + InDetDD::SiDetectorElementCollection::const_iterator iter,itermin,itermax; + for (int idet=1;idet<3;++idet) { + if (idet==1) { + itermin=m_pixman->getDetectorElementBegin(); + itermax=m_pixman->getDetectorElementEnd(); + } else { + itermin=m_sctman->getDetectorElementBegin(); + itermax=m_sctman->getDetectorElementEnd(); + } + for (iter=itermin;iter!=itermax;++iter) { + const InDetDD::SiDetectorElement* element=*iter; + if (element!=0) { + const Identifier ident=element->identify(); + std::string key=dirkey(ident,3); + // do not produce AlignableTrasnforms for SCT side 1 if option set + if (!(m_sctid->is_sct(ident) && m_sctid->side(ident)==1) || + par_scttwoside) { + if ((pat=getTransPtr(key))) { + pat->add(ident,Amg::EigenTransformToCLHEP( Amg::Transform3D::Identity() ) ); + } else { + ATH_MSG_ERROR( "Cannot retrieve AlignableTransform for key " << key ); + } + // add level 2 transform if needed - do this the first time a module + // for this level 3 key is seen + std::vector<std::string>::const_iterator ix= + find(level2.begin(),level2.end(),key); + if (ix==level2.end()) { + level2.push_back(key); + // construct identifier of level 2 transform + Identifier ident2; + if (m_pixid->is_pixel(ident)) { + ident2=m_pixid->wafer_id(m_pixid->barrel_ec(ident), + m_pixid->layer_disk(ident),m_pixid->phi_module(ident),0); // needed to be extended to phi-module due to DBM + } else if (m_sctid->is_sct(ident)) { + ident2=m_sctid->wafer_id(m_sctid->barrel_ec(ident), + m_sctid->layer_disk(ident),0,0,0); + } + std::string key2=dirkey(ident,2); + if ((pat=getTransPtr(key2))) { + pat->add(ident2,Amg::EigenTransformToCLHEP( Amg::Transform3D::Identity() ) ); + } else { + ATH_MSG_ERROR( "Cannot retrieve AlignableTransform for key " << key2 ); + } + } + } + } + } + } + // create the global ID object with positions for the pixels (one unit) + // and SCT barrel/endcap + Identifier ident1; + std::string key1=dirkey(ident1,1); + if ((pat=getTransPtr(key1))) { + Amg::Transform3D globshift; + globshift.setIdentity(); + // pixel - barrel and endcap as one, treated as barrel layer 0 module 0 + ident1=m_pixid->wafer_id(0,0,0,0); + pat->add(ident1,Amg::EigenTransformToCLHEP(globshift)); + // SCT barrel - barrel 0 module 0 + ident1=m_sctid->wafer_id(0,0,0,0,0); + pat->add(ident1,Amg::EigenTransformToCLHEP(globshift)); + // SCT endcaps A and C + ident1=m_sctid->wafer_id(-2,0,0,0,0); + pat->add(ident1,Amg::EigenTransformToCLHEP(globshift)); + ident1=m_sctid->wafer_id(2,0,0,0,0); + pat->add(ident1,Amg::EigenTransformToCLHEP(globshift)); + } else { + ATH_MSG_ERROR( "Cannot retrieve AlignableTransform for key " << key1 ); + } + // sort the created objects (in case, usually come out sorted from GeoModel) + sortTrans(); + // list out size of all created objects + ATH_MSG_DEBUG( "Dumping size of created AlignableTransform objects"); + for (unsigned int i=0;i<m_alignobjs.size();++i) + if ((pat=getTransPtr(m_alignobjs[i]))) pat->print(); +} + +bool InDetAlignDBTool::idToDetSet(const Identifier ident, int& det, int& bec, + int& layer, int& ring, int& sector, int& side) const { + // transform Identifier to list of integers specifiying dettype,bec,layer + // ring, sector, side + // note bec is +-1 or 0, not +-2 as returned by idenfitiers + + bool resok=false; + if (m_pixid->is_pixel(ident)) { + det=1; + bec=m_pixid->barrel_ec(ident)/2; + layer=m_pixid->layer_disk(ident); + ring=m_pixid->eta_module(ident); + sector=m_pixid->phi_module(ident); + side=0; + resok=true; + } else if (m_sctid->is_sct(ident)) { + det=2; + bec=m_sctid->barrel_ec(ident)/2; + layer=m_sctid->layer_disk(ident); + ring=m_sctid->eta_module(ident); + sector=m_sctid->phi_module(ident); + side=m_sctid->side(ident); + resok=true; + } + return resok; +} + +std::string InDetAlignDBTool::dirkey(const Identifier& ident, + const int level) const { + // given SCT or pixel identifier, and level (1,2 or 3) return + // directory key name for associated alignment data + int det,bec,layer,ring,sector,side; + idToDetSet(ident,det,bec,layer,ring,sector,side); + return dirkey(det,bec,layer,level,sector); +} + +// This function is redundant for the main code. +// Kept for now as I did not want to touch functions like dispCSC() etc. +std::string InDetAlignDBTool::dirkey(const int det,const int bec,const + int layer, const int level) const { + // given SCT/pixel det/bec/layer, and level (1,2 or 3) return + // directory key name for associated alignment data + std::ostringstream result; + result << par_dbkey << "/" ; + if (level==1) { + result << "ID"; + } else { + if (det==1) result << "PIX"; + if (det==2) result << "SCT"; + if (level==3) { + if (bec==1) result << "EA"; + if (bec==0) result << "B"; + if (bec==-1) result << "EC"; + result << 1+layer; + } + } + return result.str(); +} + +std::string InDetAlignDBTool::dirkey(const int det,const int bec,const + int layer, const int level, const int sector) const { + // given SCT/pixel det/bec/layer/sector, and level (1,2 or 3) return + // directory key name for associated alignment data + std::ostringstream result; + result << par_dbkey << "/" ; + if (level==1) { + result << "ID"; + } else { + if (det==1) result << "PIX"; + if (det==2) result << "SCT"; + if (level==3) { + if (det==1 && abs(bec)==2) result << DBMkey(det,bec,level,sector); + else { + if (bec==1) result << "EA"; + if (bec==0) result << "B"; + if (bec==-1) result << "EC"; + result << 1+layer; + } + } + } + return result.str(); +} + +std::string InDetAlignDBTool::DBMkey(const int det,const int bec, + const int level, const int sector) const { + // given SCT/pixel det/bec/layer/sector, and level (1,2 or 3) return + // additional directory key name for associated DBM alignment + std::ostringstream result; + if (det==1 && level==3 && abs(bec)==2) { // slightly unnecessary check + if (bec==2) result << "EADBM"; + if (bec==-2) result << "ECDBM"; + result << 1+sector; + } + return result.str(); +} + + +void InDetAlignDBTool::dispGroup(const int dettype, const int bec, + const int layer,const int ring, const int sector, + const float rphidisp, const float rdisp, const float zdisp, + const int syst, const int level, const int skip) const { + + ATH_MSG_DEBUG( "dispGroup called: level " << level << " syst " << syst); + int nmod=0; + // random number service + IRndmGenSvc* randsvc; + if (StatusCode::SUCCESS!=service("RndmGenSvc",randsvc,true)) + ATH_MSG_ERROR("Cannot find RndmGenSvc" ); + Rndm::Numbers gauss(randsvc,Rndm::Gauss(0.,1.)); + if (skip>0) { + ATH_MSG_DEBUG("Skip random numbers " << skip ); + for (int i=0;i<skip;++i) gauss(); + } + // for syst 5, choose random shifts based on the input numbers + float rpd=0,rd=0,zd=0; + if (syst==5) { + rpd=rphidisp*gauss(); + rd=rdisp*gauss(); + zd=zdisp*gauss(); + } + // keep a list of level1/2 transform IDs to make sure they are only set once + std::vector<Identifier> lvl12id; + // loop over all pixel and SCT modules + AlignableTransform* pat; + InDetDD::SiDetectorElementCollection::const_iterator iter,itermin,itermax; + for (int idet=1;idet<3;++idet) { + if (idet==1) { + itermin=m_pixman->getDetectorElementBegin(); + itermax=m_pixman->getDetectorElementEnd(); + } else { + itermin=m_sctman->getDetectorElementBegin(); + itermax=m_sctman->getDetectorElementEnd(); + } + for (iter=itermin;iter!=itermax;++iter) { + const InDetDD::SiDetectorElement* element=*iter; + if (element!=0) { + const Identifier ident=element->identify(); + int mdet,mbec,mlayer,mring,msector,mside; + idToDetSet(ident,mdet,mbec,mlayer,mring,msector,mside); + // find matching modules - note side=1 modules never touched + if ((dettype==-1 || mdet==dettype) && (bec==-1 || fabs(2*mbec)==bec) && + (layer==-1 || mlayer==layer) && (ring==-1 || mring==ring) && + (sector== -1 || msector==sector) && mside==0) { + // displace this module - first choose displacement type + // dont choose new displacements if seeing second side of SCT module + // ensures they both move together + // depends on the side1 module immediatly following side 0 in list + // which is currently the case - fragile + // also for syst 6 choose number only for new ring (eta) slice + if (dettype!=2 || mside!=1) { + if (syst==2 || syst==4 || (syst==6 && mring==-6)) { + rpd=rphidisp*gauss(); + rd=rdisp*gauss(); + zd=zdisp*gauss(); + if (syst==6) ATH_MSG_DEBUG("New rndm at layer/ring " << + mlayer << " " << mring << " z " << zd ); + } else if (syst<5) { + rpd=rphidisp; + rd=rdisp; + zd=zdisp; + } + } + // interpretation as rphi/r or x/y + float xd,yd; + if (syst<=2 || syst==6) { + // rphi displacement - calculate from module position in x/y + const Amg::Vector3D modcent=element->center(); + float dx=modcent.x(); + float dy=modcent.y(); + float dr=sqrt(dx*dx+dy*dy); + xd=(rd*dx-rpd*dy)/dr; + yd=(rd*dy+rpd*dx)/dr; + } else { + xd=rpd; + yd=rd; + } + // find the corresponding AlignableTransform object + std::string key=dirkey(mdet,mbec,mlayer,level); + // first get as const as transforms might have been read in + const AlignableTransform* cpat=cgetTransPtr(key); + pat=const_cast<AlignableTransform*>(cpat); + if (pat) { + Identifier ident2; + bool update=true; + if (level==3) { + ident2=ident; + } else if (level==2) { + // identifier for layer in level 2 transform + if (mdet==1) { + ident2=m_pixid->wafer_id(m_pixid->barrel_ec(ident), + m_pixid->layer_disk(ident),0,0); + } else { + ident2=m_sctid->wafer_id(m_sctid->barrel_ec(ident), + m_sctid->layer_disk(ident),0,0,0); + } + // check this identifier has not been updated before + std::vector<Identifier>::const_iterator ix= + find(lvl12id.begin(),lvl12id.end(),ident2); + if (ix==lvl12id.end()) { + lvl12id.push_back(ident2); + } else { + update=false; + } + } else { + // identifier for ID + if (mdet==1) { + ident2=m_pixid->wafer_id(0,0,0,0); + } else { + ident2=m_sctid->wafer_id(0,0,0,0,0); + } + // check this identifier has not been updated before + std::vector<Identifier>::const_iterator ix= + find(lvl12id.begin(),lvl12id.end(),ident2); + if (ix==lvl12id.end()) { + lvl12id.push_back(ident2); + } else { + update=false; + } + } + // update, adding to any existing shift + if (update) { + + Amg::Transform3D shift = Amg::Translation3D(xd,yd,zd) * Amg::RotationMatrix3D::Identity(); + pat->tweak(ident2,Amg::EigenTransformToCLHEP(shift)); + ATH_MSG_VERBOSE( "Updated module " << mdet << "," << mbec + << "," << mlayer << "," << mring << "," << msector << " to xyz" << + xd << "," << yd << "," << zd ); + ++nmod; + } + } else { + ATH_MSG_ERROR("Cannot find AlignableTransform for key" << key ); + } + } + } + } + } + ATH_MSG_DEBUG( "Added displacement to " << nmod << " modules " << dettype << "," + << bec << "," << layer << " [" << rphidisp << "," << rdisp + << "," << zdisp << "]" + << " type " << syst ); +} + +void InDetAlignDBTool::writeFile(const bool ntuple, const std::string file) + const { + std::ofstream* outfile=0; + INTupleSvc* ntsvc; + if (StatusCode::SUCCESS!=service("NTupleSvc",ntsvc,true)) + ATH_MSG_ERROR("Cannot find NTupleSvc" ); + const std::string path=file+"/9002"; + NTuplePtr nt(ntsvc,path); + + if (ntuple) { + ATH_MSG_DEBUG( "writeFile: Write AlignableTransforms on ntuple 9002, path: " << file ); + const int ntid=9002; + if (nt) { + ATH_MSG_DEBUG( "Ntuple " << path << " is already booked" ); + } else { + ATH_MSG_DEBUG("Attempting to book ntuple " << path ); + nt=ntsvc->book(file,ntid,CLID_ColumnWiseTuple,"AlignDB"); + if (!nt) ATH_MSG_ERROR("Ntuple booking fails" ); + } + StatusCode sc; + sc=nt->addItem("MODPROP/DetType",nt_dettype); + sc=nt->addItem("MODPROP/Bec",nt_bec); + sc=nt->addItem("MODPROP/Layer",nt_layer); + sc=nt->addItem("MODPROP/Ring",nt_ring); + sc=nt->addItem("MODPROP/Sector",nt_sector); + sc=nt->addItem("MODPROP/Side",nt_side); + sc=nt->addItem("MODPROP/Level",nt_level); + sc=nt->addItem("MODPROP/Xofs",nt_xofs); + sc=nt->addItem("MODPROP/Yofs",nt_yofs); + sc=nt->addItem("MODPROP/Zofs",nt_zofs); + sc=nt->addItem("MODPROP/Phi",nt_phi); + sc=nt->addItem("MODPROP/Theta",nt_theta); + sc=nt->addItem("MODPROP/Psi",nt_psi); + sc=nt->addItem("MODPROP/Phi",nt_phi); + sc=nt->addItem("MODPROP/Theta",nt_theta); + sc=nt->addItem("MODPROP/Psi",nt_psi); + if (sc!=StatusCode::SUCCESS) ATH_MSG_ERROR( + "Error booking ntuple 9002 contents" ); + } else { + ATH_MSG_DEBUG( "writeFile: Write AlignableTransforms on text file: " << file ); + outfile=new std::ofstream(file.c_str()); + } + int nobj=0; + int ntrans=0; + for (std::vector<std::string>::const_iterator iobj=m_alignobjs.begin(); + iobj!=m_alignobjs.end();++iobj) { + const AlignableTransform* pat; + if ((pat=cgetTransPtr(*iobj))) { + ++nobj; + if (!ntuple) *outfile << *iobj << std::endl; + for (AlignableTransform::AlignTransMem_citr cit=pat->begin(); + cit!=pat->end();++cit) { + const Identifier& ident=cit->identify(); + const Amg::Transform3D& trans=Amg::CLHEPTransformToEigen( cit->transform() ); + int det,bec,layer,ring,sector,side; + float dx,dy,dz,phi,theta,psi; + if (!idToDetSet(ident,det,bec,layer,ring,sector,side)) { + // can fail for testbeam whe identifier with all layer + // and wafer indices set to zero is not valid in the dictionary + // these transforms are not actually used by testbeam GeoModel anyway + ATH_MSG_WARNING( "Ident for unknown detector type in " << *iobj ); + det=1;bec=0;layer=0;ring=0;sector=0;side=0; + } + Amg::Vector3D shift=trans.translation(); + Amg::RotationMatrix3D rot=trans.rotation(); + dx=shift.x(); + dy=shift.y(); + dz=shift.z(); + double alpha, beta, gamma; + extractAlphaBetaGamma(trans, alpha, beta, gamma); + + ATH_MSG_WARNING("THIS NEEDS TO BE CHECKED --- InDetAlignDBTool.cxx:647"); + Amg::Vector3D ea = rot.eulerAngles(2, 0, 2); + + phi= ea[0]; + theta=ea[1]; + psi=ea[2]; + ++ntrans; + if (ntuple) { + nt_dettype=det; + nt_bec=bec; + nt_layer=layer; + nt_ring=ring; + nt_sector=sector; + nt_level=3; + // derive level by looking for objects ending in /ID, SCT or PIX + std::string r3=iobj->substr(iobj->size()-3,3); + if (r3=="/ID") nt_level=1; + if (r3=="SCT" || r3=="PIX") nt_level=2; + nt_side=side; + nt_xofs=dx; + nt_yofs=dy; + nt_zofs=dz; + nt_phi=phi; + nt_theta=theta; + nt_psi=psi; + nt_alpha=alpha; + nt_beta=beta; + nt_gamma=gamma; + if (StatusCode::SUCCESS!=nt->write()) ATH_MSG_ERROR("Problem filling ntuple 9002" ); + } else { + *outfile << "2 " << det << " " << 2*bec << " " << layer << " " << sector << + " " << ring << " " << side << " " << dx << " " << dy << " " + << dz << " " << alpha/CLHEP::mrad << " " << beta/CLHEP::mrad << " " << gamma/CLHEP::mrad << std::endl; + } + } + } else { + ATH_MSG_ERROR("Cannot find AlignableTransform for key " + << *iobj ); + } + //if (!ntuple) + // end of transform marked by line of zeros + // *outfile << "0 0 0 0 0 0 0 0 0 0 0 0" << std::endl; + } + if (ntuple) { + } else { + outfile->close(); + delete outfile; + } + ATH_MSG_DEBUG("Written " << nobj << " AlignableTransform objects" << " with " << ntrans << " transforms to text file" ); +} + + +void InDetAlignDBTool::readTextFile(const std::string file) const { + if (par_oldTextFile) return readOldTextFile(file); + + ATH_MSG_DEBUG("readTextFile - set alignment constants from text file: " << file ); + std::ifstream infile; + infile.open(file.c_str()); + if (!infile) { + ATH_MSG_ERROR("Error opening file " << file ); + return; + } + + // loop over lines in file + int nobj=0; + int ntrans=0; + + std::string channelName; // Channel name + const AlignableTransform* pat = 0; + + while (infile) { + std::string tmpline; + std::getline(infile, tmpline); + if (!infile) break; + + // Skip comment line + if ((tmpline.substr(0,2) == "//") || (tmpline.substr(0,1) == "#")) continue; + + std::istringstream instring(tmpline); + std::string tmpstr; + instring >> tmpstr; + + // Skip blank line + if (tmpstr.empty()) continue; + + if (tmpstr[0] == '/') { + // Its a valid channel name + channelName = tmpstr; + ATH_MSG_DEBUG("Read in AlignableTransform data, key " << channelName ); + // find the AlignableTransform with this key + pat = 0; + if (!(pat=cgetTransPtr(channelName))) { + ATH_MSG_ERROR("Cannot find AlignableTransform object for key" + << channelName ); + } else { + nobj++; + } + } else { + // Its a data line + + if (!pat) { + // If pat = 0, then either no channel name was specified or it could not be found. + ATH_MSG_ERROR("No channel specified. Skipping input " ); + + } else { + // normal data + std::istringstream datastream(tmpline); + + int subsystem,dettype,bec,layer,phiModule,etaModule,side; + float dx,dy,dz,alpha,beta,gamma; + datastream >> subsystem >> dettype >> bec >> layer >> phiModule >> etaModule >> side >> dx + >> dy >> dz >> alpha >> beta >> gamma; + + if (datastream.fail()) { + ATH_MSG_ERROR("Error in input" ); + } else { + alpha *= CLHEP::mrad; + beta *= CLHEP::mrad; + gamma *= CLHEP::mrad; + + // construct identifier + Identifier ident=Identifier(); + if (dettype==1) { + ident=m_pixid->wafer_id(bec,layer,phiModule,etaModule); + } else if (dettype==2) { + ident=m_sctid->wafer_id(bec,layer,phiModule,etaModule,side); + } else { + ATH_MSG_ERROR("Cannot construct identifier for dettype " + << dettype ); + } + if (!ident.is_valid()) { + ATH_MSG_ERROR("Error in identifier : " << + " [" << subsystem << "," << dettype << "," << bec << "," << layer << "," << + phiModule << "," << etaModule << "," << side << "] key " << channelName << + " shift [" << dx << "," << dy << "," << dz << "]" ); + } else { + + // construct new transform + // Order of rotations is defined as around z, then y, then x. + // For small angles it doesn't really matter though. + Amg::Translation3D newtranslation(dx,dy,dz); + Amg::Transform3D newtrans = newtranslation * Amg::RotationMatrix3D::Identity(); + newtrans *= Amg::AngleAxis3D(gamma, Amg::Vector3D(0.,0.,1.)); + newtrans *= Amg::AngleAxis3D(beta, Amg::Vector3D(0.,1.,0.)); + newtrans *= Amg::AngleAxis3D(alpha, Amg::Vector3D(1.,0.,0.)); + + + + // find pointer to existing transform, currently missing write access + // via findIdent, so have to search manually + AlignableTransform* pat2; + pat2=const_cast<AlignableTransform*>(pat); + AlignableTransform::AlignTransMem_itr itr=pat2->mbegin(); + while ((itr->identify()!=ident) && (itr!=pat2->mend())) ++itr; + if (itr!=pat2->mend()) { + ++ntrans; + itr->setTransform( Amg::EigenTransformToCLHEP(newtrans) ); + ATH_MSG_VERBOSE ( "Set transform done"); + } else { + ATH_MSG_WARNING("Cannot find existing transform for"); + } + // Can uses either id helper + ATH_MSG_DEBUG(m_pixid->show_to_string(ident) << " key " << channelName << + " shift [" << dx << "," << dy << "," << dz << "]" ); + } + } // end if (datastream.fail()) + } // end if (!pat) + } // end if (tmpstr[0] == '/') + } // end while (infile) + + infile.close(); + ATH_MSG_DEBUG( "Read " << nobj << " objects from file with " << ntrans << " transforms" ); +} + +// Old text format: +// No subsystem number (always 2 for indet and not really needed, but for completeness its in the new format). +// ring and sector number swapped. +// 0 0 0 ... lines to terminate +// No comments lines. +void InDetAlignDBTool::readOldTextFile(const std::string file) const { + ATH_MSG_DEBUG( "readTextFile - set alignment constants from text file: " << file ); + std::ifstream infile; + infile.open(file.c_str()); + // loop over lines in file + int nobj=0; + int ntrans=0; + std::string atname; + while (infile >> atname) { + ++nobj; + + ATH_MSG_DEBUG("Read in AlignableTransform data, key " << atname ); + ATH_MSG_DEBUG( "(You are using the the old format)" << file ); + + // find the AlignableTransform with this key + const AlignableTransform* pat; + if (!(pat=cgetTransPtr(atname))) + ATH_MSG_ERROR("Cannot find AlignableTransform object for key" + << atname ); + int dettype,bec,layer,ring,sector,side; + float dx,dy,dz,phi,theta,psi; + while ((infile >> dettype >> bec >> layer >> ring >> sector >> side >> dx + >> dy >> dz >> phi >> theta >> psi) && (dettype!=0)) { + // construct identifier + if (pat!=0) { + Identifier ident=Identifier(); + if (dettype==1) { + ident=m_pixid->wafer_id(bec,layer,sector,ring); + } else if (dettype==2) { + ident=m_sctid->wafer_id(bec,layer,sector,ring,side); + } else { + ATH_MSG_ERROR("Cannot construct identifier for dettype " + << dettype ); + } + // construct new transform + Amg::Translation3D newtranslation(dx,dy,dz); + Amg::Transform3D newtrans = newtranslation * Amg::RotationMatrix3D::Identity(); + newtrans *= Amg::AngleAxis3D(psi, Amg::Vector3D(0.,0.,1.)); + newtrans *= Amg::AngleAxis3D(theta, Amg::Vector3D(0.,1.,0.)); + newtrans *= Amg::AngleAxis3D(phi, Amg::Vector3D(1.,0.,0.)); + + + // find pointer to existing transform, currently missing write access + // via findIdent, so have to search manually + AlignableTransform* pat2; + pat2=const_cast<AlignableTransform*>(pat); + AlignableTransform::AlignTransMem_itr itr=pat2->mbegin(); + while ((itr->identify()!=ident) && (itr!=pat2->mend())) ++itr; + if (itr!=pat2->mend()) { + ++ntrans; + itr->setTransform( Amg::EigenTransformToCLHEP(newtrans) ); + ATH_MSG_VERBOSE ("Set transform done"); + } else { + ATH_MSG_WARNING("Cannot find existing transform for"); + } + ATH_MSG_DEBUG( " [" << dettype << "," << bec << "," << layer << "," << + ring << "," << sector << "," << side << "] key " << atname << + " shift [" << dx << "," << dy << "," << dz << "]" ); + } + } //end tranform loop + } // end input loop + infile.close(); + ATH_MSG_DEBUG("Read " << nobj << " objects from file with " << ntrans << " transforms" ); +} + +void InDetAlignDBTool::readNtuple(const std::string file) const { + ATH_MSG_DEBUG("readNtuple - set alignment constants from ntuple path: " << file ); + INTupleSvc* ntsvc; + if (StatusCode::SUCCESS!=service("NTupleSvc",ntsvc,true)) + ATH_MSG_ERROR("Cannot find NTupleSvc" ); + const std::string path=file+"/9002"; + NTuplePtr nt(ntsvc,path); + if (nt) { + StatusCode sc; + sc=nt->item( "MODPROP/DetType",nt_dettype); + sc=nt->item("MODPROP/Bec",nt_bec); + sc=nt->item("MODPROP/Layer",nt_layer); + sc=nt->item("MODPROP/Ring",nt_ring); + sc=nt->item("MODPROP/Sector",nt_sector); + sc=nt->item("MODPROP/Side",nt_side); + sc=nt->item("MODPROP/Level",nt_level); + sc=nt->item("MODPROP/Xofs",nt_xofs); + sc=nt->item("MODPROP/Yofs",nt_yofs); + sc=nt->item("MODPROP/Zofs",nt_zofs); + sc=nt->item("MODPROP/Phi",nt_phi); + sc=nt->item("MODPROP/Theta",nt_theta); + sc=nt->item("MODPROP/Psi",nt_psi); + //sc=nt->item("MODPROP/Alpha",nt_alpha); + //sc=nt->item("MODPROP/Beta",nt_beta); + //sc=nt->item("MODPROP/Gamma",nt_gamma); + if (sc!=StatusCode::SUCCESS) ATH_MSG_ERROR( + "Error booking ntuple 9002 contents" ); + int ntrans=0; + while (nt->read().isSuccess()) { + Identifier ident=Identifier(); + if (nt_dettype==1) { + ident=m_pixid->wafer_id(2*nt_bec,nt_layer,nt_sector,nt_ring); + } else if (nt_dettype==2) { + ident=m_sctid->wafer_id(2*nt_bec,nt_layer,nt_sector,nt_ring,nt_side); + } else { + ATH_MSG_ERROR("Cannot construct identifier for dettype " + << nt_dettype ); + } + + Amg::Translation3D newtranslation(nt_xofs,nt_yofs,nt_zofs); + Amg::Transform3D newtrans = newtranslation * Amg::RotationMatrix3D::Identity(); + newtrans *= Amg::AngleAxis3D(nt_psi, Amg::Vector3D(0.,0.,1.)); + newtrans *= Amg::AngleAxis3D(nt_theta, Amg::Vector3D(0.,1.,0.)); + newtrans *= Amg::AngleAxis3D(nt_phi, Amg::Vector3D(1.,0.,0.)); + + setTrans(ident,nt_level,newtrans); + ++ntrans; + } + ATH_MSG_DEBUG( "Read " << ntrans << " transforms from ntuple"); + } else { + ATH_MSG_ERROR( "Problem opening ntuple at path " << path ); + } +} + +bool InDetAlignDBTool::setTrans(const Identifier& ident, const int level, + const Amg::Transform3D& trans) const { + // find transform key, then set appropriate transform + // do storegate const retrieve, then cast to allow update of locked data + std::string key=dirkey(ident,level); + const AlignableTransform* pat; + AlignableTransform* pat2; + bool result=false; + if ((pat=cgetTransPtr(key))) { + pat2=const_cast<AlignableTransform*>(pat); + if (pat2!=0) { + result=pat2->update(ident, Amg::EigenTransformToCLHEP(trans) ); + if (!result) ATH_MSG_ERROR( + "Attempt to set non-existant transform" ); + } else { + ATH_MSG_ERROR("setTrans: cast fails for key " << key ); + } + } else { + ATH_MSG_ERROR( + "setTrans: cannot retrieve AlignableTransform for key" << key ); + } + return result; +} + +bool InDetAlignDBTool::setTrans(const Identifier& ident, const int level, + const Amg::Vector3D& translate, double alpha, double beta, double gamma) const +{ + + Amg::Translation3D newtranslation(translate); + Amg::Transform3D newtrans = newtranslation * Amg::RotationMatrix3D::Identity(); + newtrans *= Amg::AngleAxis3D(gamma, Amg::Vector3D(0.,0.,1.)); + newtrans *= Amg::AngleAxis3D(beta, Amg::Vector3D(0.,1.,0.)); + newtrans *= Amg::AngleAxis3D(alpha, Amg::Vector3D(1.,0.,0.)); + + return setTrans(ident, level, newtrans); +} + + +bool InDetAlignDBTool::tweakTrans(const Identifier& ident, const int level, + const Amg::Transform3D& trans) const { + + // find transform key, then set appropriate transform + std::string key=dirkey(ident,level); + const AlignableTransform* pat; + AlignableTransform* pat2; + bool result=false; + if ((pat=cgetTransPtr(key))) { + pat2=const_cast<AlignableTransform*>(pat); + if (pat2!=0) { + result=pat2->tweak(ident,Amg::EigenTransformToCLHEP(trans)); + if (!result) ATH_MSG_ERROR( + "Attempt to tweak non-existant transform" ); + } else { + ATH_MSG_ERROR("tweakTrans: cast fails for key " << key ); + } + } else { + ATH_MSG_ERROR( + "tweakTrans: cannot retrieve AlignableTransform for key" << key ); + } + return result; +} + +bool InDetAlignDBTool::tweakTrans(const Identifier& ident, const int level, + const Amg::Vector3D& translate, double alpha, double beta, double gamma) const +{ + + Amg::Translation3D newtranslation(translate); + Amg::Transform3D newtrans = newtranslation * Amg::RotationMatrix3D::Identity(); + newtrans *= Amg::AngleAxis3D(gamma, Amg::Vector3D(0.,0.,1.)); + newtrans *= Amg::AngleAxis3D(beta, Amg::Vector3D(0.,1.,0.)); + newtrans *= Amg::AngleAxis3D(alpha, Amg::Vector3D(1.,0.,0.)); + + return tweakTrans(ident, level, newtrans); +} + +/** convert L3 module identifier to L1 or L2 */ +Identifier InDetAlignDBTool::getL1L2fromL3Identifier( const Identifier& ident + , const int& level + ) const { + if( level == 3 ) return ident ; //!< no translation needed + /** check whether PIX */ + if( m_pixid->is_pixel(ident) ) { + if( level == 1 ) { + return m_pixid->wafer_id( 0, 0, 0, 0 ) ; //!< Whole pixel det. at L1 + } + if( level == 2 ) { + int barrel_ec = m_pixid->barrel_ec( ident ) ; + int layer_disk = m_pixid->layer_disk( ident ) ; + return m_pixid->wafer_id( barrel_ec, layer_disk, 0, 0 ) ; + } + } + /** check whether SCT */ + if( m_sctid->is_sct(ident) ) { + if( level == 1 ) { + int barrel_ec = m_sctid->barrel_ec( ident ) ; + return m_sctid->wafer_id( barrel_ec, 0, 0, 0, 0 ) ; //!< barrel + 2 x EC at L1 + } + if( level == 2 ) { + int barrel_ec = m_sctid->barrel_ec( ident ) ; + int layer_disk = m_sctid->layer_disk( ident ) ; + return m_sctid->wafer_id( barrel_ec, layer_disk, 0, 0, 0 ) ; + } + } + return ident ; //!< take care of the case where level != 1,2,3 or ident neither pix nor sct +} + +/** get cumulative L1, L2, L3 trafo for (L3-) module */ +Amg::Transform3D InDetAlignDBTool::getTransL123( const Identifier& ident ) const { + + Amg::Transform3D result ; + InDetDD::SiDetectorElement* element = m_pixman->getDetectorElement( ident ) ; + if( !element ) { + element = m_sctman->getDetectorElement( ident ) ; + } + if( !element ) { + ATH_MSG_ERROR("getTransL123(): Module not found in PIX or SCT!" ); + return result ; + } + Amg::Transform3D trfL1 = getTrans( ident, 1 ) ; + Amg::Transform3D trfL2 = getTrans( ident, 2 ) ; + Amg::Transform3D trfL3 = getTrans( ident, 3 ) ; + ATH_MSG_FATAL("Code needs to corrected otherwise you will get nonsensical results-- IndetAlignDBTool:2060"); + const Amg::Transform3D trfNominal ; //= element->defModuleTransform() ; + result = trfNominal.inverse() * trfL1 * trfL2 * trfNominal * trfL3 ; + return result ; +} + +/** return value of particular transform specified by identifier and level + calculates L1 and L2 identifiers automatically by getL1L2fromL3Identifier + if L3 identifier passed */ +Amg::Transform3D InDetAlignDBTool::getTrans(const Identifier& ident, + const int level) const { + const Identifier identifier = getL1L2fromL3Identifier( ident, level ) ; + Amg::Transform3D result; + const std::string key=dirkey(identifier,level); + const AlignableTransform* pat; + if ((pat=cgetTransPtr(key))) { + AlignableTransform::AlignTransMem_citr itr=pat->findIdent(identifier); + if (itr!=pat->end()) result= Amg::CLHEPTransformToEigen(itr->transform()); + } + return result; +} + +StatusCode InDetAlignDBTool::outputObjs() const { + + ATH_MSG_DEBUG( "Output AlignableTranform objects to stream" << par_condstream ); + // get the AthenaOutputStream tool + IAthenaOutputStreamTool* optool; + + if (StatusCode::SUCCESS!=p_toolsvc->retrieveTool("AthenaPoolOutputStreamTool",par_condstream,optool)) { + ATH_MSG_ERROR("Cannot get AthenaPoolOutputStream tool" ); + return StatusCode::FAILURE; + } + + if (StatusCode::SUCCESS!=optool->connectOutput()) { + ATH_MSG_ERROR("Could not connect stream to output" ); + return StatusCode::FAILURE; + } + // construct list of objects to be written out, either + // AlignableTransformContainer or several of AlignableTransforms + int npairs=m_alignobjs.size(); + if (par_newdb) npairs=1; + IAthenaOutputStreamTool::TypeKeyPairs typekeys(npairs); + if (par_newdb) { + typekeys[0]= + IAthenaOutputStreamTool::TypeKeyPair("AlignableTransformContainer", + par_dbroot); + if (!(detStore()->contains<AlignableTransformContainer>(par_dbroot))) + ATH_MSG_ERROR( + "Expected " << par_dbroot << " object not found" ); + } else { + for (unsigned int i=0;i<m_alignobjs.size();++i) { + typekeys[i]=IAthenaOutputStreamTool::TypeKeyPair("AlignableTransform", + m_alignobjs[i]); + if (!(detStore()->contains<AlignableTransform>(m_alignobjs[i]))) + ATH_MSG_ERROR("Expected " << m_alignobjs[i] << " object not found" ); + } + } + // write objects to stream + if (StatusCode::SUCCESS!=optool->streamObjects(typekeys)) { + ATH_MSG_ERROR("Could not stream output objects" ); + return StatusCode::FAILURE; + } + // commit output + if (StatusCode::SUCCESS!=optool->commitOutput()) { + ATH_MSG_ERROR("Could not commit output" ); + } + ATH_MSG_DEBUG( "Written " << typekeys.size() << " objects to stream " << par_condstream); + return StatusCode::SUCCESS; +} + +void InDetAlignDBTool::fillDB(const std::string tag, + const unsigned int run1, const unsigned int event1, + const unsigned int run2, const unsigned int event2) const { + + ATH_MSG_DEBUG( "fillDB: Data tag " << tag ); + ATH_MSG_DEBUG( "Run/evt1 [" << run1 << "," << event1 << "]" ); + ATH_MSG_DEBUG("Run/evt2 [" << run2 << "," << event2 << "]" ); + + // get pointer to registration svc + IIOVRegistrationSvc* regsvc; + if (StatusCode::SUCCESS!= + service("IOVRegistrationSvc",regsvc)) { + ATH_MSG_FATAL( "IOVRegistrationSvc not found" ); + return; + } + // loop over all AlignableTransform objects created earlier and save them + int nobj=0; + if (par_newdb) { + if (StatusCode::SUCCESS==regsvc->registerIOV( + "AlignableTransformContainer",par_dbroot,tag,run1,run2,event1,event2)) { + ATH_MSG_DEBUG( "Stored AlignableTransform object " << par_dbroot ); + ++nobj; + } else { + ATH_MSG_ERROR("Failed (registerIOV) to store object " << par_dbroot ); + } + } else { + // old way - register all objects separately + for (std::vector<std::string>::const_iterator iobj=m_alignobjs.begin(); + iobj!=m_alignobjs.end();++iobj) { + if (StatusCode::SUCCESS==regsvc->registerIOV("AlignableTransform", + *iobj,tag,run1,run2,event1,event2)) { + ATH_MSG_DEBUG( "Stored AlignableTransform object " << *iobj ); + ++nobj; + } else { + ATH_MSG_ERROR("Failed (registerIOV) to store object " << *iobj ); + } + } + } + ATH_MSG_DEBUG( " Written " << nobj << " AlignableTransform objects to conditions database" ); +} + +void InDetAlignDBTool::printDB(const int level) const { + + ATH_MSG_DEBUG("Printout InDetAlign database contents, detail level" << level ); + + for (std::vector<std::string>::const_iterator iobj=m_alignobjs.begin(); + iobj!=m_alignobjs.end();++iobj) { + const AlignableTransform* pat; + if ((pat=cgetTransPtr(*iobj))) { + ATH_MSG_DEBUG( "AlignableTransform object " << *iobj ); + int nobj=0; + for (AlignableTransform::AlignTransMem_citr cit=pat->begin(); + cit!=pat->end();++cit) { + const Identifier& ident=cit->identify(); + const Amg::Transform3D& trans= Amg::CLHEPTransformToEigen( cit->transform() ); + Amg::Vector3D shift=trans.translation(); + //Amg::RotationMatrix3D rot=trans.rotation(); + int det,bec,layer,ring,sector,side; + if (idToDetSet(ident,det,bec,layer,ring,sector,side)) { + if (level>1) { + double alpha, beta, gamma; + extractAlphaBetaGamma(trans, alpha, beta, gamma); + ATH_MSG_DEBUG( "ID [" << det << "," << bec << "," << layer << + "," << ring << "," << sector << "," << side << "] Trans:(" << + shift.x() << "," << shift.y() << "," << shift.z() << ") Rot:{" + << alpha << "," << beta << "," << gamma << "}"); + } + ++nobj; + } else { + ATH_MSG_ERROR("Unknown identifier in AlignableTransform" ); + } + } + ATH_MSG_DEBUG( "Object contains " << nobj << " transforms" ); + } else { + ATH_MSG_ERROR( "AlignableTransform " << *iobj << " not found" ); + } + } +} + +// ========================================== + +AlignableTransform* InDetAlignDBTool::getTransPtr(const std::string key) + const { + // look in collection to retrieve pointer to AlignableTransform object of + // given key and return it, return 0 if not collection or key value not found + AlignableTransformContainer* patc; + AlignableTransform* pat=0; + if (par_newdb) { + if (StatusCode::SUCCESS==detStore()->retrieve(patc,par_dbroot )) { + for (DataVector<AlignableTransform>::iterator dva=patc->begin(); + dva!=patc->end();++dva) { + if ((*dva)->tag()==key) { + pat=*dva; + break; + } + } + } + } else { + if (StatusCode::SUCCESS!=detStore()->retrieve(pat,key)) pat=0; + } + return pat; +} + +const AlignableTransform* InDetAlignDBTool::cgetTransPtr(const std::string key) + const { + // look in collection to retrieve pointer to AlignableTransform object of + // given key and return it, return 0 if not collection or key value not found + // const version + const AlignableTransformContainer* patc; + const AlignableTransform* pat=0; + if (par_newdb) { + if (StatusCode::SUCCESS==detStore()->retrieve(patc,par_dbroot )) { + for (DataVector<AlignableTransform>::const_iterator dva=patc->begin(); + dva!=patc->end();++dva) { + if ((*dva)->tag()==key) { + pat=*dva; + break; + } + } + } + } else { + if (StatusCode::SUCCESS!=detStore()->retrieve(pat,key)) pat=0; + } + return pat; +} + +void InDetAlignDBTool::fakeGeom(const int nbpix, const int necpix, + const int nbsct, const int necsct) { + // set alignment keys for fake geometry with given numbers of + // barrel/endcap PIX/SCT layers + // this code is somewhat fragile, trying to reproduce the order of + // keys in the same way that GeoModel returns them + // will not work for layouts with missing middle pixel layer + + int ichan3=200; + // level 1 object - ID + m_alignobjs.push_back(dirkey(1,0,0,1)); + m_alignchans.push_back(0); + // level 2 objects - pixel + if (nbpix!=0 || necpix!=0) { + m_alignobjs.push_back(dirkey(1,0,0,2)); + m_alignchans.push_back(100); + } + // level 3 objects - pixel + // objects done in this order to get COOL channels correct + // endcap A pixel + for (int i=0;i<necpix;++i) { + m_alignobjs.push_back(dirkey(1,-1,i,3)); + m_alignchans.push_back(ichan3++); + } + // barrel pixel + for (int i=0;i<nbpix;++i) { + m_alignobjs.push_back(dirkey(1,0,i,3)); + m_alignchans.push_back(ichan3++); + } + // endcap C pixel + for (int i=0;i<necpix;++i) { + m_alignobjs.push_back(dirkey(1,1,i,3)); + m_alignchans.push_back(ichan3++); + } + // level 2 objects - SCT + if (nbsct!=0 || necsct!=0) { + m_alignobjs.push_back(dirkey(2,0,0,2)); + m_alignchans.push_back(101); + } + // level 3 objects - SCT + // endcap A SCT + for (int i=0;i<necsct;++i) { + m_alignobjs.push_back(dirkey(2,-1,i,3)); + m_alignchans.push_back(ichan3++); + } + // barrel SCT + for (int i=0;i<nbsct;++i) { + m_alignobjs.push_back(dirkey(2,0,i,3)); + m_alignchans.push_back(ichan3++); + } + // endcap C SCT + for (int i=0;i<necsct;++i) { + m_alignobjs.push_back(dirkey(2,1,i,3)); + m_alignchans.push_back(ichan3++); + } +} + +void InDetAlignDBTool::sortTrans() const { + // loop through all the AlignableTransform objects and sort them + + ATH_MSG_DEBUG( "Sorting all AlignableTransforms in TDS" ); + AlignableTransform* pat; + // use cget and a const cast to allow containers that have been read in + // (and hence are locked by StoreGate) to be sorted + for (unsigned int i=0;i<m_alignobjs.size();++i) + if ((pat=const_cast<AlignableTransform*>(cgetTransPtr(m_alignobjs[i])))) pat->sortv(); +} + +void InDetAlignDBTool::extractAlphaBetaGamma(const Amg::Transform3D & trans, + double& alpha, double& beta, double &gamma) const +{ + double siny = trans(0,2); + beta = asin(siny); + // Check if cosy = 0. This requires special treatment. + // can check either element (1,2),(2,2) both equal zero + // or (0,1) and (0,0) + // Probably not likely it will be exactly 0 and may still + // have some problems when very close to zero. We mostly + // deal with small rotations so its not too important. + if ((trans(1,2) == 0) && (trans(2,2) == 0)) { + // alpha and gamma are degenerate. We arbitrarily choose + // gamma = 0. + gamma = 0; + alpha = atan2(trans(1,1),trans(2,1)); + } else { + alpha = atan2(-trans(1,2),trans(2,2)); + gamma = atan2(-trans(0,1),trans(0,0)); + if (alpha == 0) alpha = 0; // convert -0 to 0 + if (gamma == 0) gamma = 0; // convert -0 to 0 + } +} diff --git a/InnerDetector/InDetAlignTools/InDetAlignGenTools/src/InDetAlignFillSiCluster.cxx b/InnerDetector/InDetAlignTools/InDetAlignGenTools/src/InDetAlignFillSiCluster.cxx new file mode 100755 index 0000000000000000000000000000000000000000..6fbd2d10cd267c99a85efdf3125747c9b394e4b9 --- /dev/null +++ b/InnerDetector/InDetAlignTools/InDetAlignGenTools/src/InDetAlignFillSiCluster.cxx @@ -0,0 +1,387 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +////////////////////////////////////////////////////////////////////////// +// ================================================ +// InDetAlignFillSiCluster +// ================================================ +// +// InDetAlignFillSiCluster.cxx +// Source file for InDetAlignFillSiCluster +// +// Carlos Escobar, started 08/03/2008 +// +// AthAlgTool to fill silicon cluster information in a ntuple + +#include "GaudiKernel/NTuple.h" +#include "GaudiKernel/INTupleSvc.h" +#include "GaudiKernel/SmartDataPtr.h" + +#include "TrkEventPrimitives/ParamDefs.h" + +#include "InDetIdentifier/PixelID.h" +#include "InDetIdentifier/SCT_ID.h" +#include "InDetPrepRawData/PixelClusterContainer.h" +#include "InDetPrepRawData/SCT_ClusterContainer.h" +#include "InDetReadoutGeometry/SCT_DetectorManager.h" +#include "InDetReadoutGeometry/PixelDetectorManager.h" + +#include "InDetAlignGenTools/InDetAlignFillSiCluster.h" + +#include <string> + +static const int maxclusters = 24000; +static const int maxclsize = 10; + +//===================================================================== +// InDetAlignFillSiCluster() +//===================================================================== +InDetAlignFillSiCluster::InDetAlignFillSiCluster(const std::string& type, + const std::string& name, + const IInterface* parent) + : AthAlgTool(type,name,parent), + m_ntupleSvc(0), + m_pixelid(0), + m_sctID(0), + m_Pixel_clcontainer(0), + m_Sct_clcontainer(0) +{ + declareInterface<IInDetAlignFillSiCluster>(this); + declareProperty("Pixel_SiClusterContainerName", m_Pixel_SiClustersName="PixelClusters"); + declareProperty("SCT_SiClusterContainerName", m_Sct_SiClustersName="SCT_Clusters"); + + // Ntuple + declareProperty("NtupleName", m_ntupleName="/NTUPLES/GLOBFILE"); +} + + +//===================================================================== +// ~FillSiCluster() +//===================================================================== +InDetAlignFillSiCluster::~InDetAlignFillSiCluster() {} + + +//===================================================================== +// initialize() +//===================================================================== +StatusCode InDetAlignFillSiCluster::initialize() { + + ATH_MSG_DEBUG("In Initialize() of FillSiCluster()"); + + // ID Helper + + if (detStore()->retrieve(m_sctID, "SCT_ID").isFailure()){ + msg(MSG::FATAL) << "Could not get SCT ID helper" << endreq; + return StatusCode::FAILURE; + } + else if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "SCT ID is : "<< m_sctID <<endreq ; + + if (detStore()->retrieve(m_pixelid, "PixelID").isFailure()){ + msg(MSG::FATAL) << "Could not get PIXEL ID helper" << endreq; + return StatusCode::FAILURE; + } + else if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "Pixel ID is : " << m_pixelid << endreq; + + // get SCTDetectorManager + const InDetDD::SCT_DetectorManager* mgr; + if (detStore()->retrieve(mgr, "SCT").isFailure()) { + ATH_MSG_FATAL("Could not get SCT_DetectorManager!"); + return StatusCode::FAILURE; + } + else ATH_MSG_DEBUG ("Manager found!"); + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) <<"SCT ID is : "<< m_sctID << endreq; + + // get PixelDetectorManager + const InDetDD::PixelDetectorManager* pixelmgr; + if (detStore()->retrieve(pixelmgr, "Pixel").isFailure()) { + ATH_MSG_FATAL("Could not get PixelDetectorManager!"); + return StatusCode::FAILURE; + } + else ATH_MSG_DEBUG ("Pixel_DetectorManager found!"); + + + // retrieve the NTuple Service + if (StatusCode::SUCCESS != service("NTupleSvc", m_ntupleSvc)) { + ATH_MSG_FATAL ("NTupleSvc service not found!"); + return StatusCode::FAILURE; + } + + bookNtuple(); + + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "Initialize() of FillSiCluster successful" << endreq; + return StatusCode::SUCCESS; +} + + +//===================================================================== +// finalize() +//===================================================================== +StatusCode InDetAlignFillSiCluster::finalize() { + + if (msgLvl(MSG::DEBUG)) { + msg(MSG::DEBUG) << "Finalize() of FillSiCluster" << endreq; + msg(MSG::DEBUG) << "________________________________________________________" << endreq; + msg(MSG::DEBUG) << endreq; + msg(MSG::DEBUG) << " InDetAlignFillSiCluster Summary: " << endreq; + msg(MSG::DEBUG) << "________________________________________________________" << endreq; + msg(MSG::DEBUG) << endreq; + } + + return StatusCode::SUCCESS; +} + + +//===================================================================== +// FillSiCluster() +//===================================================================== +StatusCode InDetAlignFillSiCluster::FillSiCluster() { + + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "In FillSiCluster()" << endreq; + + StatusCode sc; + + if(RetrieveSCTSiClusters() != StatusCode::FAILURE) + FillSCTSiNtuple(); + + if(RetrievePixelSiClusters() != StatusCode::FAILURE) + FillPixelSiNtuple(); + + std::string nt0id = m_ntupleName + "/SiCluster"; + sc = m_ntupleSvc->writeRecord(nt0id); + if (sc.isFailure()) if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "Could not write " << nt0id << "!" << endreq; + + return StatusCode::SUCCESS; +} + + +//===================================================================== +// bookNtuple() +//===================================================================== +void InDetAlignFillSiCluster::bookNtuple() { + + ATH_MSG_DEBUG ("Booking Trk::Track Info..."); + + NTupleFilePtr file1(m_ntupleSvc, m_ntupleName); + std::string nt0id = m_ntupleName + "/SiCluster"; + std::string comments = "Silicon Cluster Information"; + + StatusCode sc; + + NTuplePtr nt0(m_ntupleSvc, nt0id); + if (nt0) ATH_MSG_DEBUG ("Ntuple is already booked"); + else { + ATH_MSG_DEBUG ("Attempting to book general ntuple"); + nt0 = m_ntupleSvc->book(nt0id,CLID_ColumnWiseTuple,comments); + + if (nt0) { + + // --------------------------------------------------------------------- + // Pixel Clusters + sc = nt0->addItem("pixel_clus_nclusters", m_pixel_nclusters, 0, maxclusters); + sc = nt0->addItem("pixel_clus_clx", m_pixel_nclusters, m_pixel_clx); + sc = nt0->addItem("pixel_clus_cly", m_pixel_nclusters, m_pixel_cly); + sc = nt0->addItem("pixel_clus_clz", m_pixel_nclusters, m_pixel_clz); + sc = nt0->addItem("pixel_clus_LocX", m_pixel_nclusters, m_pixel_LocX); + sc = nt0->addItem("pixel_clus_LocY", m_pixel_nclusters, m_pixel_LocY); + sc = nt0->addItem("pixel_clus_groupsize", m_pixel_nclusters, m_pixel_groupsize); + sc = nt0->addItem("pixel_clus_layer", m_pixel_nclusters, m_pixel_layer); + sc = nt0->addItem("pixel_clus_phi", m_pixel_nclusters, m_pixel_phi); + sc = nt0->addItem("pixel/pixel_clus_row",m_pixel_nclusters, m_pixel_clrow, maxclsize); + sc = nt0->addItem("pixel/pixel_clus_col",m_pixel_nclusters, m_pixel_clcol, maxclsize); + m_pixel_nclusters = 0; + // --------------------------------------------------------------------- + + // ---------------------------------------------------------------------- + // SCT Clusters + sc = nt0->addItem("sct_clus_nclusters", m_sct_nclusters, 0, maxclusters); + sc = nt0->addItem("sct_clus_x", m_sct_nclusters, m_sct_clx); + sc = nt0->addItem("sct_clus_y", m_sct_nclusters, m_sct_cly); + sc = nt0->addItem("sct_clus_z", m_sct_nclusters, m_sct_clz); + sc = nt0->addItem("sct_clus_groupsize", m_sct_nclusters, m_sct_groupsize); + sc = nt0->addItem("sct_clus_layer", m_sct_nclusters, m_sct_layer); + sc = nt0->addItem("sct_clus_eta", m_sct_nclusters, m_sct_eta); + sc = nt0->addItem("sct_clus_phi", m_sct_nclusters, m_sct_phi); + sc = nt0->addItem("sct_clus_side", m_sct_nclusters, m_sct_side); + m_sct_nclusters = 0; + // ---------------------------------------------------------------------- + + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "Ntuple " << nt0id << " has been booked successfully! " << endreq; + } + else + ATH_MSG_ERROR ("Error booking ntuple"); + } + +} + + +//===================================================================== +// RetrieveSCTSiClusters() +//===================================================================== +StatusCode InDetAlignFillSiCluster::RetrieveSCTSiClusters() { + + ATH_MSG_DEBUG ("In RetrieveSCTSiClusters()"); + StatusCode sc; + + // get clusters from TDS + sc = evtStore()->retrieve(m_Sct_clcontainer, m_Sct_SiClustersName); + if (sc.isFailure() || !m_Sct_clcontainer){ + ATH_MSG_DEBUG ("SCT Cluster container for SCT not found"); + return StatusCode::FAILURE; + } + else ATH_MSG_DEBUG ("SCT Cluster container for SCT found"); + + return sc; +} + + +//===================================================================== +// RetrievePixelSiClusters() +//===================================================================== +StatusCode InDetAlignFillSiCluster::RetrievePixelSiClusters() { + + ATH_MSG_DEBUG ("In RetrievePixelSiClusters()"); + StatusCode sc; + + // get pixel clusters from TDS + sc = evtStore()->retrieve(m_Pixel_clcontainer, m_Pixel_SiClustersName); + if (sc.isFailure() || !m_Pixel_clcontainer){ + ATH_MSG_DEBUG ("Pixel Cluster container for Pixels not found"); + return StatusCode::FAILURE; + } + else ATH_MSG_DEBUG ("Pixel Cluster container for Pixels found"); + + return sc; +} + +//===================================================================== +// FillSCTSiNtuple() +//===================================================================== +void InDetAlignFillSiCluster::FillSCTSiNtuple() { + + ATH_MSG_DEBUG ("In FillSCTSiNtuple()"); + + m_sct_nclusters = 0; + + // loop over SCT clusters collections + for(SCT_ClusterContainer::const_iterator it=m_Sct_clcontainer->begin(); + it!=m_Sct_clcontainer->end(); it++) { + + const InDet::SCT_ClusterCollection *colNext=&(**it); + if (!colNext) continue; + + // loop over Clusters + DataVector<InDet::SCT_Cluster>::const_iterator p_clus; + for(p_clus=colNext->begin(); p_clus!=colNext->end(); ++p_clus) { + + Identifier clId = (*p_clus)->identify(); + const InDet::SCT_Cluster& cluster = **p_clus; + int GroupSize = cluster.rdoList().size(); + + if ( m_sct_nclusters < maxclusters ) { + m_sct_clx[m_sct_nclusters] = cluster.globalPosition().x(); + m_sct_cly[m_sct_nclusters] = cluster.globalPosition().y(); + m_sct_clz[m_sct_nclusters] = cluster.globalPosition().z(); + m_sct_layer[m_sct_nclusters] = m_sctID->layer_disk(clId); + m_sct_eta[m_sct_nclusters] = m_sctID->eta_module(clId); + m_sct_phi[m_sct_nclusters] = m_sctID->phi_module(clId); + m_sct_side[m_sct_nclusters] = m_sctID->side(clId); + m_sct_groupsize[m_sct_nclusters] = GroupSize; + + // Cluster Information + if (msgLvl(MSG::DEBUG)) { + msg(MSG::DEBUG) << "SCT Cluster: [" + << m_sct_layer[m_sct_nclusters] << "/" + << m_sct_eta[m_sct_nclusters] << "/" + << m_sct_phi[m_sct_nclusters] << "/" + << m_sct_side[m_sct_nclusters] << "] - (" + << m_sct_clx[m_sct_nclusters] << "," + << m_sct_cly[m_sct_nclusters] << "," + << m_sct_clz[m_sct_nclusters] << ")"; + } + + if (GroupSize>1) { + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << " strips: (" << m_sctID->strip(clId) - GroupSize + << "," << m_sctID->strip(clId) + << ") - GroupSize: " << GroupSize + << endreq; + } + else { + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << " strip: (" << m_sctID->strip(clId) + << ") - GroupSize: " << GroupSize + << endreq; + } + + } + m_sct_nclusters++; + } + // if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "number of clusters: " << m_sct_nclusters << endreq; + } + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "SCT Clusters: " << m_sct_nclusters << endreq; + return; +} + +//===================================================================== +// FillPixelSiNtuple() +//===================================================================== +void InDetAlignFillSiCluster::FillPixelSiNtuple() { + + ATH_MSG_DEBUG ("In FillPixelSiNtuple()"); + + m_pixel_nclusters = 0; + + // loop over Pixel clusters collections + for(PixelClusterContainer::const_iterator it=m_Pixel_clcontainer->begin(); + it!=m_Pixel_clcontainer->end(); it++) { + + const InDet::PixelClusterCollection *colNext=&(**it); + if (!colNext) continue; + + // loop over Clusters + DataVector<InDet::PixelCluster>::const_iterator p_clus; + for(p_clus=colNext->begin(); p_clus!=colNext->end(); ++p_clus) { + Identifier clId = (*p_clus)->identify(); + + const InDet::PixelCluster& cluster = **p_clus; + + int GroupSize = cluster.rdoList().size(); + + if (m_pixel_nclusters < maxclusters ){ + m_pixel_LocX[m_pixel_nclusters] = cluster.localPosition().x(); + m_pixel_LocY[m_pixel_nclusters] = cluster.localPosition().y(); + m_pixel_layer[m_pixel_nclusters] = m_pixelid->layer_disk(clId); + m_pixel_phi[m_pixel_nclusters] = m_pixelid->phi_module(clId); + m_pixel_groupsize[m_pixel_nclusters] = GroupSize; + m_pixel_clx[m_pixel_nclusters] = cluster.globalPosition().x(); + m_pixel_cly[m_pixel_nclusters] = cluster.globalPosition().y(); + m_pixel_clz[m_pixel_nclusters] = cluster.globalPosition().z(); + const std::vector<Identifier>& rdoList =cluster.rdoList(); + std::vector<Identifier>::const_iterator nextRDO; + int i(0); + for(nextRDO=rdoList.begin(); nextRDO !=rdoList.end(); ++nextRDO){ + Identifier rdoId = (*nextRDO); + if(i<maxclsize){ + m_pixel_clrow[m_pixel_nclusters][i] = m_pixelid->phi_index(rdoId); + m_pixel_clcol[m_pixel_nclusters][i] = m_pixelid->eta_index(rdoId); + } + ++i; + } + + // Cluster Information + if (msgLvl(MSG::DEBUG)) { + msg(MSG::DEBUG) << "Pixel Cluster (Local Pos): [" + << m_pixel_layer[m_pixel_nclusters] << "/" + << m_pixel_phi[m_pixel_nclusters] << "] - (" + << m_pixel_LocX[m_pixel_nclusters] << "," + << m_pixel_LocY[m_pixel_nclusters] << ")" + << " Groupsize: " << GroupSize << endreq; + } + + } + m_pixel_nclusters++; + } + } + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "Pixel Clusters: " << m_pixel_nclusters << endreq; + + return; +} diff --git a/InnerDetector/InDetAlignTools/InDetAlignGenTools/src/InDetAlignFillTrack.cxx b/InnerDetector/InDetAlignTools/InDetAlignGenTools/src/InDetAlignFillTrack.cxx new file mode 100755 index 0000000000000000000000000000000000000000..e7daf5ccf1b9dd1f51fc59075cf2380331729e9f --- /dev/null +++ b/InnerDetector/InDetAlignTools/InDetAlignGenTools/src/InDetAlignFillTrack.cxx @@ -0,0 +1,1290 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +////////////////////////////////////////////////////////////////////////// +// ================================================ +// InDetAlignFillTrack +// ================================================ +// +// InDetAlignFillTrack.cxx +// Source file for InDetAlignFillTrack +// +// Carlos Escobar, started 27/12/2007 +// +// AlgTool to fill track information (including truth) in a ntuple + + +#include "GaudiKernel/NTuple.h" +#include "GaudiKernel/INTupleSvc.h" +#include "GaudiKernel/SmartDataPtr.h" + +#include "GaudiKernel/IPartPropSvc.h" + +#include "DataModel/DataVector.h" + +#include "TrkTrack/Track.h" + +#include "TrkTruthData/TrackTruth.h" +#include "TrkTruthData/TrackTruthCollection.h" + +#include "TrkParameters/TrackParameters.h" + +#include "TrackRecord/TrackRecord.h" +#include "TrackRecord/TrackRecordCollection.h" + + +#include "TrkToolInterfaces/ITruthToTrack.h" + +#include "TrkTrackSummary/TrackSummary.h" +#include "TrkToolInterfaces/ITrackSummaryTool.h" +#include "TrkEventPrimitives/FitQuality.h" +#include "CLHEP/GenericFunctions/CumulativeChiSquare.hh" + +#include "InDetAlignGenTools/InDetAlignFillTrack.h" + +#include "HepMC/GenParticle.h" +#include "HepMC/GenVertex.h" + +#include <string> + +static const int maxTracks = 10000; // maximal number of tracks per event + +//===================================================================== +// InDetAlignFillTrack() +//===================================================================== +InDetAlignFillTrack::InDetAlignFillTrack (const std::string& type, + const std::string& name, + const IInterface* parent) + : AthAlgTool(type,name,parent), + m_ntupleSvc(0), + m_totaltrks(0), + m_totalhits(0), + m_totalPixhits(0), + m_totalSCThits(0), + m_totalTRThits(0), + m_totalUptrks(0), + m_totalUphits(0), + m_totalUpPixhits(0), + m_totalUpSCThits(0), + m_totalUpTRThits(0), + m_totalLowtrks(0), + m_totalLowhits(0), + m_totalLowPixhits(0), + m_totalLowSCThits(0), + m_totalLowTRThits(0), + m_events(0), + m_truthToTrack("Trk::TruthToTrack", this), + m_extrapolator("Trk::Extrapolator/CosmicsExtrapolator", this), + m_trackSumTool("Trk::TrackSummaryTool", this) +{ + declareInterface<IInDetAlignFillTrack>(this); + declareProperty("InputTrkCol", m_inputCol="Tracks"); + declareProperty("InputUpTrkCol", m_inputUpCol=""); + declareProperty("InputLowTrkCol", m_inputLowCol=""); + + // cosmic segments matching + declareProperty("doMatching", m_doMatching = true); + declareProperty("dRCut", m_matchedRcut = 100.); + declareProperty("minimumdR", m_mindR = 10000.); + + // Truth information + declareProperty("doTruth", m_doTruth = false); + declareProperty("TruthTrkCol", m_TruthTrkCol = "TrackTruthCollection"); + + // Ntuple + declareProperty("NtupleName", m_ntupleName="/NTUPLES/GLOBFILE"); + + // Tools + declareProperty("TruthToTrackTool", m_truthToTrack, + "tool to produce perigee track parameters from generated parameters"); + declareProperty("ExtrapolationTool",m_extrapolator, + "tool to extrapolate tracks"); + declareProperty("TrackSummaryTool",m_trackSumTool, + "tool to extract track info"); +} + + +//===================================================================== +// ~FillTrack() +//===================================================================== +InDetAlignFillTrack::~InDetAlignFillTrack() {} + + +//===================================================================== +// initialize() +//===================================================================== +StatusCode InDetAlignFillTrack::initialize() { + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "Initialize() of FillTrack" << endreq; + + // retrieve the NTuple Service + if (StatusCode::SUCCESS != service("NTupleSvc", m_ntupleSvc)) { + if (msgLvl(MSG::FATAL)) msg(MSG::FATAL) << "NTupleSvc service not found" << endreq; + return StatusCode::FAILURE; + } + + // get TrackSummaryTool + if ( m_trackSumTool.retrieve().isFailure() ) { + if (msgLvl(MSG::FATAL)) msg(MSG::FATAL) << "Failed to retrieve tool " << m_trackSumTool << endreq; + return StatusCode::FAILURE; + } + else { if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "Retrieved tool " << m_trackSumTool << endreq; } + + // if Truth... + if (m_doTruth) { + // Get TruthToTrack + if ( m_truthToTrack.retrieve().isFailure() ) { + if (msgLvl(MSG::FATAL)) msg(MSG::FATAL) << "Failed to retrieve tool " << m_truthToTrack << endreq; + return StatusCode::FAILURE; + } + else { if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "Retrieved tool " << m_truthToTrack << endreq; } + + // Get Extrapolator Tool + if ( m_extrapolator.retrieve().isFailure() ) { + if (msgLvl(MSG::FATAL)) msg(MSG::FATAL) << "Failed to retrieve tool " << m_extrapolator << endreq; + return StatusCode::FAILURE; + } + else { if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "Retrieved tool " << m_extrapolator << endreq; } + + // retrieve the PartPropSvc service (need for Cosmics) + IPartPropSvc* p_PartPropSvc; + static const bool CREATEIFNOTTHERE(true); + StatusCode PartPropStatus = svcLoc()->service( "PartPropSvc", p_PartPropSvc, CREATEIFNOTTHERE ); + + if (!PartPropStatus.isSuccess() || 0 == p_PartPropSvc) { + if (msgLvl(MSG::FATAL)) msg(MSG::FATAL) << "FillTrack: could not find PartPropSvc" << endreq; + return StatusCode::FAILURE; + } + m_mctable = const_cast<HepPDT::ParticleDataTable*>(p_PartPropSvc->PDT()); + if(m_mctable == 0){ + if (msgLvl(MSG::FATAL)) msg(MSG::FATAL) << "FillTrack: PDG table not found" << endreq; + return StatusCode::FAILURE; + } + } + + // Book Ntuple + bookNtuple(); + if (m_inputUpCol!="") bookUpNtuple(); + if (m_inputLowCol!="") bookLowNtuple(); + if (m_doMatching && m_inputUpCol!="" && m_inputLowCol!="") bookMatchingNtuple(); + + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "Initialize() of FillTrack successful" << endreq; + return StatusCode::SUCCESS; +} + + +//===================================================================== +// finalize() +//===================================================================== +StatusCode InDetAlignFillTrack::finalize() { + + if (msgLvl(MSG::DEBUG)) { + + msg(MSG::DEBUG) << "Finalize() of FillTrack" << endreq; + + msg(MSG::DEBUG) << "________________________________________________________" << endreq; + msg(MSG::DEBUG) << endreq; + msg(MSG::DEBUG) << " InDetAlignFillTrack Summary: " << endreq; + msg(MSG::DEBUG) << " - " << m_events << " events" << endreq; + msg(MSG::DEBUG) << " - " << m_totaltrks << " tracks" << endreq; + msg(MSG::DEBUG) << " - " << m_totalhits << " hits" << endreq; + msg(MSG::DEBUG) << " - " << m_totalPixhits << " Pixel hits" << endreq; + msg(MSG::DEBUG) << " - " << m_totalSCThits << " SCT hits" << endreq; + msg(MSG::DEBUG) << " - " << m_totalTRThits << " TRT hits" << endreq; + msg(MSG::DEBUG) << "________________________________________________________" << endreq; + + if (m_inputUpCol!="") { + msg(MSG::DEBUG) << " Up Track Summary: " << endreq; + msg(MSG::DEBUG) << " - " << m_totalUptrks << " tracks" << endreq; + msg(MSG::DEBUG) << " - " << m_totalUphits << " hits" << endreq; + msg(MSG::DEBUG) << " - " << m_totalUpPixhits << " Pixel hits" << endreq; + msg(MSG::DEBUG) << " - " << m_totalUpSCThits << " SCT hits" << endreq; + msg(MSG::DEBUG) << " - " << m_totalUpTRThits << " TRT hits" << endreq; + msg(MSG::DEBUG) << "________________________________________________________" << endreq; + } + + if (m_inputLowCol!="") { + msg(MSG::DEBUG) << " Low Track Summary: " << endreq; + msg(MSG::DEBUG) << " - " << m_totalLowtrks << " tracks" << endreq; + msg(MSG::DEBUG) << " - " << m_totalLowhits << " hits" << endreq; + msg(MSG::DEBUG) << " - " << m_totalLowPixhits << " Pixel hits" << endreq; + msg(MSG::DEBUG) << " - " << m_totalLowSCThits << " SCT hits" << endreq; + msg(MSG::DEBUG) << " - " << m_totalLowTRThits << " TRT hits" << endreq; + msg(MSG::DEBUG) << "________________________________________________________" << endreq; + } + + msg(MSG::DEBUG) << endreq; + + } + + return StatusCode::SUCCESS; +} + + +//===================================================================== +// FillTrack() +//===================================================================== +StatusCode InDetAlignFillTrack::FillTrack() { + + if (msgLvl(MSG::DEBUG)) { + msg(MSG::DEBUG) << "In FillTrack()" << endreq; + msg(MSG::DEBUG) << " event " << m_events << endreq; + } + + StatusCode sc; + + const TrackCollection* tracks = new TrackCollection; + const TrackCollection* Uptracks = new TrackCollection; + const TrackCollection* Lowtracks = new TrackCollection; + + // retrieve all tracks from TDS + if (StatusCode::SUCCESS!=evtStore()->retrieve(tracks,m_inputCol)) { + if (msgLvl(MSG::ERROR)) msg(MSG::ERROR) << "Cannot find " << m_inputCol << endreq; + return StatusCode::FAILURE; + } + + // retrieve all Up tracks from TDS + if (m_inputUpCol!="") { + if (StatusCode::SUCCESS!=evtStore()->retrieve(Uptracks,m_inputUpCol)) { + if (msgLvl(MSG::ERROR)) msg(MSG::ERROR) << "Cannot find " << m_inputUpCol << endreq; + return StatusCode::FAILURE; + } + else ATH_MSG_DEBUG( "Collection with name "<< m_inputUpCol <<" found in StoreGate" ); + } + + // retrieve all Low tracks from TDS + if (m_inputLowCol!="") { + if (StatusCode::SUCCESS!=evtStore()->retrieve(Lowtracks,m_inputLowCol)) { + if (msgLvl(MSG::ERROR)) msg(MSG::ERROR) << "Cannot find " << m_inputLowCol << endreq; + return StatusCode::FAILURE; + } + else ATH_MSG_DEBUG( "Collection with name "<< m_inputLowCol <<" found in StoreGate" ); + } + + nt_ntracks = tracks->size(); + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "Retrieved Track Collection size: " << tracks->size() << endreq; + + if (m_inputUpCol!="") { + nt_nUptracks = Uptracks->size(); + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "Retrieved Up Track Collection size: " << Uptracks->size() << endreq; + } + + if (m_inputLowCol!="") { + nt_nLowtracks = Lowtracks->size(); + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "Retrieved Low Track Collection size: " << Lowtracks->size() << endreq; + } + + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "Printing input track collection" << endreq; + m_totaltrks += dumpTrackCol(tracks); + if (m_inputUpCol!="") m_totalUptrks += dumpTrackCol(Uptracks,"Up"); + if (m_inputLowCol!="") m_totalLowtrks += dumpTrackCol(Lowtracks,"Low"); + + // matching + if (m_doMatching && m_inputUpCol!="" && m_inputLowCol!="") + if(StatusCode::SUCCESS!=dumpMatching(Uptracks,Lowtracks)) { + if (msgLvl(MSG::ERROR)) msg(MSG::ERROR) << "dumpMatching failure" << endreq; + return StatusCode::FAILURE; + } + + + // if truth is available... + if (m_doTruth) { + + int nTracks = 0; + + const TrackTruthCollection *truthCol = NULL; + + // retrieve all track truth collection from TDS + if (StatusCode::SUCCESS!=evtStore()->retrieve(truthCol, m_TruthTrkCol)) { + if (msgLvl(MSG::ERROR)) msg(MSG::ERROR) << "Cannot find " << m_inputCol << endreq; + return StatusCode::FAILURE; + } + else { + if (msgLvl(MSG::DEBUG)) { + msg(MSG::DEBUG) << "Collection with name "<< m_TruthTrkCol <<" found in StoreGate" << endreq; + msg(MSG::DEBUG) << "Retrieved "<< truthCol->size() <<" truth tracks from StoreGate" << endreq; + } + + nt_nmctracks = truthCol->size(); + + TrackCollection::const_iterator trackItr = tracks->begin(); + TrackCollection::const_iterator trackItrE = tracks->end(); + + //looping over tracks + for (; trackItr != trackItrE && nTracks < maxTracks; ++trackItr) { + + const Trk::Track* track = *trackItr; + if(track == NULL){ + if (msgLvl(MSG::WARNING)) msg(MSG::WARNING) << "No associated Trk::Track object found for track " + << nTracks << endreq; + continue; + } + + if (truthCol) { + + // the key for the truth std::map is an ElementLink<TrackCollection> object + // comprises a pointer to the track and reconstructed track collection + ElementLink<TrackCollection> trackLink; + trackLink.setElement(const_cast<Trk::Track*>(track)); + trackLink.setStorableObject(*tracks); + const ElementLink<TrackCollection> trackLink2 = trackLink; + + // trying to find the std::map entry for this reconstructed track + TrackTruthCollection::const_iterator found = truthCol->find(trackLink2); + + if (found != truthCol->end()) { + + // getting the TrackTruth object - the map element + TrackTruth trkTruth = found->second; + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "got trkTruth" << endreq; + + // probability of the reco<->truth match + float trkTruthProb = trkTruth.probability(); + + HepMcParticleLink HMPL = trkTruth.particleLink(); + + if ( HMPL.isValid()) { + const HepMC::GenParticle *genParticle = HMPL.cptr(); + + double charge = 1.0; + if (genParticle->pdg_id()<0) charge=-charge; + + Amg::Vector3D productionVertex(genParticle->production_vertex()->position().x(), + genParticle->production_vertex()->position().y(), + genParticle->production_vertex()->position().z() ); + + if (msgLvl(MSG::DEBUG)) { + msg(MSG::DEBUG) << nTracks << ". Generated Particle "<< endreq; + msg(MSG::DEBUG) << " * PDG " << genParticle->pdg_id() + << ", Status " << genParticle->status() + << ", mass " << genParticle->momentum().m() << " CLHEP::MeV/c" + << endreq; + } + + float genPt = sqrt((genParticle->momentum().x())*(genParticle->momentum().x()) + + (genParticle->momentum().y())*(genParticle->momentum().y())); + + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << " * pt " << genPt/CLHEP::GeV << " CLHEP::GeV/c" + << ", p " << genParticle->momentum().e()/CLHEP::GeV << " CLHEP::GeV/c" + << ", eta " << genParticle->momentum().eta() + << ", phi "<< genParticle->momentum().phi() << " CLHEP::rad" + << endreq; + + nt_mc_trkistruth[nTracks] = 1; + nt_mc_Trk_pdg[nTracks] = genParticle->pdg_id(); + nt_mc_Trk_prob[nTracks] = trkTruthProb; + float pX = genParticle->momentum().px(); float pY = genParticle->momentum().py(); + float genParticlePt = sqrt((pX*pX)+(pY*pY)); + nt_mc_Trk_genParticlePt[nTracks] = genParticlePt; + nt_mc_Trk_genParticleEta[nTracks] = genParticle->momentum().eta(); + nt_mc_Trk_genParticlePhi[nTracks] = genParticle->momentum().phi(); + + if(genParticle->pdg_id()==0) + ATH_MSG_WARNING("Particle with PDG 0!"); + else if(!genParticle->production_vertex()) + ATH_MSG_WARNING( "No GenVertex (generator level) production vertex found!"); + else { + // currently cannot configure the TruthToTrack tool properly + + const Trk::TrackParameters* generatedTrackPerigee = 0; + + //using a tool to produce perigee track parameters from generated parameters + generatedTrackPerigee = m_truthToTrack->makePerigeeParameters(genParticle); + + if (!generatedTrackPerigee) { + if (msgLvl(MSG::DEBUG)) { + msg(MSG::DEBUG) << "Unable to extrapolate genParticle to perigee!" << endreq; + msg(MSG::DEBUG) << "Trying to extrapolate directly to exclude material effects!" << endreq; + } + + const TrackRecordCollection* recordCollection; + + sc = evtStore()->retrieve(recordCollection, "CaloEntryLayer"); + if (sc==StatusCode::FAILURE) { + if (msgLvl(MSG::ERROR)) msg(MSG::ERROR) << "Could not get track record!" << endreq; + return sc; + } + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "reading from track record, size = " + << recordCollection->size() << endreq; + + int m_nmctracks = 0; + + for (TrackRecordCollection::const_iterator record = recordCollection->begin(); + record != recordCollection->end();++record) { + + const HepPDT::ParticleData* particle = m_mctable->particle(abs((**record).GetPDGCode())); + if (!particle) continue; + + double charge=particle->charge(); + if (std::abs(charge)<0.01) continue; + + HepGeom::Point3D<double> productionVertex = (**record).GetPosition(); + if ((**record).GetPDGCode()<0) charge=-charge; + if (fabs((**record).GetPDGCode())!=13) continue; + + Amg::Vector3D direction((**record).GetMomentum().x(), + (**record).GetMomentum().y(), + (**record).GetMomentum().z()); + + double momentum = direction.mag(); + if (momentum<500) continue; + double genPar_qOverP = 1./direction.mag(); + direction *= genPar_qOverP; + if (charge<0) genPar_qOverP = -genPar_qOverP; + + + double genPar_phi = direction.phi(); + if(genPar_phi<0.) genPar_phi = genPar_phi+(4*asin(1.)); + + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "Production vertex (x,y,z): (" + << productionVertex.x() << ", " + << productionVertex.y() << ", " + << productionVertex.z() << ")" + << endreq; + + double genPar_theta = direction.theta(); + + // Create a planar surface and transform the vertex information to a TrackParameters object + Amg::Transform3D* globalSurfaceCentre = new Amg::Transform3D(); + globalSurfaceCentre->setIdentity(); + *globalSurfaceCentre *= Amg::Translation3D(productionVertex.x(), productionVertex.y(), productionVertex.z() ); + + Trk::PlaneSurface planeSurface( globalSurfaceCentre, 5., 5. ); + + const Amg::Vector3D productionVertexAsGlobalPosition( productionVertex.x(), + productionVertex.y(), + productionVertex.z() ); + + const Trk::AtaPlane* productionVertexTrackParams + = new Trk::AtaPlane( productionVertexAsGlobalPosition, + genPar_phi, genPar_theta, genPar_qOverP, + planeSurface ); + + // Create a new perigee surface + Trk::PerigeeSurface perigeeSurface; + + if (!tracks->empty()) + perigeeSurface=((**tracks->begin()).perigeeParameters()->associatedSurface()); + + const Amg::Vector3D& perigeeGlobalPosition = perigeeSurface.center(); + + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "Surface global centre (x,y,z): (" + << perigeeGlobalPosition.x() << ", " + << perigeeGlobalPosition.y() << ", " + << perigeeGlobalPosition.z() << ")" + << endreq; + + // Extrapolate the TrackParameters object to the perigee + + // ( Establish the distance between perigee and generated vertex. + // If less than tolerance don't bother with the propagation ) + const Amg::Vector3D difference = productionVertexAsGlobalPosition - perigeeGlobalPosition; + + double distance = sqrt( difference.x() * difference.x() + difference.y() * difference.y() ); + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "Distance between perigee point and generated vertex: " + << distance/CLHEP::m << " m" << endreq; + + const Trk::TrackParameters* generatedTrackPerigee = 0; + + // Extrapolate directly to exclude material effects! + if ( distance > 1.e-4 ) { + + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "Distance between perigee and generated vertex exceeds tolerance (" + << 1.e-4 << " mm)... Extrapolating!" << endreq; + + generatedTrackPerigee = m_extrapolator->extrapolateDirectly( *productionVertexTrackParams, + perigeeSurface, + Trk::anyDirection, + false, + Trk::nonInteracting ); + if (!generatedTrackPerigee) continue; + } + + else { + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "Distance between perigee and generated vertex is less than tolerance (" + << 1.e-4 << " CLHEP::mm)... " << "No propagation to perigee required" << endreq; + + // Clone the parameters from the AtaPlane object on to perigee + generatedTrackPerigee = new Trk::Perigee( 0., 0., + genPar_phi, genPar_theta, genPar_qOverP, + Trk::PerigeeSurface(productionVertexAsGlobalPosition) ); + + } + + dumpPerigee(generatedTrackPerigee, m_nmctracks); + nt_mc_Trk_vtxX[m_nmctracks] = productionVertex.x(); + nt_mc_Trk_vtxY[m_nmctracks] = productionVertex.y(); + nt_mc_Trk_vtxZ[m_nmctracks] = productionVertex.z(); + + delete productionVertexTrackParams; + delete generatedTrackPerigee; + + m_nmctracks++; + } + } + + if (generatedTrackPerigee) { + dumpPerigee(generatedTrackPerigee, nTracks); + nt_mc_Trk_vtxX[nTracks] = genParticle->production_vertex()->position().x(); + nt_mc_Trk_vtxY[nTracks] = genParticle->production_vertex()->position().y(); + nt_mc_Trk_vtxZ[nTracks] = genParticle->production_vertex()->position().z(); + + delete generatedTrackPerigee; + + } + + } + } + nt_mc_trkistruth[nTracks] = 0; + } + } + + nTracks++; + + } + } + } // end if m_doTruth + + // Store TrkTrack branch + std::string nt0id = m_ntupleName + "/TrkTrack"; + sc = m_ntupleSvc->writeRecord(nt0id); + if (sc.isFailure()) if (msgLvl(MSG::ERROR)) msg(MSG::ERROR) << "Could not write " << nt0id << "!" << endreq; + + if (m_inputUpCol!="") { + // Store TrkTrack_Up branch + std::string nt1id = m_ntupleName + "/TrkTrack_Up"; + sc = m_ntupleSvc->writeRecord(nt1id); + if (sc.isFailure()) if (msgLvl(MSG::ERROR)) msg(MSG::ERROR) << "Could not write " << nt1id << "!" << endreq; + } + + if (m_inputLowCol!="") { + // Store TrkTrack_Low branch + std::string nt2id = m_ntupleName + "/TrkTrack_Low"; + sc = m_ntupleSvc->writeRecord(nt2id); + if (sc.isFailure()) if (msgLvl(MSG::ERROR)) msg(MSG::ERROR) << "Could not write " << nt2id << "!" << endreq; + } + + if (m_doMatching && m_inputUpCol!="" && m_inputLowCol!="") { + // Store Matching Up/Low TrkTracks branch + std::string nt3id = m_ntupleName + "/TrkTrack_Matching"; + sc = m_ntupleSvc->writeRecord(nt3id); + if (sc.isFailure()) if (msgLvl(MSG::ERROR)) msg(MSG::ERROR) << "Could not write " << nt3id << "!" << endreq; + } + + + m_events++; + + // check for negative values + if (m_totaltrks<0) m_totaltrks = 0; + if (m_totalhits<0) m_totalhits = 0; + if (m_totalPixhits<0) m_totalPixhits = 0; + if (m_totalSCThits<0) m_totalSCThits = 0; + if (m_totalTRThits<0) m_totalTRThits = 0; + + if (m_inputUpCol!="") { + if (m_totalUptrks<0) m_totalUptrks = 0; + if (m_totalUphits<0) m_totalUphits = 0; + if (m_totalUpPixhits<0) m_totalUpPixhits = 0; + if (m_totalUpSCThits<0) m_totalUpSCThits = 0; + if (m_totalUpTRThits<0) m_totalUpTRThits = 0; + } + + if (m_inputLowCol!="") { + if (m_totalLowtrks<0) m_totalLowtrks = 0; + if (m_totalLowhits<0) m_totalLowhits = 0; + if (m_totalLowPixhits<0) m_totalLowPixhits = 0; + if (m_totalLowSCThits<0) m_totalLowSCThits = 0; + if (m_totalLowTRThits<0) m_totalLowTRThits = 0; + } + + if (m_events<0) m_events = 0; + + return StatusCode::SUCCESS; +} + + +//===================================================================== +// bookNtuple() +//===================================================================== +void InDetAlignFillTrack::bookNtuple() { + + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "Booking Trk::Track Info..." << endreq; + + NTupleFilePtr file1(m_ntupleSvc, m_ntupleName); + std::string nt0id = m_ntupleName + "/TrkTrack"; + std::string comments = "Trk::Track Information"; + + NTuplePtr nt0(m_ntupleSvc, nt0id); + if (nt0) {if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "Ntuple is already booked" << endreq;} + else { + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "Attempting to book general ntuple" << endreq; + nt0 = m_ntupleSvc->book(nt0id,CLID_ColumnWiseTuple,comments); + + if (nt0) { + StatusCode sc; + // nt0->addItem("event", nt_event); // event number + sc = nt0->addItem("nTracks", nt_ntracks, 0, maxTracks); // number of tracks + if (m_doTruth) sc =nt0->addItem("mc_nTracks", nt_nmctracks, 0, maxTracks); // number of mc tracks + + // ---------------------------------------------------------------------- + // Trk::Track parameters + sc = nt0->addItem("Trk_d0", nt_ntracks, nt_Trk_d0); + sc = nt0->addItem("Trk_z0", nt_ntracks, nt_Trk_z0); + sc = nt0->addItem("Trk_phi0", nt_ntracks, nt_Trk_phi0); + sc = nt0->addItem("Trk_theta0", nt_ntracks, nt_Trk_theta0); + sc = nt0->addItem("Trk_qoverp", nt_ntracks, nt_Trk_qoverp); + sc = nt0->addItem("Trk_pt", nt_ntracks, nt_Trk_pt); + // ---------------------------------------------------------------------- + + // ---------------------------------------------------------------------- + // Trk::Track hits... + sc = nt0->addItem("Trk_nHits", nt_ntracks, nt_Trk_nHits); + sc = nt0->addItem("Trk_nhitsPixels", nt_ntracks, nt_Trk_nhitspix); + sc = nt0->addItem("Trk_nhitsSCT", nt_ntracks, nt_Trk_nhitssct); + sc = nt0->addItem("Trk_nhitsTRT", nt_ntracks, nt_Trk_nhitstrt); + + sc = nt0->addItem("Trk_nsharedPixels", nt_ntracks, nt_Trk_nsharedPixels); + sc = nt0->addItem("Trk_nsharedSCT", nt_ntracks, nt_Trk_nsharedSCT); + sc = nt0->addItem("Trk_nshared", nt_ntracks, nt_Trk_nshared); + + sc = nt0->addItem("Trk_nholesPixels", nt_ntracks, nt_Trk_nholesPixels); + sc = nt0->addItem("Trk_nholesSCT", nt_ntracks, nt_Trk_nholesSCT); + sc = nt0->addItem("Trk_nholes", nt_ntracks, nt_Trk_nholes); + + sc = nt0->addItem("Trk_chi2", nt_ntracks, nt_Trk_chi2); + sc = nt0->addItem("Trk_ndof", nt_ntracks, nt_Trk_ndof); + sc = nt0->addItem("Trk_chi2Prob", nt_ntracks, nt_Trk_chi2Prob); + // ---------------------------------------------------------------------- + + if (m_doTruth) { + // ---------------------------------------------------------------------- + sc = nt0->addItem("mc_TrkIsTruth", nt_ntracks, nt_mc_trkistruth); + + // generated particle parameters + sc = nt0->addItem("mc_Trk_genParticlePt", nt_ntracks, nt_mc_Trk_genParticlePt); + sc = nt0->addItem("mc_Trk_genParticleEta", nt_ntracks, nt_mc_Trk_genParticleEta); + sc = nt0->addItem("mc_Trk_genParticlePhi", nt_ntracks, nt_mc_Trk_genParticlePhi); + + // MonteCarlo Track parameters + sc = nt0->addItem("mc_Trk_d0", nt_ntracks, nt_mc_Trk_d0); + sc = nt0->addItem("mc_Trk_z0", nt_ntracks, nt_mc_Trk_z0); + sc = nt0->addItem("mc_Trk_phi0", nt_ntracks, nt_mc_Trk_phi0); + sc = nt0->addItem("mc_Trk_theta", nt_ntracks, nt_mc_Trk_theta0); + sc = nt0->addItem("mc_Trk_eta", nt_ntracks, nt_mc_Trk_eta); + sc = nt0->addItem("mc_Trk_qoverp", nt_ntracks, nt_mc_Trk_qoverp); + sc = nt0->addItem("mc_Trk_qoverpt", nt_ntracks, nt_mc_Trk_qoverpt); + sc = nt0->addItem("mc_Trk_pt", nt_ntracks, nt_mc_Trk_pt); + sc = nt0->addItem("mc_Trk_charge", nt_ntracks, nt_mc_Trk_charge); + sc = nt0->addItem("mc_Trk_prob", nt_ntracks, nt_mc_Trk_prob); + sc = nt0->addItem("mc_Trk_pdg", nt_ntracks, nt_mc_Trk_pdg); + + sc = nt0->addItem("mc_Trk_vtxX", nt_ntracks, nt_mc_Trk_vtxX); + sc = nt0->addItem("mc_Trk_vtxY", nt_ntracks, nt_mc_Trk_vtxY); + sc = nt0->addItem("mc_Trk_vtxZ", nt_ntracks, nt_mc_Trk_vtxZ); + // ---------------------------------------------------------------------- + } + + if (sc.isFailure()) msg(MSG::FATAL) << "Failed ntupleSvc()" << endreq; + else msg(MSG::DEBUG) << "Ntuple " << nt0id << " has been booked successfully! " << endreq; + //if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "Ntuple " << nt0id << " has been booked successfully! " << endreq; + } + + else { if (msgLvl(MSG::ERROR)) msg(MSG::ERROR) << "Error booking ntuple" << endreq; } + } + + // return; +} + + +//===================================================================== +// bookUpNtuple() +//===================================================================== +void InDetAlignFillTrack::bookUpNtuple() { + + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "Booking Up Trk::Track Info..." << endreq; + + std::string nt1id = m_ntupleName + "/TrkTrack_Up"; + std::string comments = "Trk::UpTrack Information"; + + NTuplePtr nt1(m_ntupleSvc, nt1id); + if (nt1) {if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "Ntuple is already booked" << endreq;} + else { + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "Attempting to book general ntuple" << endreq; + nt1 = m_ntupleSvc->book(nt1id,CLID_ColumnWiseTuple,comments); + + if (nt1) { + StatusCode sc; + // nt1->addItem("event", nt_event); // event number + sc = nt1->addItem("nTracks_Up", nt_nUptracks, 0, maxTracks); // number of tracks + + // ---------------------------------------------------------------------- + // Trk::Track parameters + sc = nt1->addItem("Trk_d0_Up", nt_nUptracks, nt_Trk_d0_Up); + sc = nt1->addItem("Trk_z0_Up", nt_nUptracks, nt_Trk_z0_Up); + sc = nt1->addItem("Trk_phi0_Up", nt_nUptracks, nt_Trk_phi0_Up); + sc = nt1->addItem("Trk_theta0_Up", nt_nUptracks, nt_Trk_theta0_Up); + sc = nt1->addItem("Trk_qoverp_Up", nt_nUptracks, nt_Trk_qoverp_Up); + sc = nt1->addItem("Trk_pt_Up", nt_nUptracks, nt_Trk_pt_Up); + // ---------------------------------------------------------------------- + + // ---------------------------------------------------------------------- + // Trk::Track hits... + sc = nt1->addItem("Trk_nHits_Up", nt_nUptracks, nt_Trk_nHits_Up); + sc = nt1->addItem("Trk_nhitsPixels_Up", nt_nUptracks, nt_Trk_nhitspix_Up); + sc = nt1->addItem("Trk_nhitsSCT_Up", nt_nUptracks, nt_Trk_nhitssct_Up); + sc = nt1->addItem("Trk_nhitsTRT_Up", nt_nUptracks, nt_Trk_nhitstrt_Up); + + sc = nt1->addItem("Trk_nsharedPixels_Up", nt_nUptracks, nt_Trk_nsharedPixels_Up); + sc = nt1->addItem("Trk_nsharedSCT_Up", nt_nUptracks, nt_Trk_nsharedSCT_Up); + sc = nt1->addItem("Trk_nshared_Up", nt_nUptracks, nt_Trk_nshared_Up); + + sc = nt1->addItem("Trk_nholesPixels_Up", nt_nUptracks, nt_Trk_nholesPixels_Up); + sc = nt1->addItem("Trk_nholesSCT_Up", nt_nUptracks, nt_Trk_nholesSCT_Up); + sc = nt1->addItem("Trk_nholes_Up", nt_nUptracks, nt_Trk_nholes_Up); + + sc = nt1->addItem("Trk_chi2_Up", nt_nUptracks, nt_Trk_chi2_Up); + sc = nt1->addItem("Trk_ndof_Up", nt_nUptracks, nt_Trk_ndof_Up); + sc = nt1->addItem("Trk_chi2Prob_Up", nt_nUptracks, nt_Trk_chi2Prob_Up); + // ---------------------------------------------------------------------- + + if (sc.isFailure()) msg(MSG::FATAL) << "Failed ntupleSvc()" << endreq; + else msg(MSG::DEBUG) << "Ntuple " << nt1id << " has been booked successfully! " << endreq; + } + else { if (msgLvl(MSG::ERROR)) msg(MSG::ERROR) << "Error booking ntuple" << endreq; } + } + +} + + +//===================================================================== +// bookLowNtuple() +//===================================================================== +void InDetAlignFillTrack::bookLowNtuple() { + + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "Booking Low Trk::Track Info..." << endreq; + + std::string nt2id = m_ntupleName + "/TrkTrack_Low"; + std::string comments = "Trk::LowTrack Information"; + + NTuplePtr nt2(m_ntupleSvc, nt2id); + if (nt2) {if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "Ntuple is already booked" << endreq;} + else { + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "Attempting to book general ntuple" << endreq; + nt2 = m_ntupleSvc->book(nt2id,CLID_ColumnWiseTuple,comments); + + if (nt2) { + StatusCode sc; + // sc = nt2->addItem("event", nt_event); // event number + sc = nt2->addItem("nTracks_Low", nt_nLowtracks, 0, maxTracks); // number of tracks + + // ---------------------------------------------------------------------- + // Trk::Track parameters + sc = nt2->addItem("Trk_d0_Low", nt_nLowtracks, nt_Trk_d0_Low); + sc = nt2->addItem("Trk_z0_Low", nt_nLowtracks, nt_Trk_z0_Low); + sc = nt2->addItem("Trk_phi0_Low", nt_nLowtracks, nt_Trk_phi0_Low); + sc = nt2->addItem("Trk_theta0_Low", nt_nLowtracks, nt_Trk_theta0_Low); + sc = nt2->addItem("Trk_qoverp_Low", nt_nLowtracks, nt_Trk_qoverp_Low); + sc = nt2->addItem("Trk_pt_Low", nt_nLowtracks, nt_Trk_pt_Low); + // ---------------------------------------------------------------------- + + // ---------------------------------------------------------------------- + // Trk::Track hits... + sc = nt2->addItem("Trk_nHits_Low", nt_nLowtracks, nt_Trk_nHits_Low); + sc = nt2->addItem("Trk_nhitsPixels_Low", nt_nLowtracks, nt_Trk_nhitspix_Low); + sc = nt2->addItem("Trk_nhitsSCT_Low", nt_nLowtracks, nt_Trk_nhitssct_Low); + sc = nt2->addItem("Trk_nhitsTRT_Low", nt_nLowtracks, nt_Trk_nhitstrt_Low); + + sc = nt2->addItem("Trk_nsharedPixels_Low", nt_nLowtracks, nt_Trk_nsharedPixels_Low); + sc = nt2->addItem("Trk_nsharedSCT_Low", nt_nLowtracks, nt_Trk_nsharedSCT_Low); + sc = nt2->addItem("Trk_nshared_Low", nt_nLowtracks, nt_Trk_nshared_Low); + + sc = nt2->addItem("Trk_nholesPixels_Low", nt_nLowtracks, nt_Trk_nholesPixels_Low); + sc = nt2->addItem("Trk_nholesSCT_Low", nt_nLowtracks, nt_Trk_nholesSCT_Low); + sc = nt2->addItem("Trk_nholes_Low", nt_nLowtracks, nt_Trk_nholes_Low); + + sc = nt2->addItem("Trk_chi2_Low", nt_nLowtracks, nt_Trk_chi2_Low); + sc = nt2->addItem("Trk_ndof_Low", nt_nLowtracks, nt_Trk_ndof_Low); + sc = nt2->addItem("Trk_chi2Prob_Low", nt_nLowtracks, nt_Trk_chi2Prob_Low); + // ---------------------------------------------------------------------- + if (sc.isFailure()) msg(MSG::FATAL) << "Failed ntupleSvc()" << endreq; + else msg(MSG::DEBUG) << "Ntuple " << nt2id << " has been booked successfully! " << endreq; + + } + else { if (msgLvl(MSG::ERROR)) msg(MSG::ERROR) << "Error booking ntuple" << endreq; } + } +} + + +//===================================================================== +// bookMatchingNtuple() +//===================================================================== +void InDetAlignFillTrack::bookMatchingNtuple() { + + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "Booking Matching between Up and Low Trk::Track Collections..." << endreq; + + std::string nt3id = m_ntupleName + "/TrkTrack_Matching"; + std::string comments = "Matching between Up and Low Trk::Track Collections"; + + NTuplePtr nt3(m_ntupleSvc, nt3id); + if (nt3) {if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "Ntuple is already booked" << endreq;} + else { + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "Attempting to book general ntuple" << endreq; + m_ntupleSvc->book(nt3id,CLID_ColumnWiseTuple,comments); + + if (nt3) { + StatusCode sc; + // sc = nt3->addItem("event", nt_event); // event number + + sc = nt3->addItem("nTracks_Match",nt_matchingTrk, 0, maxTracks); // number of tracks + + // ---------------------------------------------------------------------- + // Matching for the usual Trk::Track parameters + sc = nt3->addItem("Trk_delta_d0", nt_matchingTrk, nt_Trk_delta_d0); + sc = nt3->addItem("Trk_delta_phi0", nt_matchingTrk, nt_Trk_delta_phi0); + sc = nt3->addItem("Trk_delta_theta0", nt_matchingTrk, nt_Trk_delta_theta0); + sc = nt3->addItem("Trk_delta_eta", nt_matchingTrk, nt_Trk_delta_eta); + sc = nt3->addItem("Trk_delta_z0", nt_matchingTrk, nt_Trk_delta_z0); + sc = nt3->addItem("Trk_delta_qoverpt", nt_matchingTrk, nt_Trk_delta_qoverpt); + sc = nt3->addItem("Trk_delta_pt", nt_matchingTrk, nt_Trk_delta_pt); + sc = nt3->addItem("Trk_delta_charge", nt_matchingTrk, nt_Trk_delta_charge); + // ---------------------------------------------------------------------- + + if (sc.isFailure()) msg(MSG::FATAL) << "Failed ntupleSvc()" << endreq; + else msg(MSG::DEBUG) << "Ntuple " << nt3id << " has been booked successfully! " << endreq; + + } + else { if (msgLvl(MSG::ERROR)) msg(MSG::DEBUG) << "Error booking ntuple" << endreq; } + } + + return; +} + + +//===================================================================== +// InDetAlignFillTrack::dumpTrackCol() +//===================================================================== +int InDetAlignFillTrack::dumpTrackCol(const TrackCollection* tracks) { + return dumpTrackCol(tracks,""); +} + +//===================================================================== +// InDetAlignFillTrack::dumpTrackCol() +//===================================================================== +int InDetAlignFillTrack::dumpTrackCol(const TrackCollection* tracks, + const std::string TrkColName) { + + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "In dump" << TrkColName << "TrackCol()" << endreq; + + int itrk = 0; + + TrackCollection::const_iterator trackItr = tracks->begin(); + TrackCollection::const_iterator trackItrE = tracks->end(); + + //looping over tracks + for (; trackItr != trackItrE && itrk < maxTracks; ++trackItr) { + + if(*trackItr!=0) + dumpTrack(itrk, (*trackItr), TrkColName); + + itrk++; + } + + return itrk; +} + + +//===================================================================== +// InDetAlignFillTrack::dumpTrack() +//===================================================================== +void InDetAlignFillTrack::dumpTrack(int itrk, const Trk::Track* trk, + std::string TrkColName) { + if (msgLvl(MSG::VERBOSE)) msg(MSG::VERBOSE) << "In dump" << TrkColName << "Track()" << endreq; + + const Trk::Perigee *aMeasPer =trk->perigeeParameters(); + + if (aMeasPer==0){ + if (msgLvl(MSG::ERROR)) msg(MSG::ERROR) << "Could not get Trk::MeasuredPerigee" << endreq;} + else { + double d0 = aMeasPer->parameters()[Trk::d0]; + double z0 = aMeasPer->parameters()[Trk::z0]; + double phi0 = aMeasPer->parameters()[Trk::phi0]; + double theta = aMeasPer->parameters()[Trk::theta]; + double qOverP = aMeasPer->parameters()[Trk::qOverP]; + + if (msgLvl(MSG::DEBUG)) { + msg(MSG::DEBUG) << itrk << ". " << TrkColName + << " Track Parameters (d0, z0, phi0, theta, q/p)" << endreq; + msg(MSG::DEBUG) << " " << d0 << ", " << z0 << ", " + << phi0 << ", " << theta << ", " << qOverP << endreq; + } + + float transverseMomentum = sqrt((aMeasPer->momentum().x())*(aMeasPer->momentum().x()) + + (aMeasPer->momentum().y())*(aMeasPer->momentum().y())); + + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << " p = " << aMeasPer->momentum().mag()/CLHEP::GeV << " CLHEP::GeV/c, " + << " pT = " << transverseMomentum/CLHEP::GeV << " CLHEP::GeV/c" << endreq; + + // number of hits + int nHits = (*trk).measurementsOnTrack()->size(); + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << " - number of hits : " << nHits << endreq; + if (nHits!=0) { + if (TrkColName=="Up") m_totalUphits += nHits; + else if (TrkColName=="Low") m_totalLowhits += nHits; + else m_totalhits += nHits; + } + + // number of hits, shared hits and holes + int nhitspix=0, nhitssct=0, nhitstrt=0; + int nshared=0, nshpix=0, nshsct=0; + int nholes=0, nhpix=0, nhsct=0; + + const Trk::TrackSummary* summary = m_trackSumTool->createSummary((*trk)); + if(summary==0) + ATH_MSG_ERROR("Could not create TrackSummary"); + + else { + + // hits + nhitspix = summary->get(Trk::numberOfPixelHits); + nhitssct = summary->get(Trk::numberOfSCTHits); + nhitstrt = summary->get(Trk::numberOfTRTHits); + + if (msgLvl(MSG::DEBUG)) { + msg(MSG::DEBUG) << " -- number of Pixel hits : " << nhitspix << endreq; + msg(MSG::DEBUG) << " -- number of SCT hits : " << nhitssct << endreq; + msg(MSG::DEBUG) << " -- number of TRT hits : " << nhitstrt << endreq; + } + + if (nhitspix!=0) { + if (TrkColName=="Up") m_totalUpPixhits += nhitspix; + else if (TrkColName=="Low") m_totalLowPixhits += nhitspix; + else m_totalPixhits += nhitspix; + } + if (nhitssct!=0) { + if (TrkColName=="Up") m_totalUpSCThits += nhitssct; + else if (TrkColName=="Low") m_totalLowSCThits += nhitssct; + else m_totalSCThits += nhitssct; + } + if (nhitstrt!=0) { + if (TrkColName=="Up") m_totalUpTRThits += nhitstrt; + else if (TrkColName=="Low") m_totalLowTRThits += nhitstrt; + else m_totalTRThits += nhitstrt; + } + + // shared hits + nshpix = summary->get(Trk::numberOfPixelSharedHits); + nshsct = summary->get(Trk::numberOfSCTSharedHits); + nshared = nshpix + nshsct; + + if (nshpix<0) nshpix=0; + if (nshsct<0) nshsct=0; + if (nshared<0) nshared=0; + + if (msgLvl(MSG::DEBUG)) { + msg(MSG::DEBUG) << " - number of shared hits : " << nshared << endreq; + msg(MSG::DEBUG) << " -- number of Pixel shared hits : " << nshpix << endreq; + msg(MSG::DEBUG) << " -- number of SCT shared hits : " << nshsct << endreq; + } + + // holes + nhpix = summary->get(Trk::numberOfPixelHoles); + nhsct = summary->get(Trk::numberOfSCTHoles); + nholes = nhpix + nhsct; + + if (nhpix<0) nhpix=0; + if (nhsct<0) nhsct=0; + if (nholes<0) nholes=0; + + if (msgLvl(MSG::DEBUG)) { + msg(MSG::DEBUG) << " - number of Pixel holes : " << nhpix << endreq; + msg(MSG::DEBUG) << " -- number of SCT holes : " << nhsct << endreq; + msg(MSG::DEBUG) << " -- number of holes : " << nholes << endreq; + } + + } + + delete summary; + + // get fit quality and chi2 probability of track + // chi2Prob = TMath::Prob(chi2,DoF) ROOT function + double chi2Prob=0.; + const Trk::FitQuality* fitQual = (*trk).fitQuality(); + if (fitQual==0) { + if (msgLvl(MSG::ERROR)) msg(MSG::ERROR) << "No fit quality assigned to the track" << endreq; + chi2Prob = -1e12; // return big value + } + else { + if (fitQual->chiSquared() > 0. && fitQual->numberDoF() > 0) { + Genfun::CumulativeChiSquare probabilityFunction( fitQual->numberDoF() ); + chi2Prob = 1 - probabilityFunction( fitQual->chiSquared() ); + + if (msgLvl(MSG::DEBUG)) { + msg(MSG::DEBUG) << " - chi2 : " << fitQual->chiSquared() << endreq; + msg(MSG::DEBUG) << " - ndof : " << fitQual->numberDoF() << endreq; + msg(MSG::DEBUG) << " - chi2 propability : " << chi2Prob << endreq; + } + + } + } + // Fill ntuple + if (TrkColName=="Up") { + nt_Trk_d0_Up[itrk] = d0; + nt_Trk_z0_Up[itrk] = z0; + nt_Trk_phi0_Up[itrk] = phi0; + nt_Trk_theta0_Up[itrk] = theta; + nt_Trk_qoverp_Up[itrk] = qOverP; + nt_Trk_pt_Up[itrk] = transverseMomentum; + + nt_Trk_nHits_Up[itrk] = nHits; + nt_Trk_nhitspix_Up[itrk] = nhitspix; + nt_Trk_nhitssct_Up[itrk] = nhitssct; + nt_Trk_nhitstrt_Up[itrk] = nhitstrt; + + nt_Trk_nsharedPixels_Up[itrk] = nshpix; + nt_Trk_nsharedSCT_Up[itrk] = nshsct; + nt_Trk_nshared_Up[itrk] = nshared; + + nt_Trk_nholesPixels_Up[itrk] = nhpix; + nt_Trk_nholesSCT_Up[itrk] = nhsct; + nt_Trk_nholes_Up[itrk] = nholes; + + nt_Trk_chi2_Up[itrk] = fitQual->chiSquared(); + nt_Trk_ndof_Up[itrk] = fitQual->numberDoF(); + nt_Trk_chi2Prob_Up[itrk] = chi2Prob; + } + else if (TrkColName=="Low") { + nt_Trk_d0_Low[itrk] = d0; + nt_Trk_z0_Low[itrk] = z0; + nt_Trk_phi0_Low[itrk] = phi0; + nt_Trk_theta0_Low[itrk] = theta; + nt_Trk_qoverp_Low[itrk] = qOverP; + nt_Trk_pt_Low[itrk] = transverseMomentum; + + nt_Trk_nHits_Low[itrk] = nHits; + nt_Trk_nhitspix_Low[itrk] = nhitspix; + nt_Trk_nhitssct_Low[itrk] = nhitssct; + nt_Trk_nhitstrt_Low[itrk] = nhitstrt; + + nt_Trk_nsharedPixels_Low[itrk] = nshpix; + nt_Trk_nsharedSCT_Low[itrk] = nshsct; + nt_Trk_nshared_Low[itrk] = nshared; + + nt_Trk_nholesPixels_Low[itrk] = nhpix; + nt_Trk_nholesSCT_Low[itrk] = nhsct; + nt_Trk_nholes_Low[itrk] = nholes; + + nt_Trk_chi2_Low[itrk] = fitQual->chiSquared(); + nt_Trk_ndof_Low[itrk] = fitQual->numberDoF(); + nt_Trk_chi2Prob_Low[itrk] = chi2Prob; + } + else { + nt_Trk_d0[itrk] = d0; + nt_Trk_z0[itrk] = z0; + nt_Trk_phi0[itrk] = phi0; + nt_Trk_theta0[itrk] = theta; + nt_Trk_qoverp[itrk] = qOverP; + nt_Trk_pt[itrk] = transverseMomentum; + + nt_Trk_nHits[itrk] = nHits; + nt_Trk_nhitspix[itrk] = nhitspix; + nt_Trk_nhitssct[itrk] = nhitssct; + nt_Trk_nhitstrt[itrk] = nhitstrt; + + nt_Trk_nsharedPixels[itrk] = nshpix; + nt_Trk_nsharedSCT[itrk] = nshsct; + nt_Trk_nshared[itrk] = nshared; + + nt_Trk_nholesPixels[itrk] = nhpix; + nt_Trk_nholesSCT[itrk] = nhsct; + nt_Trk_nholes[itrk] = nholes; + + nt_Trk_chi2[itrk] = fitQual->chiSquared(); + nt_Trk_ndof[itrk] = fitQual->numberDoF(); + nt_Trk_chi2Prob[itrk] = chi2Prob; + } + } + + return; +} + + +//===================================================================== +// InDetAlignFillTrack::dumpPerigee() +//===================================================================== +void InDetAlignFillTrack::dumpPerigee(const Trk::TrackParameters* generatedTrackPerigee, + int index) { + + if (msgLvl(MSG::VERBOSE)) msg(MSG::VERBOSE) << "In dumpPerigee()" << endreq; + + float d0 = generatedTrackPerigee->parameters()[Trk::d0]; + float z0 = generatedTrackPerigee->parameters()[Trk::z0]; + float phi0 = generatedTrackPerigee->parameters()[Trk::phi0]; + float theta = generatedTrackPerigee->parameters()[Trk::theta]; + float eta = generatedTrackPerigee->eta(); + float charge = generatedTrackPerigee->charge(); + float qoverp = generatedTrackPerigee->parameters()[Trk::qOverP]; + float qoverpt = generatedTrackPerigee->parameters()[Trk::qOverP]/(sin(theta)); + float pt = (1/qoverpt)*(charge); + // int pdg = genParticle->pdg_id(); + + if (msgLvl(MSG::DEBUG)) { + msg(MSG::DEBUG) << " - Extrapolated genParticle perigee parameters " + << "(d0, z0, phi0, theta, q/p)" << endreq; + msg(MSG::DEBUG) << " " << d0 << ", " << z0 << ", " + << phi0 << ", " << theta << ", " << qoverp << endreq; + + msg(MSG::DEBUG) << " p = " << fabs(1/qoverp)/CLHEP::GeV << " CLHEP::GeV/c, " + << " pT = " << pt/CLHEP::GeV << " CLHEP::GeV/c" << endreq; + } + + nt_mc_Trk_d0[index] = d0; + nt_mc_Trk_z0[index] = z0; + nt_mc_Trk_phi0[index] = phi0; + nt_mc_Trk_theta0[index] = theta; + nt_mc_Trk_eta[index] = eta; + nt_mc_Trk_qoverp[index] = qoverp; + nt_mc_Trk_qoverpt[index] = qoverpt; + nt_mc_Trk_pt[index] = pt; + nt_mc_Trk_charge[index] = charge; + + return; +} + + +//===================================================================== +// InDetAlignFillTrack::dumpMatching() +//===================================================================== +StatusCode InDetAlignFillTrack::dumpMatching(const TrackCollection* tracksUpper, + const TrackCollection* tracksLower) { + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "In dumpMatching()" << endreq; + + int nTracksUpper = 0; + + // looping over the upper barrel tracks + TrackCollection::const_iterator trackItrUpper = tracksUpper->begin(); + TrackCollection::const_iterator trackItrUpperE = tracksUpper->end(); + for (; trackItrUpper != trackItrUpperE; ++trackItrUpper) { + + const Trk::Track* trackUpper = *trackItrUpper; + if(trackUpper == NULL) { + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "No associated Trk::Track object found for track "<< nTracksUpper << endreq; + continue; + } + + const Trk::Perigee* measUpperPer = trackUpper->perigeeParameters(); + + // Get the track parameters from the Upper Track + float d0Up = measUpperPer->parameters()[Trk::d0]; + float phi0Up = measUpperPer->parameters()[Trk::phi0]; + float eta0Up = measUpperPer->eta(); + float z0Up = measUpperPer->parameters()[Trk::z0]; + float thetaUp = measUpperPer->parameters()[Trk::theta]; + float qOverPtUp = measUpperPer->parameters()[Trk::qOverP]*1000/sin(thetaUp); + float chargeUp = measUpperPer->charge(); + float ptUp = measUpperPer->pT()/1000.; + + if (msgLvl(MSG::DEBUG)) { + msg(MSG::DEBUG) << nTracksUpper << ". UpTrack -> Track Parameters:" << endreq; + msg(MSG::DEBUG) << " d0, z0, phi0, theta, q/p" << d0Up << ", " << z0Up << ", " + << phi0Up << ", " << thetaUp << ", " << qOverPtUp << endreq; + } + + int nTracksLower = 0; + + bool matchFound = false; + float Matched_Low_d0 = -999; + float Matched_Low_phi0 = -999; + //** + float Matched_Low_theta = -999; + float Matched_Low_eta0 = -999; + float Matched_Low_z0 = -999; + float Matched_Low_qOverPt = -999; + float Matched_Low_charge = -999; + float Matched_Low_pt = -999; + + // looping over the lower barrel tracks + DataVector<Trk::Track>::const_iterator trackItrLower = tracksLower->begin(); + DataVector<Trk::Track>::const_iterator trackItrLowerE = tracksLower->end(); + for (; trackItrLower != trackItrLowerE; ++trackItrLower) { //looping over Lower tracks + + const Trk::Track* trackLower = *trackItrLower; + if(trackLower == NULL){ + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "No associated Trk::Track object found for track "<< nTracksLower << endreq; + continue; + } + + const Trk::Perigee* measLowerPer = trackLower->perigeeParameters(); + + float d0Low = measLowerPer->parameters()[Trk::d0]; + float phi0Low = measLowerPer->parameters()[Trk::phi0]; + float eta0Low = measLowerPer->eta(); + float z0Low = measLowerPer->parameters()[Trk::z0]; + float thetaLow = measLowerPer->parameters()[Trk::theta]; + float qOverPtLow = measLowerPer->parameters()[Trk::qOverP]*1000/sin(thetaLow); + float chargeLow = measLowerPer->charge(); + float ptLow = measLowerPer->pT()/1000.; + + //selecting Lower track that is closest to Upper in eta-phi + float dphi2 = (phi0Up - phi0Low)*(phi0Up - phi0Low); + + //For TRT only tracks we will ignore the delta eta + // and just require a delta phi match + float deta2 = (eta0Up - eta0Low)*(eta0Up - eta0Low); + + // look for a matching track within a cone + float dR = sqrt(dphi2 + deta2); + ATH_MSG_DEBUG("dR (sqrt(dphi+deta2): " << dR); + ATH_MSG_DEBUG("minmdR (sqrt(dphi+deta2): " << m_mindR); + if(dR < m_mindR) { + // m_mindR = dR; + Matched_Low_d0 = d0Low; + Matched_Low_phi0 = phi0Low; + //** + Matched_Low_theta = thetaLow; + Matched_Low_eta0 = eta0Low; + Matched_Low_z0 = z0Low; + Matched_Low_qOverPt = qOverPtLow; + Matched_Low_charge = chargeLow; + Matched_Low_pt = ptLow; + + if(dR < m_matchedRcut) { + + if (msgLvl(MSG::DEBUG)) { + msg(MSG::DEBUG) << nTracksLower << ". LowTrack -> Track Parameters:" << endreq; + msg(MSG::DEBUG) << " d0, z0, phi0, theta, q/p: " << d0Low << ", " << z0Low << ", " + << phi0Low << ", " << thetaLow << ", " << qOverPtLow << endreq; + } + + matchFound = true; + } + } + else { + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "No matching low track found!!" << endreq; + } + nTracksLower++; + + } //looping over lower tracks + + if(matchFound) { + + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "Match found!" <<endreq; + nt_matchingTrk = 1; + nt_Trk_delta_d0[nTracksUpper] = d0Up - Matched_Low_d0; + nt_Trk_delta_phi0[nTracksUpper] = phi0Up - Matched_Low_phi0; + //** + nt_Trk_delta_theta0[nTracksUpper] = thetaUp - Matched_Low_theta; + nt_Trk_delta_eta[nTracksUpper] = eta0Up - Matched_Low_eta0; + nt_Trk_delta_z0[nTracksUpper] = z0Up - Matched_Low_z0; + nt_Trk_delta_qoverpt[nTracksUpper] = qOverPtUp - Matched_Low_qOverPt; + nt_Trk_delta_pt[nTracksUpper] = ptUp - Matched_Low_pt; + nt_Trk_delta_charge[nTracksUpper] = chargeUp - Matched_Low_charge; + + } // end match found + + + nTracksUpper++; + + } //looping over upper tracks + + return StatusCode::SUCCESS; +} diff --git a/InnerDetector/InDetAlignTools/InDetAlignGenTools/src/InDetAlignHitQualSelTool.cxx b/InnerDetector/InDetAlignTools/InDetAlignGenTools/src/InDetAlignHitQualSelTool.cxx new file mode 100644 index 0000000000000000000000000000000000000000..e090fd9696d7996dc75e15b073d3b0f5af224ea0 --- /dev/null +++ b/InnerDetector/InDetAlignTools/InDetAlignGenTools/src/InDetAlignHitQualSelTool.cxx @@ -0,0 +1,244 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +#include "TrkTrack/TrackStateOnSurface.h" +#include "TrkMeasurementBase/MeasurementBase.h" +#include "TrkRIO_OnTrack/RIO_OnTrack.h" +#include "TrkPrepRawData/PrepRawData.h" +#include "InDetPrepRawData/SiCluster.h" +// geometry +#include "Identifier/Identifier.h" +#include "InDetReadoutGeometry/SiDetectorElement.h" +#include "AtlasDetDescr/AtlasDetectorID.h" +#include "InDetReadoutGeometry/PixelDetectorManager.h" +#include "InDetReadoutGeometry/SCT_DetectorManager.h" + +#include "InDetAlignGenTools/InDetAlignHitQualSelTool.h" + +using namespace std ; + +InDetAlignHitQualSelTool::InDetAlignHitQualSelTool( const std::string& t + , const std::string& n + , const IInterface* p + ) + : AthAlgTool(t,n,p) + , m_rejectOutliers( true ) + , m_maxClusterSize( 5 ) + , m_rejectEdgeChannels( true ) + , m_rejectGangedPixels( false ) + , m_maxIncidAngle( 0.8 ) //!< corresponds to 45 deg. -- recomm. by pixel ppl +{ + declareInterface<IInDetAlignHitQualSelTool>(this) ; + declareProperty( "RejectOutliers", m_rejectOutliers ) ; + declareProperty( "MaxClusterSize", m_maxClusterSize ) ; + declareProperty( "RejectEdgeChannels", m_rejectEdgeChannels ) ; + declareProperty( "RejectGangedPixels", m_rejectGangedPixels ) ; + declareProperty( "MaxIncidAngle", m_maxIncidAngle ) ; +} + + +InDetAlignHitQualSelTool::~InDetAlignHitQualSelTool() {} + + +StatusCode InDetAlignHitQualSelTool::initialize() { + StatusCode sc = AlgTool::initialize() ; + if( sc.isFailure() ) return sc ; + // get DetectorStore service + sc = detStore().retrieve() ; + if( sc.isFailure() ) { + ATH_MSG_ERROR( "DetectorStore service not found !" ); + return sc ; + } else { + ATH_MSG_DEBUG( "DetectorStore retrieved!" ); + } + if (detStore()->retrieve(m_sctID, "SCT_ID").isFailure()){ + msg(MSG::FATAL) << "Could not get SCT ID helper" << endreq; + return StatusCode::FAILURE; + } + else if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "SCT ID is : "<< m_sctID <<endreq ; + + if (detStore()->retrieve(m_pixelid, "PixelID").isFailure()){ + msg(MSG::FATAL) << "Could not get PIXEL ID helper" << endreq; + return StatusCode::FAILURE; + } + else if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "Pixel ID is : " << m_pixelid << endreq; + + // get pixel manager + sc = detStore()->retrieve( m_PIXManager, "Pixel" ) ; + if( sc.isFailure() ) { + ATH_MSG_ERROR( "Could not get PIXManager !" ) ; + return sc; + } + // get SCT manager + sc = detStore()->retrieve( m_SCTManager, "SCT" ) ; + if( sc.isFailure() ) { + ATH_MSG_ERROR( "Could not get SCTManager !" ) ; + return sc; + } + ATH_MSG_DEBUG( "initialize() successful in " << name() ) ; + + + return StatusCode::SUCCESS ; +} + + +StatusCode InDetAlignHitQualSelTool::finalize() { + ATH_MSG_DEBUG( "finalize() successful in " << name() ) ; + return AlgTool::finalize() ; +} + + +const Trk::RIO_OnTrack* InDetAlignHitQualSelTool::getGoodHit( const Trk::TrackStateOnSurface* tsos ) const { + if( tsos == NULL ) { + ATH_MSG_ERROR( "0 pointer passed for TSOS!" ) ; + return 0 ; + } + if( !tsos->type(Trk::TrackStateOnSurface::Measurement) ) { + ATH_MSG_DEBUG( "not a hit, cast to MeasurementBase will fail, so reject" ) ; + return 0 ; + } + if( m_rejectOutliers && tsos->type(Trk::TrackStateOnSurface::Outlier) ) { + ATH_MSG_DEBUG( "outlier, reject" ) ; + return 0 ; + } + const Trk::MeasurementBase* measBase = tsos->measurementOnTrack() ; + if( measBase == NULL ) { + ATH_MSG_DEBUG( "tsos->measurementOnTrack() returned 0 pointer" ) ; + return 0 ; + } + const Trk::RIO_OnTrack* hit = dynamic_cast <const Trk::RIO_OnTrack*>( measBase ) ; + if( hit == NULL ) { + ATH_MSG_DEBUG( "dynamic_cast <const Trk::RIO_OnTrack*>( measBase ) returned 0 pointer" ) ; + return 0 ; + } + // cut on the cluster size + const Trk::PrepRawData* prd = hit->prepRawData() ; + if( prd == NULL ) { + ATH_MSG_WARNING( "hit->prepRawData() method failed" ) ; + return 0 ; + } + if( m_rejectGangedPixels && isGangedPixel( prd ) ) return 0 ; + const vector<Identifier> idVec = prd->rdoList() ; + if( m_maxClusterSize > 0 && !isGoodClusterSize( idVec ) ) return 0 ; + // cut on edge channels + if( m_rejectEdgeChannels && isEdgeChannel( idVec ) ) return 0 ; + // cut on the track incidence angle alpha + const Trk::TrackParameters* trkPar = tsos->trackParameters() ; + if( trkPar == NULL ) { + ATH_MSG_WARNING( "tsos->trackParameters() returned 0 pointer" ) ; + return 0 ; + } + const InDetDD::SiDetectorElement *detEle + = dynamic_cast<const InDetDD::SiDetectorElement*>( hit->detectorElement() ) ; + if( detEle == NULL ) { + ATH_MSG_WARNING( "hit cast to SiDetectorElement returned 0 pointer" ) ; + return 0 ; + } + if( !isGoodAngle( trkPar, detEle ) ) return 0 ; + return hit ; +} + + +bool InDetAlignHitQualSelTool::getGoodHole( const Trk::TrackStateOnSurface* tsos ) const { + if( tsos == NULL ) { + ATH_MSG_ERROR( "0 pointer passed for TSOS!" ) ; + return false ; + } + if( !tsos->type(Trk::TrackStateOnSurface::Hole) ) { + ATH_MSG_DEBUG( "This is not a hole, reject" ) ; + return false ; + } + // for holes only cut on the track incidence angle alpha + const Trk::TrackParameters* trkPar = tsos->trackParameters() ; + if( trkPar == NULL ) { + ATH_MSG_WARNING( "tsos->trackParameters() returned 0 pointer" ) ; + return false ; + } + const InDetDD::SiDetectorElement *detEle = dynamic_cast<const InDetDD::SiDetectorElement*>( + tsos->trackParameters()->associatedSurface().associatedDetectorElement() ) ; + if( detEle == NULL ) { + ATH_MSG_WARNING( "hole cast to SiDetectorElement returned 0 pointer" ) ; + return false ; + } + if( !isGoodAngle( trkPar, detEle ) ) return false ; + return true; +} + + +bool InDetAlignHitQualSelTool::isGangedPixel( const Trk::PrepRawData* prd ) const { + const InDet::SiCluster* cluster = dynamic_cast<const InDet::SiCluster*>( prd ) ; + if( cluster == NULL ) { + ATH_MSG_WARNING( "dynamic_cast<const InDet::SiCluster*>( prd ) failed!" ) ; + return false ; + } + if( cluster->gangedPixel() ) { //!< cut only if m_maxClusterSize set + ATH_MSG_DEBUG( "cluster contains a ganged pixel, reject" ) ; + return true ; + } + return false ; +} + + +bool InDetAlignHitQualSelTool::isGoodClusterSize( const std::vector<Identifier>& idVec ) const { + int clusterSize = idVec.size() ; + ATH_MSG_DEBUG( "clusterSize = " << clusterSize ) ; + if( clusterSize > m_maxClusterSize ) { //!< cut only if m_maxClusterSize set + ATH_MSG_DEBUG( "clusterSize = " << clusterSize << " > " << m_maxClusterSize << ", reject" ) ; + return false ; + } + return true ; +} + + +bool InDetAlignHitQualSelTool::isEdgeChannel( const vector<Identifier>& idVec ) const { + for( unsigned int i=0, i_max=idVec.size() ; i!=i_max ; ++i ) { + if( m_SCTManager->identifierBelongs(idVec[i]) ) { + int stripId = m_sctID->strip(idVec[i]) ; + if( stripId == 0 || stripId == 767 ) { + ATH_MSG_DEBUG( " SCT strip " << i << " with id " << stripId << " is an edge channel " ) ; + return true ; + } + if( stripId < 0 || stripId > 767 ) { + ATH_MSG_FATAL( " WRONG DETECTOR INFORMATION " ) ; + } + } + if( m_PIXManager->identifierBelongs(idVec[i]) ) { + int pixelIdPhi = m_pixelid->phi_index(idVec[i]) ; + int pixelIdEta = m_pixelid->eta_index(idVec[i]) ; + if( pixelIdPhi == 0 || pixelIdPhi == 327 || pixelIdEta == 0 || pixelIdEta == 143 ) { + ATH_MSG_DEBUG( " pixel hit " << i << " with idPhi " << pixelIdPhi << " and idEta " << pixelIdEta << " is an edge channel " ) ; + return true ; + } + if( pixelIdPhi < 0 || pixelIdPhi > 327 || pixelIdEta < 0 || pixelIdEta > 143 ) { + ATH_MSG_FATAL( " WRONG DETECTOR INFORMATION " ) ; + } + } + } + return false ; +} + +bool InDetAlignHitQualSelTool::isGoodAngle( const Trk::TrackParameters* trkPar + , const InDetDD::SiDetectorElement* detEle + ) const { + const double trkIncidAngle = incidAngle( trkPar, detEle ) ; + if( fabs(trkIncidAngle) > m_maxIncidAngle ) { + ATH_MSG_DEBUG( "trkIncidAngle = |" << trkIncidAngle << "| > " << m_maxIncidAngle << ", reject" ) ; + return false ; + } + return true; +} + + +double InDetAlignHitQualSelTool::incidAngle( const Trk::TrackParameters* trkPar + , const InDetDD::SiDetectorElement* detEle + ) const { + Amg::Vector3D trkDir = trkPar->momentum() ; + Amg::Vector3D detElePhi = detEle->phiAxis() ; //!< local x axis in global frame + Amg::Vector3D detEleNormal = detEle->normal() ; //!< local z axis in global frame + double trkDotPhi = trkDir.dot( detElePhi ) ; //!< scalar product + double trkDotNormal = trkDir.dot( detEleNormal ) ; + double trkIncidAngle = atan( trkDotPhi/trkDotNormal ) ; + ATH_MSG_DEBUG( "trkIncidAngle = " << trkIncidAngle ) ; + return trkIncidAngle ; +} diff --git a/InnerDetector/InDetAlignTools/InDetAlignGenTools/src/InDetAlignOverlapTool.cxx b/InnerDetector/InDetAlignTools/InDetAlignGenTools/src/InDetAlignOverlapTool.cxx new file mode 100644 index 0000000000000000000000000000000000000000..e0e34f93595afddc4eae7a5a1ea3f034185f77fc --- /dev/null +++ b/InnerDetector/InDetAlignTools/InDetAlignGenTools/src/InDetAlignOverlapTool.cxx @@ -0,0 +1,468 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +#include <iostream> +#include <iomanip> +#include <fstream> +#include <string> +#include <math.h> +#include <vector> + +#include "DataModel/DataVector.h" +#include "InDetAlignTrkInfo/AlignSiModule.h" +#include "InDetAlignTrkInfo/AlignTrk.h" +#include "InDetAlignTrkInfo/AlignTrkContainer.h" +#include "TrkRIO_OnTrack/RIO_OnTrack.h" +#include "InDetAlignTrkInfo/AlignSiHit.h" +#include "InDetAlignGenTools/InDetAlignOverlapTool.h" +//#include "InDetAlignGenAlgs/InDetAlignNt.h" +#include "InDetAlignGenTools/AlignSiModuleList.h" + + +using namespace std ; + + +InDetAlignOverlapTool::InDetAlignOverlapTool( const std::string& t , const std::string& n , const IInterface* p) + : AthAlgTool(t,n,p) + // ,m_minPixelOverlapPerTrack(1) + // ,m_minSCTOverlapPerTrack(2) + +{ + declareInterface<IInDetAlignOverlapTool>(this) ; + //declareProperty( "MinimumPixelOverlapPerTrack", m_minPixelOverlapPerTrack ) ; + //declareProperty( "MinimumSCTOverlapPerTrack", m_minSCTOverlapPerTrack ) ; +} + +InDetAlignOverlapTool::~InDetAlignOverlapTool() {} + +StatusCode InDetAlignOverlapTool::initialize() { + + StatusCode sc = AlgTool::initialize() ; + if( sc.isFailure() ) return sc ; + + // get DetectorStore service + sc = detStore().retrieve() ; + if( sc.isFailure() ) { + ATH_MSG_ERROR( "DetectorStore service not found !" ); + return sc ; + } else { + ATH_MSG_DEBUG( "DetectorStore retrieved!" ); + } + + // get pointer to AlignSiModuleList; after setup has created it + if (StatusCode::SUCCESS!=detStore()->retrieve(p_modlist,"InDetAlignNt")) + msg(MSG::WARNING) << "Could not retrieve Modlist in execute" << endreq; + + + ATH_MSG_DEBUG( "initialize() successful in " << name() ) ; + return StatusCode::SUCCESS ; +} + + +StatusCode InDetAlignOverlapTool::finalize() { + ATH_MSG_DEBUG( "finalize() successful in " << name() ) ; + return AlgTool::finalize() ; +} + +//** This function give us the number of overlaps in pixel detector + +int InDetAlignOverlapTool::getNumberOverlapPIX( const AlignTrk& trk ) const{ + ATH_MSG_DEBUG(" InDetAlignOverlapTool::get getNumberOverlapPIX inizialized"); + + int nPixelBarrelOverlapEta=0; + int nPixelBarrelOverlapPhi=0; + int nPixelECOverlapPhi=0; + int PixelSector[3] = {21,37,51}; + int ECPixelSector = 48; + int nOverlapPIX=0; + int nhits2=0; + + //-------------------- + // first loop: + //-------------------- + + nhits2=trk.nhits(); + ATH_MSG_DEBUG(" Num de Hits "<< nhits2); + ATH_MSG_DEBUG(" Anem a entrar en el loop"); + + for (std::vector<AlignSiHit>::const_iterator hit=trk.hitlist_cbegin(); + hit!=trk.hitlist_cend();++hit) { + + int index = hit->index()-1; // note index starts at 1 in ntuple! + + if(msgLvl(MSG::VERBOSE)){ + if((p_modlist->vec[index]).bec()==0 && (p_modlist->vec[index]).dettype()==1 ){ + msg().setColor(MSG::GREEN); + msg(MSG::VERBOSE)<<"Index "<< index<<" Detector: "<< (p_modlist->vec[index]).dettype()<< " layer : " <<(p_modlist->vec[index]).layer() << " ring : " + <<(p_modlist->vec[index]).ring() <<" sector : "<<(p_modlist->vec[index]).sector()<< endreq; + } + } + + //-------------------- + // second loop: + //-------------------- + for(std::vector<AlignSiHit>::const_iterator Allhit=hit; + Allhit!=trk.hitlist_cend();++Allhit) { + + int Allindex = Allhit->index()-1; // note index starts at 1 in ntuple! + + if(msgLvl(MSG::VERBOSE)){ + if((p_modlist->vec[Allindex]).bec()==0 && (p_modlist->vec[Allindex]).dettype()==1 ){ + msg().setColor(MSG::BLUE); + msg(MSG::VERBOSE)<<"Index "<< Allindex << endreq; + } + } + + //------------------------// + // PIXEL SYSTEM // + //------------------------// + + if ((p_modlist->vec[index]).dettype()==1 && (p_modlist->vec[Allindex]).dettype()==1){ + if ((p_modlist->vec[index]).bec()==0 && (p_modlist->vec[Allindex]).bec()==0){ // BARREL + + if ((p_modlist->vec[Allindex]).layer()==(p_modlist->vec[index]).layer() && (p_modlist->vec[Allindex]).ring()==(p_modlist->vec[index]).ring() && + (fabs((p_modlist->vec[Allindex]).sector()-(p_modlist->vec[index]).sector())==1. || + fabs((p_modlist->vec[Allindex]).sector()-((p_modlist->vec[index]).sector()+ PixelSector[(p_modlist->vec[index]).layer()])) ==1.)){ + + // msg() + if(msgLvl(MSG::DEBUG)){ + msg().setColor(MSG::RED); + msg(MSG::DEBUG)<<" Index:"<< index << " layer " << (p_modlist->vec[index]).layer()<<" ring "<<(p_modlist->vec[index]).ring() << + " sector "<<(p_modlist->vec[index]).sector()<<endreq;//<<"rPhi Residuals: " << hit->rphi_resid() <<" eta Residuals: " << hit->z_resid()); + msg().setColor(MSG::RED); + msg(MSG::DEBUG)<< " AllIndex:"<< Allindex <<" layer "<<(p_modlist->vec[Allindex]).layer()<<" ring "<<(p_modlist->vec[Allindex]).ring() << + " sector "<<(p_modlist->vec[Allindex]).sector() <<endreq;//<<" rPhi Residuals: " << Allhit->rphi_resid() <<" eta Residuals: " << Allhit->z_resid()); + msg().setColor(MSG::RED); + msg(MSG::DEBUG)<<"Pixel Overlap in Phi found"<<endreq; + } //msg() end + + nPixelBarrelOverlapPhi++; + } + + if ((p_modlist->vec[Allindex]).layer()==(p_modlist->vec[index]).layer() && fabs((p_modlist->vec[Allindex]).ring()-(p_modlist->vec[index]).ring())==1 && + (p_modlist->vec[Allindex]).sector()== (p_modlist->vec[index]).sector()){ + + //msg() + if(msgLvl(MSG::DEBUG)){ + msg().setColor(MSG::RED); + msg(MSG::DEBUG)<<" Index:"<< index << " layer " << (p_modlist->vec[index]).layer()<<" ring "<<(p_modlist->vec[index]).ring() << + " sector "<<(p_modlist->vec[index]).sector()<<endreq;//<<"rPhi Residuals: " << hit->rphi_resid() <<" eta Residuals: " << hit->z_resid()); + msg().setColor(MSG::RED); + msg(MSG::DEBUG)<<" AllIndex:"<< Allindex <<" layer "<<(p_modlist->vec[Allindex]).layer()<<" ring "<<(p_modlist->vec[Allindex]).ring() << + " sector "<<(p_modlist->vec[Allindex]).sector()<<endreq;// <<" rPhi Residuals: " << Allhit->rphi_resid() <<" eta Residuals: " << Allhit->z_resid()); + msg().setColor(MSG::RED); + msg(MSG::DEBUG)<<"Pixel Overlap in Eta found"<<endreq; + } //msg() end + + nPixelBarrelOverlapEta++; + } + } // END OF BARREL + + if ((p_modlist->vec[index]).bec()!=0 && (p_modlist->vec[index]).bec()==(p_modlist->vec[Allindex]).bec() ){ // ENDCAP + if ((p_modlist->vec[Allindex]).layer()==(p_modlist->vec[index]).layer() && (p_modlist->vec[Allindex]).ring()==(p_modlist->vec[index]).ring() && + (fabs((p_modlist->vec[Allindex]).sector()-(p_modlist->vec[index]).sector())==1. || + fabs((p_modlist->vec[Allindex]).sector()-((p_modlist->vec[index]).sector()+ ECPixelSector)) ==1.)){ + + //msg() + if(msgLvl(MSG::DEBUG)){ + msg().setColor(MSG::RED); + msg(MSG::DEBUG)<<" Index:"<< index << " layer " << (p_modlist->vec[index]).layer()<<" ring "<<(p_modlist->vec[index]).ring() << + " sector "<<(p_modlist->vec[index]).sector()<<endreq; + msg(MSG::DEBUG)<<" AllIndex:"<< Allindex <<" layer "<<(p_modlist->vec[Allindex]).layer()<<" ring "<<(p_modlist->vec[Allindex]).ring() << + " sector "<<(p_modlist->vec[Allindex]).sector()<<endreq; + msg(MSG::DEBUG)<<"EC Pixel Overlap in Phi found"<<endreq; + }// msg() end + + nPixelECOverlapPhi++; + } + } // END OF ENDCAP + + } + } + } + + nOverlapPIX = nPixelBarrelOverlapPhi+nPixelBarrelOverlapEta+nPixelECOverlapPhi; + ATH_MSG_DEBUG(" Total number of Pixel overlaps in eta " << nPixelBarrelOverlapEta ); + ATH_MSG_DEBUG(" Total number of Pixel overlaps in phi " << nPixelBarrelOverlapPhi ); + ATH_MSG_DEBUG(" Total number of Pixel EC overlaps in phi " << nPixelECOverlapPhi ); + + return nOverlapPIX ; + +} + +//** This function give us the number of overlaps in pixel detector. +//** The structure is the same that in the pixel case ( differs in dettype and number of sectors, etc... + +int InDetAlignOverlapTool::getNumberOverlapSCT( const AlignTrk& trk ) const{ + + ATH_MSG_DEBUG(" InDetAlignOverlapTool::get getNumberOverlapSCT inizialized"); + + // counters for overlaps eta and overlaps Phi + int nSCTBarrelOverlapEta=0; + int nSCTBarrelOverlapPhi=0; + int nSCTECOverlapPhi=0; + int nSCTECOverlapEta=0; + + int SCTSector[4] = {31,39,47,55}; + int ECSCTSector= 52; + + int nOverlapSCT=0; + int nhits2=0; + + nhits2=trk.nhits(); + ATH_MSG_DEBUG(" Num de Hits "<< nhits2); + ATH_MSG_DEBUG(" Anem a entrar en el loop"); + + //------------- + // first loop: + //------------- + + for (std::vector<AlignSiHit>::const_iterator hit=trk.hitlist_cbegin(); + hit!=trk.hitlist_cend();++hit) { + + int index = hit->index()-1; // note index starts at 1 in ntuple! + + if(msgLvl(MSG::VERBOSE)){ + if((p_modlist->vec[index]).bec()==0 && (p_modlist->vec[index]).dettype()==1 ){ + msg().setColor(MSG::GREEN); + msg(MSG::VERBOSE)<<"Index "<< index<<" Detector: "<< (p_modlist->vec[index]).dettype()<< " layer : " <<(p_modlist->vec[index]).layer() << " ring : " + <<(p_modlist->vec[index]).ring() <<" sector : "<<(p_modlist->vec[index]).sector()<< endreq; + } + } + + //------------- + // second loop: + //------------- + + for(std::vector<AlignSiHit>::const_iterator Allhit=hit; + Allhit!=trk.hitlist_cend();++Allhit) { + + int Allindex = Allhit->index()-1; // note index starts at 1 in ntuple! + + if(msgLvl(MSG::VERBOSE)){ + if((p_modlist->vec[Allindex]).bec()==0 && (p_modlist->vec[Allindex]).dettype()==1 ){ + msg().setColor(MSG::BLUE); + msg(MSG::VERBOSE)<<"Index "<< Allindex << endreq; + } + } + + //------------------------// + // SCT SYSTEM // + //------------------------// + + if ((p_modlist->vec[index]).dettype()==2 &&(p_modlist->vec[Allindex]).dettype()==2 ){ + if((p_modlist->vec[index]).bec()==0 && (p_modlist->vec[Allindex]).bec()==0){ // BARREL + + if ((p_modlist->vec[Allindex]).layer()==(p_modlist->vec[index]).layer() && (p_modlist->vec[Allindex]).ring()==(p_modlist->vec[index]).ring() && + (fabs((p_modlist->vec[Allindex]).sector()-(p_modlist->vec[index]).sector())==1.|| + fabs((p_modlist->vec[Allindex]).sector()-((p_modlist->vec[index]).sector()+ SCTSector[(p_modlist->vec[index]).layer()])) ==1.)){ + + if(msgLvl(MSG::DEBUG)){ + msg().setColor(MSG::RED); + msg(MSG::DEBUG)<<" Index:"<< index << " layer " << (p_modlist->vec[index]).layer()<<" ring "<<(p_modlist->vec[index]).ring() << + " sector "<<(p_modlist->vec[index]).sector()<<endreq; + msg().setColor(MSG::RED); + msg(MSG::DEBUG)<<" AllIndex:"<< Allindex <<" layer "<<(p_modlist->vec[Allindex]).layer()<<" ring "<<(p_modlist->vec[Allindex]).ring() << + " sector "<<(p_modlist->vec[Allindex]).sector()<< endreq; + msg().setColor(MSG::RED); + msg(MSG::DEBUG)<<"SCT Overlap in Phi found"<< endreq; + } + nSCTBarrelOverlapPhi++; + } + + if ((p_modlist->vec[Allindex]).layer()==(p_modlist->vec[index]).layer() && fabs((p_modlist->vec[Allindex]).ring()-(p_modlist->vec[index]).ring())==1 && + (p_modlist->vec[Allindex]).sector()== (p_modlist->vec[index]).sector()){ + if(msgLvl(MSG::DEBUG)){ + msg().setColor(MSG::RED); + msg(MSG::DEBUG)<<" Index:"<< index << " layer " << (p_modlist->vec[index]).layer()<<" ring "<<(p_modlist->vec[index]).ring() << + " sector "<<(p_modlist->vec[index]).sector()<<endreq; + msg().setColor(MSG::RED); + msg(MSG::DEBUG)<<" AllIndex:"<< Allindex <<" layer "<<(p_modlist->vec[Allindex]).layer()<<" ring "<<(p_modlist->vec[Allindex]).ring() << + " sector "<<(p_modlist->vec[Allindex]).sector()<<endreq; + msg().setColor(MSG::RED); + msg(MSG::DEBUG)<<"SCT Overlap in Eta found"<<endreq; + } + nSCTBarrelOverlapEta++; + } + } + + if((p_modlist->vec[index]).bec()!=0 && (( p_modlist->vec[index]).bec()== (p_modlist->vec[Allindex]).bec()) ){ // EndCap ( EC values 2 i -2) + + if ((p_modlist->vec[Allindex]).layer()==(p_modlist->vec[index]).layer() && (p_modlist->vec[Allindex]).ring()==(p_modlist->vec[index]).ring() && + (fabs((p_modlist->vec[Allindex]).sector()-(p_modlist->vec[index]).sector())==1.|| + fabs((p_modlist->vec[Allindex]).sector()-((p_modlist->vec[index]).sector()+ ECSCTSector)) ==1.)){ + if(msgLvl(MSG::DEBUG)){ + msg().setColor(MSG::RED); + msg(MSG::DEBUG)<<" Index:"<< index << " layer " << (p_modlist->vec[index]).layer()<<" ring "<<(p_modlist->vec[index]).ring() << + " sector "<<(p_modlist->vec[index]).sector()<<endreq; + msg().setColor(MSG::RED); + msg(MSG::DEBUG)<<" AllIndex:"<< Allindex <<" layer "<<(p_modlist->vec[Allindex]).layer()<<" ring "<<(p_modlist->vec[Allindex]).ring() << + " sector "<<(p_modlist->vec[Allindex]).sector()<<endreq; + + msg(MSG::DEBUG)<<"EC SCT Overlap in Phi found"<<endreq; + } + nSCTECOverlapPhi++; + } + + if ((p_modlist->vec[Allindex]).layer()==(p_modlist->vec[index]).layer() && fabs((p_modlist->vec[Allindex]).ring()-(p_modlist->vec[index]).ring())==1 && + (p_modlist->vec[Allindex]).sector()== (p_modlist->vec[index]).sector()){ + if(msgLvl(MSG::DEBUG)){ + msg().setColor(MSG::GREEN); + msg(MSG::DEBUG)<<" Index:"<< index << " layer " << (p_modlist->vec[index]).layer()<<" ring "<<(p_modlist->vec[index]).ring() << + " sector "<<(p_modlist->vec[index]).sector()<<endreq; + msg().setColor(MSG::GREEN); + msg(MSG::DEBUG)<<" AllIndex:"<< Allindex <<" layer "<<(p_modlist->vec[Allindex]).layer()<<" ring "<<(p_modlist->vec[Allindex]).ring() << + " sector "<<(p_modlist->vec[Allindex]).sector()<<endreq; + + msg(MSG::DEBUG)<<"EC SCT Overlap in Eta found"<<endreq; + } + nSCTECOverlapEta++; + } + + } + } + + } + } + + nOverlapSCT = nSCTBarrelOverlapPhi + nSCTBarrelOverlapEta + nSCTECOverlapPhi + nSCTECOverlapEta; + + ATH_MSG_DEBUG(" Total number of STC overlaps in eta " << nSCTBarrelOverlapEta ); + ATH_MSG_DEBUG(" Total number of SCT overlaps in phi " << nSCTBarrelOverlapPhi ); + ATH_MSG_DEBUG(" Total number of STC EC overlaps in phi " << nSCTECOverlapPhi ); + ATH_MSG_DEBUG(" Total number of STC EC overlaps in Eta " << nSCTECOverlapEta ); + + return nOverlapSCT ; + +} + + + +//** This fuctions give a DataVector<AlignSiHit> with overlap hits in aligntrk +//** Now we take into account all hits in 2 loops + +std::vector<AlignSiHit> InDetAlignOverlapTool::getOverlapHit( const AlignTrk& trk ) { //const { + // const std::vector<AlignSiHit> InDetAlignOverlapTool::getOverlapHit( const AlignTrk& trk ){//const { + + + //m_Overlaphits.clear(); + + //const DataVector<AlignSiHit>* InDetAlignOverlapTool::getOverlapHit( const AlignTrk& trk ) const{ + //DataVector<AlignSiHit>* m_Overlaphits= new DataVector<AlignSiHit>; + + //std::vector<AlignSiHit> m_Overlaphits= new std::vector<AlignSiHit>; + + ATH_MSG_DEBUG(" InDetAlignOverlapTool::get getOverlapHit inizialized"); + + bool isPIXOverlap =false; + bool isSCTOverlap =false; + bool isOverlap=false; + + int PixelSector[3] = {21,37,51}; + int ECPixelSector = 48; + int SCTSector[4] = {31,39,47,55}; + int ECSCTSector= 52; + + //---------------------------------------------------- + // first loop: preparing info for next processing + //---------------------------------------------------- + for (std::vector<AlignSiHit>::const_iterator hit=trk.hitlist_cbegin(); + hit!=trk.hitlist_cend();++hit) { + + int index = hit->index()-1; // note index starts at 1 in ntuple! + + isPIXOverlap =false; + isSCTOverlap =false; + isOverlap=false; + + //---------------------------------------------------- + // Second loop: preparing info for next processing + //---------------------------------------------------- + for (std::vector<AlignSiHit>::const_iterator hit2=trk.hitlist_cbegin(); + hit2!=trk.hitlist_cend();++hit2) { + + int newindex = hit2->index()-1; // note index starts at 1 in ntuple! + + //msg().setColor(MSG::GREEN); + // ATH_MSG_DEBUG("Hit : "<< newindex); + + //------------------------// + // PIXEL SYSTEM // + //------------------------// + + if ((p_modlist->vec[index]).dettype()==1 && (p_modlist->vec[newindex]).dettype()==1){ + if ((p_modlist->vec[index]).bec()==0 && (p_modlist->vec[newindex]).bec()==0){ // BARREL + + if ((p_modlist->vec[newindex]).layer()==(p_modlist->vec[index]).layer() && (p_modlist->vec[newindex]).ring()==(p_modlist->vec[index]).ring() && + (fabs((p_modlist->vec[newindex]).sector()-(p_modlist->vec[index]).sector())==1. || + fabs((p_modlist->vec[newindex]).sector()-((p_modlist->vec[index]).sector()+ PixelSector[(p_modlist->vec[index]).layer()])) ==1.)){ + isPIXOverlap=true; + } + + if ((p_modlist->vec[newindex]).layer()==(p_modlist->vec[index]).layer() && fabs((p_modlist->vec[newindex]).ring()-(p_modlist->vec[index]).ring())==1 && + (p_modlist->vec[newindex]).sector()== (p_modlist->vec[index]).sector()){ + isPIXOverlap=true; + } + } + + if ((p_modlist->vec[index]).bec()!=0 && (p_modlist->vec[index]).bec()==(p_modlist->vec[newindex]).bec() ){ // ENDCAP + if ((p_modlist->vec[newindex]).layer()==(p_modlist->vec[index]).layer() && (p_modlist->vec[newindex]).ring()==(p_modlist->vec[index]).ring() && + (fabs((p_modlist->vec[newindex]).sector()-(p_modlist->vec[index]).sector())==1. || + fabs((p_modlist->vec[newindex]).sector()-((p_modlist->vec[index]).sector()+ ECPixelSector)) ==1.)){ + isPIXOverlap=true; + } + } + } // End Pixel system + + + //------------------------// + // SCT SYSTEM // + //------------------------// + + if ((p_modlist->vec[index]).dettype()==2 &&(p_modlist->vec[newindex]).dettype()==2 ){ + if((p_modlist->vec[index]).bec()==0 && (p_modlist->vec[newindex]).bec()==0){ // BARREL + + if ((p_modlist->vec[newindex]).layer()==(p_modlist->vec[index]).layer() && (p_modlist->vec[newindex]).ring()==(p_modlist->vec[index]).ring() && + (fabs((p_modlist->vec[newindex]).sector()-(p_modlist->vec[index]).sector())==1.|| + fabs((p_modlist->vec[newindex]).sector()-((p_modlist->vec[index]).sector()+ SCTSector[(p_modlist->vec[index]).layer()])) ==1.)){ + + isSCTOverlap=true; + } + + if ((p_modlist->vec[newindex]).layer()==(p_modlist->vec[index]).layer() && fabs((p_modlist->vec[newindex]).ring()-(p_modlist->vec[index]).ring())==1 && + (p_modlist->vec[newindex]).sector()== (p_modlist->vec[index]).sector()){ + isSCTOverlap=true; + } + } + + if((p_modlist->vec[index]).bec()!=0 && (( p_modlist->vec[index]).bec()== (p_modlist->vec[newindex]).bec()) ){ // // ENDCAP ( EC values 2 i -2) + + if ((p_modlist->vec[newindex]).layer()==(p_modlist->vec[index]).layer() && (p_modlist->vec[newindex]).ring()==(p_modlist->vec[index]).ring() && + (fabs((p_modlist->vec[newindex]).sector()-(p_modlist->vec[index]).sector())==1.|| + fabs((p_modlist->vec[newindex]).sector()-((p_modlist->vec[index]).sector()+ ECSCTSector)) ==1.)){ + + isSCTOverlap=true; + } + + if ((p_modlist->vec[newindex]).layer()==(p_modlist->vec[index]).layer() && fabs((p_modlist->vec[newindex]).ring()-(p_modlist->vec[index]).ring())==1 && + (p_modlist->vec[newindex]).sector()== (p_modlist->vec[index]).sector()){ + isSCTOverlap=true; + } + + + } // End SCT system + + } + + + } + if (isPIXOverlap || isSCTOverlap) { + msg().setColor(MSG::GREEN); + ATH_MSG_DEBUG("This Hit is overlap"<< index ); + m_Overlaphits.push_back(*hit); + } + } + + return m_Overlaphits; +} diff --git a/InnerDetector/InDetAlignTools/InDetAlignGenTools/src/InDetAlignTrackSelTool.cxx b/InnerDetector/InDetAlignTools/InDetAlignGenTools/src/InDetAlignTrackSelTool.cxx new file mode 100755 index 0000000000000000000000000000000000000000..1a790fe3106ed23df14e6a3c55435a783e3aa5b6 --- /dev/null +++ b/InnerDetector/InDetAlignTools/InDetAlignGenTools/src/InDetAlignTrackSelTool.cxx @@ -0,0 +1,254 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +// InDetAlignTrackSelTool.h + +// Sergio Gonzalez Sevilla, started 04/7/05 +// Miguel Olivo Gomez, extended 07/6/06 + +// AlgTool to select high quality tracks for the inner detector +// (Pixel+SCT) alignment algorithms. +// This AlgTool provides a track selection based on the following cut variables: +// Momentum, pt, number of shared hits, number of holes and chi2 probability. +// Returns 0 in case a track is not accepted, otherwise 1 + +#include "InDetAlignGenTools/InDetAlignTrackSelTool.h" + +#include "TrkTrackSummary/TrackSummary.h" +#include "TrkToolInterfaces/ITrackSummaryTool.h" + +#include "TrkParameters/TrackParameters.h" +#include "TrkEventPrimitives/FitQuality.h" +#include "CLHEP/GenericFunctions/CumulativeChiSquare.hh" + +InDetAlignTrackSelTool::InDetAlignTrackSelTool( const std::string& type + , const std::string& name + , const IInterface* parent + ) + : AthAlgTool(type,name,parent), + m_trackSumTool("Trk::TrackSummaryTool", this), + m_minMomentum(0), + m_minPt(2), + m_maxShared(0), + m_maxHoles(1), + m_minChi2Prob(0.2) +{ + declareInterface<IInDetAlignTrackSelTool>(this); + declareProperty("MinMomentum" , m_minMomentum); + declareProperty("MinPt" , m_minPt); + declareProperty("MaxShared" , m_maxShared); + declareProperty("MaxHoles" , m_maxHoles); + declareProperty("MinChi2Prob" , m_minChi2Prob); + + // Tools + declareProperty("TrackSummaryTool" , m_trackSumTool, + "tool to extract track info"); +} + +InDetAlignTrackSelTool::~InDetAlignTrackSelTool() +{} + +//////////////////////////////////////////////////////////////////////////////////////// +StatusCode InDetAlignTrackSelTool::initialize(){ + //////////////////////////////////////////////////////////////////////////////////////// + + // get TrackSummaryTool + if ( m_trackSumTool.retrieve().isFailure() ) { + ATH_MSG_FATAL( "Failed to retrieve tool " << m_trackSumTool ) ; + return StatusCode::FAILURE; + } else ATH_MSG_DEBUG( "Retrieved tool " << m_trackSumTool ) ; + + ATH_MSG_DEBUG( "Cuts selected : min_Momentum(CLHEP::GeV)=" << m_minMomentum + << " min_pt(CLHEP::GeV)=" << m_minPt + << " max_shared=" << m_maxShared + << " max_holes=" << m_maxHoles + << " min_chi2Prob=" << m_minChi2Prob ) ; + + ATH_MSG_DEBUG( "InDetAlignTrackSelTool initialize() successful" ) ; + return StatusCode::SUCCESS; +} + +//////////////////////////////////////////////////////////////////////////////////////// +StatusCode InDetAlignTrackSelTool::finalize(){ + //////////////////////////////////////////////////////////////////////////////////////// + ATH_MSG_DEBUG( "InDetAlignTrackSelTool finalize method called" ) ; + return StatusCode::SUCCESS; +} + +//////////////////////////////////////////////////////////////////////////////////////// +double InDetAlignTrackSelTool::Momentum(const Trk::Track& track) const { + //////////////////////////////////////////////////////////////////////////////////////// + ATH_MSG_DEBUG( "in Momentum() " ) ; + double mom=0.; + + // get measured perigee and momentum of track + const Trk::Perigee* perigee = track.perigeeParameters(); + + if ( !perigee->covariance()) { + ATH_MSG_ERROR( "No measured perigee parameters assigned to the track" ) ; + mom = -1e12; // return big value + } + else{ + Amg::VectorX perigeeParams = perigee->parameters(); + mom = fabs(1./perigeeParams[Trk::qOverP]); + mom /= 1000.; //mom in GeV + } + + return mom; +} + +//////////////////////////////////////////////////////////////////////////////////////// +double InDetAlignTrackSelTool::Pt(const Trk::Track& track) const { + //////////////////////////////////////////////////////////////////////////////////////// + ATH_MSG_DEBUG( "in Pt() " ) ; + double pt=0.; + + // get measured perigee and pt of track + const Trk::Perigee* perigee = track.perigeeParameters(); + + if (!perigee->covariance()) { + ATH_MSG_ERROR( "No measured perigee parameters assigned to the track" ) ; + pt = -1e12; // return big value + } + else{ + Amg::VectorX perigeeParams = perigee->parameters(); + pt = fabs(sin(perigeeParams[Trk::theta])/perigeeParams[Trk::qOverP]); + pt /= 1000.; // pt in GeV + } + + return pt; +} + +//////////////////////////////////////////////////////////////////////////////////////// +int InDetAlignTrackSelTool::nShared(const Trk::Track& track) const { + //////////////////////////////////////////////////////////////////////////////////////// + ATH_MSG_DEBUG( "in nShared()" ) ; + int nshared=0, nshpix, nshsct; + + const Trk::TrackSummary* summary = m_trackSumTool->createSummary(track); + if(summary==0){ + ATH_MSG_ERROR( "Could not create TrackSummary" ) ; + nshared = 1000; + } + else{ + nshpix = summary->get(Trk::numberOfPixelSharedHits); + nshsct = summary->get(Trk::numberOfSCTSharedHits); + + if(nshpix==-1) + nshpix=0; + + if(nshsct==-1) + nshsct=0; + + nshared = nshpix + nshsct; + } + + return nshared; +} + +//////////////////////////////////////////////////////////////////////////////////////// +int InDetAlignTrackSelTool::nHoles(const Trk::Track& track) const { + //////////////////////////////////////////////////////////////////////////////////////// + ATH_MSG_DEBUG( "in nHoles() " ) ; + int nholes=0, nhpix, nhsct; + + const Trk::TrackSummary* summary = m_trackSumTool->createSummary(track); + if(summary==0){ + ATH_MSG_ERROR( "Could not create TrackSummary" ) ; + nholes = 1000; + } + else{ + ATH_MSG_VERBOSE( "TrackSummary created. Dumping it..." ) ; + ATH_MSG_VERBOSE( (*summary) ) ; + + nhpix = summary->get(Trk::numberOfPixelHoles); + nhsct = summary->get(Trk::numberOfSCTHoles); + + if(nhpix==-1) + nhpix = 0; + + if(nhsct==-1) + nhsct = 0; + + nholes = nhpix + nhsct; + } + + return nholes; +} + +//////////////////////////////////////////////////////////////////////////////////////// +double InDetAlignTrackSelTool::chi2Prob(const Trk::Track& track) const { + //////////////////////////////////////////////////////////////////////////////////////// + ATH_MSG_DEBUG( "in chi2Prob()" ) ; + double chi2Prob=0.; + + // get fit quality and chi2 probability of track + // chi2Prob = TMath::Prob(chi2,DoF) ROOT function + const Trk::FitQuality* fitQual = track.fitQuality(); + + if (fitQual==0) { + ATH_MSG_ERROR( "No fit quality assigned to the track" ) ; + chi2Prob = -1e12; // return big value + } + else { + if (fitQual->chiSquared() > 0. && fitQual->numberDoF() > 0) { + Genfun::CumulativeChiSquare probabilityFunction( fitQual->numberDoF() ); + chi2Prob = 1 - probabilityFunction( fitQual->chiSquared() ); + } + } + + return chi2Prob; +} + +//////////////////////////////////////////////////////////////////////////////////////// +int InDetAlignTrackSelTool::getStatus(const Trk::Track& track) const { + //////////////////////////////////////////////////////////////////////////////////////// + ATH_MSG_DEBUG( "in getStatus()" ) ; + int stat=1, nholes, nshared; + double mom, pt, chi2prob; + + // momentum + mom = Momentum(track); + if (mom < m_minMomentum) { + stat=0; + } + + // transverse momentum + pt = Pt(track); + if (pt < m_minPt) { + stat=0; + } + + // number of holes + nholes = nHoles(track); + if (nholes > m_maxHoles) { + stat=0; + } + + // number of shared hits + nshared = nShared(track); + if (nshared > m_maxShared) { + stat=0; + } + + // chi2 Probability + chi2prob = chi2Prob(track); + if (chi2prob < m_minChi2Prob) { + stat=0; + } + + ATH_MSG_DEBUG( " momentum(CLHEP::GeV)=" << mom + << " pt (CLHEP::GeV)=" << pt + << " nshared=" << nshared + << " nholes=" << nholes + << " chi2Prob=" << chi2prob + ) ; + + if(!stat) + ATH_MSG_DEBUG( "Track not accepted" ) ; + else + ATH_MSG_DEBUG( "Track accepted" ) ; + + return stat; +} diff --git a/InnerDetector/InDetAlignTools/InDetAlignGenTools/src/RefitSiOnlyTool.cxx b/InnerDetector/InDetAlignTools/InDetAlignGenTools/src/RefitSiOnlyTool.cxx new file mode 100755 index 0000000000000000000000000000000000000000..bcacac3e390e095e2222e5d6c88a66665ce91c0d --- /dev/null +++ b/InnerDetector/InDetAlignTools/InDetAlignGenTools/src/RefitSiOnlyTool.cxx @@ -0,0 +1,250 @@ +/* + Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration +*/ + +////////////////////////////////////////////////////////////////////////// +// ================================================ +// RefitSiOnly +// ================================================ +// +// RefitSiOnlyTool.cxx +// Source file for RefitSiOnlyTool +// +// namespace InDetAlignment +// +// AlgTool to create refitted tracks with only silicon hits present + +#include "DataModel/DataVector.h" + +#include "Identifier/Identifier.h" +#include "AtlasDetDescr/AtlasDetectorID.h" + +//Track-things +#include "TrkTrack/Track.h" +#include "TrkRIO_OnTrack/RIO_OnTrack.h" +#include "TrkPrepRawData/PrepRawData.h" +#include "TrkMeasurementBase/MeasurementBase.h" +#include "TrkParameters/TrackParameters.h" + +#include "TrkFitterInterfaces/ITrackFitter.h" +#include "TrkEventUtils/TrkParametersComparisonFunction.h" + +#include "InDetPrepRawData/SiCluster.h" + +#include "InDetAlignGenTools/RefitSiOnlyTool.h" + +namespace InDetAlignment { + + RefitSiOnlyTool::RefitSiOnlyTool (const std::string& type, const std::string& name, const IInterface* parent): + AthAlgTool(type,name,parent), + m_ITrkFitter("Trk::KalmanFitter"), + m_comPar(0), + m_ParticleHypothesis(Trk::nonInteracting) + { + declareInterface<ICollectionProcessor>(this); + declareProperty("OutlierRemoval" , m_OutlierRemoval = true); + declareProperty("FitterTool" , m_ITrkFitter); + declareProperty("ParticleNumber" , m_ParticleNumber = 2); + declareProperty("Cosmic" , m_doCosmic = false); + declareProperty("StripPixelHits" , m_stripPixelHits = false); + declareProperty("StripSCTHits" , m_stripSCTHits = false); + declareProperty("StripTRTHits" , m_stripTRTHits = true); + } + //_________________________________________________________________________________________ + + //destructor + RefitSiOnlyTool::~RefitSiOnlyTool() {} + //_________________________________________________________________________________________ + + StatusCode RefitSiOnlyTool::initialize() + { + StatusCode sc; + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "Initialize() of RefitSiOnlyTool" << endreq; + + // make for cosmic tracks new reference point + Amg::Vector3D refGP(0.,15000.,0.); + m_comPar = new Trk::TrkParametersComparisonFunction(refGP); + + //ID Helper + if (detStore()->retrieve(m_idHelper,"AtlasID").isFailure()) { + if (msgLvl(MSG::ERROR)) msg(MSG::ERROR) << "Unable to initialize ID helper." << endreq; + return StatusCode::FAILURE; + } + + //Fitter + if (m_ITrkFitter.retrieve().isFailure()) { + if (msgLvl(MSG::FATAL)) msg(MSG::FATAL) << "Can not retrieve Fitter of type " + << m_ITrkFitter.typeAndName() << endreq; + return StatusCode::FAILURE; + } else if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "Retrieved tool " << m_ITrkFitter.typeAndName() << endreq; + + if (m_ParticleNumber == 0) m_ParticleHypothesis = Trk::nonInteracting; + else if (m_ParticleNumber == 1) m_ParticleHypothesis = Trk::electron; + else if (m_ParticleNumber == 2) m_ParticleHypothesis = Trk::muon; + else if (m_ParticleNumber == 3) m_ParticleHypothesis = Trk::pion; + else if (m_ParticleNumber == 4) m_ParticleHypothesis = Trk::kaon; + else if (m_ParticleNumber == 5) m_ParticleHypothesis = Trk::proton; + else if (m_ParticleNumber == 99) m_ParticleHypothesis = Trk::noHypothesis; + + if (msgLvl(MSG::DEBUG)) { + msg(MSG::DEBUG) << "ParticleNumber: " << m_ParticleNumber << endreq; + msg(MSG::DEBUG) << "ParticleHypothesis: " << m_ParticleHypothesis << endreq; + msg(MSG::DEBUG) << "Initialize() of RefitSiOnlyTool successful" << endreq; + } + + return StatusCode::SUCCESS; + } + //_________________________________________________________________________________________ + + StatusCode RefitSiOnlyTool::finalize() + { + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "Finalize() of RefitSiOnlyTool" << endreq; + return StatusCode::SUCCESS; + } + //_________________________________________________________________________________________ + StatusCode RefitSiOnlyTool::configure() + { + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "nothing to configure - everything has already been configured in initialize()" << endreq; + return StatusCode::SUCCESS; + } + //_________________________________________________________________________________________ + const DataVector<Trk::Track>* RefitSiOnlyTool::processTrackCollection(const DataVector<Trk::Track>* trks) { + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "processTrackCollection() of RefitSiOnlyTool" << endreq; + + m_trkindex=0; + DataVector<Trk::Track>* newTrks = new DataVector<Trk::Track>; + + for (DataVector<Trk::Track>::const_iterator t=trks->begin(); t!=trks->end();++t) { + if (*t) { + Trk::Track* refittedTrack = refit(*t); + + if (refittedTrack==0) + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << " Refit of the track " << m_trkindex << " did not work. Track skipped." << endreq; + + newTrks->push_back(refittedTrack); + } else { + newTrks->push_back(0); + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << " Track " << m_trkindex << " is empty." << endreq; + } + + m_trkindex++; + } + + if (msgLvl(MSG::DEBUG)) { + dumpTrackCol(newTrks); + } + + return newTrks; + } + + //_________________________________________________________________________________________ + Trk::Track* RefitSiOnlyTool::refit(const Trk::Track* tr) const { + // Step through all hits on a Trk::Track, pushback SiliconHits into RIO_Collection + // and Refit Track with SiliconHits + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "Entering RefitSiOnlyTool::refit() for track #" << m_trkindex << endreq; + + const Trk::Perigee* initialPerigee = tr->perigeeParameters(); + if (initialPerigee) { + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "Track # " << m_trkindex << " has momentum before refit: " + << initialPerigee->momentum().mag() / CLHEP::GeV << " CLHEP::GeV/c" << endreq; + } + else { if (msgLvl(MSG::WARNING)) msg(MSG::WARNING) << "Track # " << m_trkindex << " without perigee" << endreq; } + + Trk::Track* SiOnlyTrack = 0; + bool containsGangedPixels = false; + std::vector<const Trk::MeasurementBase*> MeasurementBase_Collection; + + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "Building stripped RIO_Collection" << endreq; + for (DataVector<const Trk::MeasurementBase>::const_iterator measBase = tr->measurementsOnTrack()->begin(); + measBase !=tr->measurementsOnTrack()->end(); ++measBase) { + + const Trk::RIO_OnTrack* rio = dynamic_cast <const Trk::RIO_OnTrack*>(*measBase); + if (rio != 0) { + const Identifier& surfaceID = (rio->identify()); + if(filterHit(surfaceID)==false) { + MeasurementBase_Collection.push_back(*measBase); + + // protect against ganged Pixel Hits + if (rio->prepRawData()) { + const InDet::SiCluster* siHit = dynamic_cast <const InDet::SiCluster*>(rio->prepRawData()); + if (siHit) { + if (siHit->gangedPixel()) { + containsGangedPixels = true; + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << " Reject track # " << m_trkindex + << " because it contains ganged pixel hits" << endreq; + } + } + } + + } + } + } + + // Do Silicon only Refit with Silicon only HitCollection and TrackParameters at StartingPoint of track + if (!containsGangedPixels) { + if (m_doCosmic) { + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "Refitting using Cosmic reference point" << endreq; + + // order the parameters according to the Cosmics reference points + // works only for Si, not for TRT! + const Trk::TrackParameters* minPar = *(std::min_element(tr->trackParameters()->begin(), + tr->trackParameters()->end(), + *m_comPar)); + + SiOnlyTrack = m_ITrkFitter->fit(MeasurementBase_Collection, *minPar, m_OutlierRemoval, m_ParticleHypothesis); + } + else + SiOnlyTrack = m_ITrkFitter->fit(MeasurementBase_Collection, *(tr->trackParameters()->front()), + m_OutlierRemoval, m_ParticleHypothesis); + } + else { + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << " No refit was done for track # " << m_trkindex << endreq; + return 0; + } + + return SiOnlyTrack; + } + + + //_________________________________________________________________________________________ + bool RefitSiOnlyTool::filterHit(const Identifier& id) const + { + if (m_idHelper->is_sct(id) && m_stripSCTHits) {return true;} + else if (m_idHelper->is_pixel(id) && m_stripPixelHits) {return true;} + else if (m_idHelper->is_trt(id) && m_stripTRTHits) {return true;} + else {return false;} + } + + + //_________________________________________________________________________________________ + void RefitSiOnlyTool::dumpTrackCol(DataVector<Trk::Track>* tracks) { + if (msgLvl(MSG::DEBUG)) msg(MSG::DEBUG) << "In dumpTrackCol() with size=" << tracks->size() << endreq; + + int itrk = 0; + for (DataVector<Trk::Track>::const_iterator it=tracks->begin(); + it!=tracks->end();++it) { + + if(*it!=0){ + + const Trk::Perigee* aMeasPer = (*it)->perigeeParameters(); + if (aMeasPer==0){ + if (msgLvl(MSG::ERROR)) msg(MSG::ERROR) << "Could not get Trk::MeasuredPerigee" << endreq;} + else { + double d0 = aMeasPer->parameters()[Trk::d0]; + double z0 = aMeasPer->parameters()[Trk::z0]; + double phi0 = aMeasPer->parameters()[Trk::phi0]; + double theta = aMeasPer->parameters()[Trk::theta]; + double qOverP = aMeasPer->parameters()[Trk::qOverP]; + if (msgLvl(MSG::DEBUG)) { + msg(MSG::DEBUG) << itrk << ". Track Parameters (d0,z0,phi0,theta,qOverP)" << endreq; + msg(MSG::DEBUG) << " " << d0 << ", " << z0 << ", " + << phi0 << ", " << theta << ", " << qOverP << endreq; + msg(MSG::DEBUG) << " - momentum: " + << aMeasPer->momentum().mag() / CLHEP::GeV << " CLHEP::GeV/c" << endreq; + } + } + } + itrk++; + } + } +} //namespace diff --git a/InnerDetector/InDetAlignTools/InDetAlignGenTools/src/components/InDetAlignGenTools_entries.cxx b/InnerDetector/InDetAlignTools/InDetAlignGenTools/src/components/InDetAlignGenTools_entries.cxx new file mode 100755 index 0000000000000000000000000000000000000000..66d5f6b463623f4731b0400c327e32373ba484ba --- /dev/null +++ b/InnerDetector/InDetAlignTools/InDetAlignGenTools/src/components/InDetAlignGenTools_entries.cxx @@ -0,0 +1,30 @@ +#include "InDetAlignGenTools/InDetAlignDBTool.h" +#include "InDetAlignGenTools/InDetAlignTrackSelTool.h" +#include "InDetAlignGenTools/InDetAlignFillTrack.h" +//** +#include "InDetAlignGenTools/InDetAlignOverlapTool.h" +#include "InDetAlignGenTools/InDetAlignFillSiCluster.h" +#include "InDetAlignGenTools/RefitSiOnlyTool.h" +#include "InDetAlignGenTools/InDetAlignHitQualSelTool.h" +#include "GaudiKernel/DeclareFactoryEntries.h" + + +DECLARE_TOOL_FACTORY( InDetAlignDBTool ) +DECLARE_TOOL_FACTORY( InDetAlignTrackSelTool ) +DECLARE_TOOL_FACTORY( InDetAlignFillTrack ) +//** +DECLARE_TOOL_FACTORY( InDetAlignOverlapTool ) +DECLARE_TOOL_FACTORY( InDetAlignFillSiCluster ) +DECLARE_TOOL_FACTORY( InDetAlignHitQualSelTool ) +DECLARE_NAMESPACE_TOOL_FACTORY( InDetAlignment, RefitSiOnlyTool ) + +DECLARE_FACTORY_ENTRIES(InDetAlignGenTools) { + DECLARE_TOOL( InDetAlignDBTool ); + DECLARE_TOOL( InDetAlignTrackSelTool ); + DECLARE_TOOL( InDetAlignFillTrack ); + //** + DECLARE_TOOL( InDetAlignOverlapTool ); + DECLARE_TOOL( InDetAlignFillSiCluster ); + DECLARE_TOOL( InDetAlignHitQualSelTool ); + DECLARE_NAMESPACE_TOOL( InDetAlignment, RefitSiOnlyTool ); +} diff --git a/InnerDetector/InDetAlignTools/InDetAlignGenTools/src/components/InDetAlignGenTools_load.cxx b/InnerDetector/InDetAlignTools/InDetAlignGenTools/src/components/InDetAlignGenTools_load.cxx new file mode 100755 index 0000000000000000000000000000000000000000..1735bb29256a304d3fa0a7509cc15340d5a8cbb3 --- /dev/null +++ b/InnerDetector/InDetAlignTools/InDetAlignGenTools/src/components/InDetAlignGenTools_load.cxx @@ -0,0 +1,3 @@ +#include "GaudiKernel/LoadFactoryEntries.h" + +LOAD_FACTORY_ENTRIES( InDetAlignGenTools )