Skip to content
Snippets Groups Projects

Draft: Ne

Merged Ho Yin Derek Yeung requested to merge ne into master
1 file
+ 73
8
Compare changes
  • Side-by-side
  • Inline
@@ -23,6 +23,7 @@
#include "VPDet/DeVP.h"
#include "fmt/format.h"
#include <Gaudi/Accumulators/Histogram.h>
#include "GaudiAlg/GaudiTupleAlg.h"
namespace LHCb::Tr::Monitor {
@@ -30,7 +31,6 @@ namespace LHCb::Tr::Monitor {
template <typename TNode>
inline const Gaudi::TrackProjectionMatrix1D& projectionMatrix( const TNode& node ) {
// const Gaudi::TrackProjectionMatrix1D& Hk;
if constexpr ( std::is_same<TNode, LHCb::Pr::Tracks::Fit::Node>::value ) { return node.projectionMatrix(); }
if constexpr ( std::is_same<TNode, LHCb::FitNode>::value ) {
return node.template visit_r<const Gaudi::TrackProjectionMatrix1D&>(
@@ -44,7 +44,6 @@ namespace LHCb::Tr::Monitor {
template <typename TNode>
inline LHCb::State& state( const TNode& node ) {
// const Gaudi::TrackProjectionMatrix1D& Hk;
if constexpr ( std::is_same<TNode, LHCb::Pr::Tracks::Fit::Node>::value ) {
return LHCb::Pr::Tracks::Fit::state( node );
}
@@ -69,10 +68,13 @@ namespace LHCb::Tr::Monitor {
} // namespace
struct TrackVPOverlapMonitor final : LHCb::Algorithm::Consumer<void( const LHCb::Track::Range& )> {
struct TrackVPOverlapMonitor final : LHCb::Algorithm::Consumer<void( const LHCb::Track::Range& , DeVP const& ) , DetDesc::usesBaseAndConditions<GaudiTupleAlg, DeVP>> {
Gaudi::Property<bool> Cuts{this, "Cuts", true, "Adding the cuts to reduce bias in the overlap residuals"};
Gaudi::Property<bool> m_expertMode{this, "VerboseMode", false, "Cuts"};
/// Standard constructor
TrackVPOverlapMonitor( const std::string& name, ISvcLocator* pSvcLocator )
: Consumer{name, pSvcLocator, {KeyValue{"TrackContainer", TrackLocation::Default}}} {}
: Consumer{name, pSvcLocator, { KeyValue{"TrackContainer", TrackLocation::Default} , KeyValue{"VPDetectorLocation", DeVPLocation::Default} } } {}
private:
// define histograms which don't need string formatting
@@ -81,6 +83,9 @@ namespace LHCb::Tr::Monitor {
mutable Gaudi::Accumulators::Histogram<2> x_y_tileoverlap_histo{
this, "x vs y tile overlap", "x vs y tile overlap", {100, -80.0, 80.0, "X"}, {100, -80.0, 80.0, "Y"}};
mutable Gaudi::Accumulators::Histogram<1> module_histo{this, "module", "module", {52, 0.0, 52.0}};
mutable Gaudi::Accumulators::Histogram<1> module_CLI_NLO_histo{this, "CLI-NLO overlap per module", "CLI-NLO overlap per module", {52, 0.0, 52.0}};
mutable Gaudi::Accumulators::Histogram<1> module_CLI_NSI_histo{this, "CLI-NSI overlap per module", "CLI-NSI overlap per module", {52, 0.0, 52.0}};
mutable Gaudi::Accumulators::Histogram<1> module_CSO_NSI_histo{this, "CSO-NSI overlap per module", "CSO-NSI overlap per module", {52, 0.0, 52.0}};
mutable Gaudi::Accumulators::Histogram<1> station_histo{this, "station", "station", {27, 0.0, 27.0}};
mutable Gaudi::Accumulators::Histogram<2> m_overlap_residual_x_CLI_NLO{
@@ -132,24 +137,24 @@ namespace LHCb::Tr::Monitor {
public:
///< Algorithm execution
void operator()( const LHCb::Track::Range& tracks ) const override {
void operator()( const LHCb::Track::Range& tracks, const DeVP& det ) const override {
for ( const Track* track : tracks ) { // start of track-loop
if ( track->hasVelo() && track->checkFitStatus( Track::FitStatus::Fitted ) ) {
// dispatch based on trackfitresult type
const auto prfr = dynamic_cast<const LHCb::PrKalmanFitResult*>( track->fitResult() );
if ( prfr )
fill( *track, *prfr );
fill( *track, *prfr, det );
else {
const auto fr = dynamic_cast<const LHCb::TrackFitResult*>( track->fitResult() );
if ( fr ) fill( *track, *fr );
if ( fr ) fill( *track, *fr, det );
}
}
}
}
template <typename TFitResult>
void fill( const LHCb::Track& track, const TFitResult& fr ) const {
void fill( const LHCb::Track& track, const TFitResult& fr, const DeVP& det ) const {
using TNode = typename TFitResult::NodeType;
// As it turns out, it doesn't
@@ -262,7 +267,9 @@ namespace LHCb::Tr::Monitor {
// residuals
auto jnode = xynodes.begin();
auto inode = jnode;
int HitCounter = 0;
for ( ++jnode; jnode < xynodes.end(); ++jnode, ++inode ) {
HitCounter += 1;
// station overlaps
if ( jnode->vpid.station() == inode->vpid.station() ) {
const auto iside = inode->vpid.sidepos();
@@ -295,18 +302,76 @@ namespace LHCb::Tr::Monitor {
const auto globalpos = state( *( node1.xnode ) ).position();
++x_y_tileoverlap_histo[{globalpos.x(), globalpos.y()}];
const DeVPSensor& sens1 = det.sensor( node1.vpid );
const DeVPSensor& sens2 = det.sensor( node2.vpid );
const auto nodeGlobal1 = state( *( node1.xnode ) ).position();
const auto nodeGlobal2 = state( *( node2.xnode ) ).position();
const auto nodeLocal1 = sens1.globalToLocal( Gaudi::XYZPoint{nodeGlobal1.x(), nodeGlobal1.y(), nodeGlobal1.z()} );
const auto nodeLocal2 = sens2.globalToLocal( Gaudi::XYZPoint{nodeGlobal2.x(), nodeGlobal2.y(), nodeGlobal2.z()} );
if ( Cuts ) {
bool Within_One = ( node1.vpid.scol() == 0 || node1.vpid.scol() == 767 || node2.vpid.scol() == 0 || node2.vpid.scol() == 767 ||
node1.vpid.row() == VPChannelID::RowID(0) || node1.vpid.row() == VPChannelID::RowID(255) || node2.vpid.row() == VPChannelID::RowID(0) || node2.vpid.row() == VPChannelID::RowID(255) );
bool Hit_Cut = ( ( track.isVeloBackward() && ( HitCounter>3 || int(module)/2+1==HitCounter ) && int(track.nHits())-HitCounter>3 ) ||
( !(track.isVeloBackward()) && ( HitCounter>3 || 26-int(module)/2==HitCounter ) && int(track.nHits())-HitCounter>3 ) );
if ( Within_One || !Hit_Cut ) ) continue;
}
if ( sensor1 == 0 && sensor2 == 1 ) { // CLI-NLO
++m_overlap_residual_x_CLI_NLO[{module, node1.residualX - node2.residualX}];
++m_overlap_residual_y_CLI_NLO[{module, node1.residualY - node2.residualY}];
++module_CLI_NLO_histo[module];
} else if ( sensor1 == 0 && sensor2 == 2 ) { // CLI-NSI
++m_overlap_residual_x_CLI_NSI[{module, node1.residualX - node2.residualX}];
++m_overlap_residual_y_CLI_NSI[{module, node1.residualY - node2.residualY}];
++module_CLI_NSI_histo[module];
} else if ( sensor1 == 2 && sensor2 == 3 ) { // NSI-CSO
++m_overlap_residual_x_CSO_NSI[{module, node1.residualX - node2.residualX}];
++m_overlap_residual_y_CSO_NSI[{module, node1.residualY - node2.residualY}];
++module_CSO_NSI_histo[module];
} else {
warning() << "How can these sensors overlap?: " << sensor1 << " " << sensor2 << endmsg;
}
if ( m_expertMode.value() ) {
Tuple nodeTuple = nTuple( "TrackVPOverlapMonitor_overlapnodes", "" );
nodeTuple->column( "sensor1", int( sensor1 ) ).ignore();
nodeTuple->column( "sensor2", int( sensor2 ) ).ignore();
nodeTuple->column( "module", int( module ) ).ignore();
nodeTuple->column( "sensor1_row", to_unsigned( node1.vpid.row() ) ).ignore();
nodeTuple->column( "sensor1_scol", node1.vpid.scol() ).ignore();
nodeTuple->column( "sensor2_row", to_unsigned( node2.vpid.row() ) ).ignore();
nodeTuple->column( "sensor2_scol", node2.vpid.scol() ).ignore();
nodeTuple->column( "sensor1_residualX", node1.residualX ).ignore();
nodeTuple->column( "sensor1_residualY", node1.residualY ).ignore();
nodeTuple->column( "sensor2_residualX", node2.residualX ).ignore();
nodeTuple->column( "sensor2_residualY", node2.residualY ).ignore();
nodeTuple->column( "overlap_residualX", node1.residualX - node2.residualX ).ignore();
nodeTuple->column( "overlap_residualY", node1.residualY - node2.residualY ).ignore();
nodeTuple->column( "theta", 2 * std::atan( std::exp( -track.pseudoRapidity() ) ) ).ignore();
nodeTuple->column( "phi", track.phi() ).ignore();
nodeTuple->column( "hits", int( track.nHits() ) ).ignore();
nodeTuple->column( "backward", int( track.isVeloBackward() ) ).ignore();
nodeTuple->column( "x", globalpos.x() ).ignore();
nodeTuple->column( "y", globalpos.y() ).ignore();
nodeTuple->column( "z", globalpos.z() ).ignore();
nodeTuple->column( "HitCounter", HitCounter ).ignore();
nodeTuple->column( "Hit_Cut", Hit_Cut ).ignore();
nodeTuple->column( "Within_One", Within_One ).ignore();
nodeTuple->column( "sensor1_NodeGlobalX", nodeGlobal1.x() ).ignore();
nodeTuple->column( "sensor1_NodeGlobalY", nodeGlobal1.y() ).ignore();
nodeTuple->column( "sensor1_NodeLocalX", nodeLocal1.x() ).ignore();
nodeTuple->column( "sensor1_NodeLocalY", nodeLocal1.y() ).ignore();
nodeTuple->column( "sensor2_NodeGlobalX", nodeGlobal2.x() ).ignore();
nodeTuple->column( "sensor2_NodeGlobalY", nodeGlobal2.y() ).ignore();
nodeTuple->column( "sensor2_NodeLocalX", nodeLocal2.x() ).ignore();
nodeTuple->column( "sensor2_NodeLocalY", nodeLocal2.y() ).ignore();
nodeTuple->write().ignore();
}
}
}
}
Loading