Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

ME0Segment update for muon tdr samples #18127

Merged
merged 7 commits into from Apr 6, 2017
Merged

ME0Segment update for muon tdr samples #18127

Show file tree
Hide file tree
Changes from 5 commits
Commits
Show all changes
7 commits
Select commit Hold shift + click to select a range
3c4c2d8
Moved the time info to float since it is ~integer data and is even ou…
nickmccoll Mar 29, 2017
ef259a5
Store deltaPhi variable in ME0Segment
nickmccoll Mar 29, 2017
d9de5c3
ME0Chamber layer functions had a bug, were checking partition numbers
nickmccoll Mar 29, 2017
fed5cc4
Add parameter to require coincidence with central BX when building se…
nickmccoll Mar 29, 2017
c2365ec
Change parameter defaults to baseline geometry: 8 eta partitions/384 …
nickmccoll Mar 29, 2017
548e308
Moved the deltaPhi function to ME0Chamebr
nickmccoll Mar 30, 2017
d74aaa1
Code cleanup
nickmccoll Apr 3, 2017
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Jump to
Jump to file
Failed to load files.
Diff view
Diff view
2 changes: 2 additions & 0 deletions DataFormats/GEMRecHit/BuildFile.xml
View file
Open in desktop
Expand Up @@ -4,6 +4,8 @@

<use name="DataFormats/GeometryVector"/>
<use name="DataFormats/CSCRecHit"/>
<use name="Geometry/GEMGeometry"/>
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@nickmccoll - DataFormats cannot depend on Geometry

<use name="root"/>
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@nickmccoll - also, brining in root dependency here is not a good idea

<use name="rootrflx"/>

<export>
Expand Down
15 changes: 10 additions & 5 deletions DataFormats/GEMRecHit/interface/ME0Segment.h
View file
Open in desktop
Expand Up @@ -15,21 +15,22 @@

#include <iosfwd>

class ME0Chamber;
class ME0DetId;

class ME0Segment final : public RecSegment {

public:

/// Default constructor
ME0Segment() : theChi2(0.){}
ME0Segment() : theChi2(0.), theTimeValue(0.), theTimeUncrt(0.), theDeltaPhi(0.){}

/// Constructor
ME0Segment(const std::vector<const ME0RecHit*>& proto_segment, const LocalPoint& origin,
const LocalVector& direction, const AlgebraicSymMatrix& errors, double chi2);

ME0Segment(const std::vector<const ME0RecHit*>& proto_segment, const LocalPoint& origin,
const LocalVector& direction, const AlgebraicSymMatrix& errors, double chi2, double time, double timeErr);
const LocalVector& direction, const AlgebraicSymMatrix& errors, double chi2, float time, float timeErr, float deltaPhi);
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Having moved the other class members to single precision, is there any reason for still keeping chi2 in double precision?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The base class, "RecSegment" returns chi2 in double precision. I think it is fine for us to change it to a float here.


/// Destructor
virtual ~ME0Segment();
Expand Down Expand Up @@ -72,7 +73,10 @@ class ME0Segment final : public RecSegment {

float time() const { return theTimeValue; }
float timeErr() const { return theTimeUncrt; }


float deltaPhi() const { return theDeltaPhi; }
static float computeDeltaPhi(const ME0Chamber * chamber, const LocalPoint& position, const LocalVector& direction );
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@nickmccoll - this function does not belong in ME0Segment, nor in DataFormats. I can see it is used in a plugin ME0SegmentAlgorithm. Shouldn't it be moved there? @slava77 - please, suggest a better place for it to avoid circular dependency.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I can put it in the ME0SegmentAlgorithm, if there is no problem that RecoMuonProducer depends on it. Otherwise I can duplicate the code and place it both in the global muon reconstruction and in the segment algorithm.

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@nickmccoll - which class represents RecoMuonProducer? Duplicating the code is not a good idea, nor introducing dependency on a plugin.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@ianna It would be used under: RecoMuon/GlobalMuonProducer/src/GlobalMuonProducer.h and the related classes when building RecoMuons. Technically it should also be used in ME0SegmentMatcher.

I couldn't really think of a good place to put this function, and now that I think about it the ME0Segment isn't a very good place for it at all. I am open to any suggestion!

Thanks!

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@nickmccoll - would RecoMuon/​GlobalTrackingTools be a good place for it?

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I believe the bending angle is something strictly related to the segment, then later to the matching and the muon ID. so the function can go into the algorithms that build the segment. why do you think it is not a good place?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@ianna If there is no problem with ME0SegmentAlgorithm depending on it, then I think it is fine.

@calabria We also need to calculate it for projected pixel tracks as we compare projected and segment bending angles when matching.

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@nickmccoll Usually Reco depends on LocalReco and not vice versa. How about placing it in Geometry/GEMGeometry?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@ianna That is a good idea, it naturally fits in as a ME0Chamber function. I updated the PR, and there are no new dependencies!


void print() const;

private:
Expand All @@ -82,8 +86,9 @@ class ME0Segment final : public RecSegment {
LocalVector theLocalDirection; // in chamber frame - the GeomDet local coordinate system
AlgebraicSymMatrix theCovMatrix; // the covariance matrix
double theChi2; // the Chi squared of the segment fit
double theTimeValue; // the best time estimate of the segment
double theTimeUncrt; // the uncertainty on the time estimation
float theTimeValue; // the best time estimate of the segment
float theTimeUncrt; // the uncertainty on the time estimation
float theDeltaPhi; // Difference in segment phi position: outer layer - inner lay

};

Expand Down
26 changes: 25 additions & 1 deletion DataFormats/GEMRecHit/src/ME0Segment.cc
View file
Open in desktop
Expand Up @@ -5,6 +5,9 @@
*/
#include "FWCore/MessageLogger/interface/MessageLogger.h"
#include "DataFormats/GEMRecHit/interface/ME0Segment.h"
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@nickmccoll - I think, the place for this class is in Geometry/GEMGeometry

#include "Geometry/GEMGeometry/interface/ME0Chamber.h"
#include "Geometry/GEMGeometry/interface/ME0Layer.h"
#include <TVector2.h>
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@nickmccoll - isn't there a CMS equivalent for it in DataFormats/GeometryVector?

#include <iostream>

namespace {
Expand Down Expand Up @@ -38,17 +41,19 @@ ME0Segment::ME0Segment(const std::vector<const ME0RecHit*>& proto_segment, const
theLocalDirection(direction), theCovMatrix(errors), theChi2(chi2){
theTimeValue = 0.0;
theTimeUncrt = 0.0;
theDeltaPhi = 0.0;
for(unsigned int i=0; i<proto_segment.size(); ++i)
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

You can profit of this PR and modify here (and in several other places in this file) in order to use range based for loops

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

also, if you are at it, you can use initialization of a base or a non-static member by constructor initializer list

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Done.

theME0RecHits.push_back(*proto_segment[i]);
}

ME0Segment::ME0Segment(const std::vector<const ME0RecHit*>& proto_segment, const LocalPoint& origin,
const LocalVector& direction, const AlgebraicSymMatrix& errors, double chi2, double time, double timeErr) :
const LocalVector& direction, const AlgebraicSymMatrix& errors, double chi2, float time, float timeErr, float deltaPhi) :
RecSegment(buildDetId(proto_segment.front()->me0Id())),
theOrigin(origin),
theLocalDirection(direction), theCovMatrix(errors), theChi2(chi2){
theTimeValue = time;
theTimeUncrt = timeErr;
theDeltaPhi = deltaPhi;

for(unsigned int i=0; i<proto_segment.size(); ++i)
theME0RecHits.push_back(*proto_segment[i]);
Expand Down Expand Up @@ -124,3 +129,22 @@ std::ostream& operator<<(std::ostream& os, const ME0Segment& seg) {
return os;
}


float ME0Segment::computeDeltaPhi(const ME0Chamber * chamber, const LocalPoint& position, const LocalVector& direction ) {
auto extrap = [] (const LocalPoint& point, const LocalVector& dir, double extZ) -> LocalPoint {
double extX = point.x()+extZ*dir.x()/dir.z();
double extY = point.y()+extZ*dir.y()/dir.z();
return LocalPoint(extX,extY,extZ);
};

const float beginOfChamber = chamber->layer(1)->position().z();
const float centerOfChamber = chamber->position().z();
const float endOfChamber = chamber->layer(chamber->nLayers())->position().z();

LocalPoint projHigh = extrap(position,direction, (centerOfChamber < 0 ? -1.0 : 1.0) * ( endOfChamber- centerOfChamber));
LocalPoint projLow = extrap(position,direction, (centerOfChamber < 0 ? -1.0 : 1.0) *( beginOfChamber- centerOfChamber));
auto globLow = chamber->toGlobal(projLow );
auto globHigh = chamber->toGlobal(projHigh);
return TVector2::Phi_mpi_pi(globHigh.phi() - globLow.phi());
}

3 changes: 2 additions & 1 deletion DataFormats/GEMRecHit/src/classes_def.xml
View file
Open in desktop
Expand Up @@ -45,7 +45,8 @@
<class name="edm::Wrapper<edm::RangeMap<GEMDetId,edm::OwnVector<GEMSegment,edm::ClonePolicy<GEMSegment> >,edm::ClonePolicy<GEMSegment> > >" splitLevel="0"/>
<class name="GEMSegmentRef" splitLevel="0"/>

<class name="ME0Segment" ClassVersion="12">
<class name="ME0Segment" ClassVersion="13">
<version ClassVersion="13" checksum="3580185462"/>
<version ClassVersion="12" checksum="759556370"/>
<version ClassVersion="11" checksum="513208975"/>
<version ClassVersion="10" checksum="2562385753"/>
Expand Down
8 changes: 4 additions & 4 deletions Geometry/GEMGeometry/src/ME0Chamber.cc
View file
Open in desktop
Expand Up @@ -41,13 +41,13 @@ int ME0Chamber::nLayers() const {

const ME0Layer* ME0Chamber::layer(ME0DetId id) const {
if (id.chamberId()!=detId_) return 0; // not in this layer!
return layer(id.roll());
return layer(id.layer());
}

const ME0Layer* ME0Chamber::layer(int isl) const {
for (auto roll : layers_){
if (roll->id().roll()==isl)
return roll;
for (auto layer : layers_){
if (layer->id().layer()==isl)
return layer;
}
return 0;
}
Expand Down
18 changes: 15 additions & 3 deletions RecoLocalMuon/GEMSegment/plugins/ME0SegAlgoRU.cc
View file
Open in desktop
Expand Up @@ -29,14 +29,16 @@ ME0SegAlgoRU::ME0SegAlgoRU(const edm::ParameterSet& ps)

allowWideSegments = ps.getParameter<bool>("allowWideSegments");
doCollisions = ps.getParameter<bool>("doCollisions");

stdParameters.maxChi2Additional = ps.getParameter<double>("maxChi2Additional");
stdParameters.maxChi2Prune = ps.getParameter<double>("maxChi2Prune" );
stdParameters.maxChi2GoodSeg = ps.getParameter<double>("maxChi2GoodSeg" ) ;
stdParameters.maxPhiSeeds = ps.getParameter<double>("maxPhiSeeds" ) ;
stdParameters.maxPhiAdditional = ps.getParameter<double>("maxPhiAdditional" );
stdParameters.maxETASeeds = ps.getParameter<double>("maxETASeeds" );
stdParameters.maxTOFDiff = ps.getParameter<double>("maxTOFDiff" );
stdParameters.minNumberOfHits = ps.getParameter<unsigned int>("minNumberOfHits" );
stdParameters.requireCentralBX = ps.getParameter<bool>("requireCentralBX" );
stdParameters.minNumberOfHits = ps.getParameter<unsigned int>("minNumberOfHits" );
stdParameters.requireBeamConstr = true;
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

ParameterSetDescription should be used for these parameters
Are you at ease at having the wideParameters and displacedParameters (i.e. their differences wrt stdParameters) hard-coded in the c++ class code? Shouldn't you better also define them via configurable parameters?

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

See #17687 (review) for a previous discussion about that

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thank you @kpedro88 for remembering me : indeed in that (already merged) PR we argued that in a next update of the code the hardcoded parameters should have been made configurable, and everything would have possibly moved to ParameterSetDescription: aren't we arrived at that "next update" yet?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@perrotta I was thinking that the next update in the algorithm that can actually use the displaced and wide parameters. As of now, they are not used in the code at all as there have been no studies to understand if they are needed/or what to do. If you would like I can make them parameters, but I feel more odd about having parameters that don't change anything. I also feel a little odd about having the unused code there, but elected to keep it in for easier future development.

I can do the ParameterSetDescription

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@perrotta I addressed the other issues brought up in my last commit but did not add the "ParameterSetDescription." I understand that this will be a rather large change as the class uses a vector of parameter sets with different parameter for different algorithms. Is there a suggested way to handle this?

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thank you @nickmccoll . One can define an unitialized ParameterSet in the base class, which is then explicitely initialized in the different algos, see https://twiki.cern.ch/twiki/bin/view/CMSPublic/SWGuideConfigurationValidationAndHelp
If you can, it wouldn't be bad to fix it since now. But if it is a big issue, you can factorize your job and plan it for a forthcoming (hopefully next) pull request

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@perrotta thanks for the information. I've looked at this page and don't see an application that will work for me. For example, I can have each algorithm have its own "ConfigurationDescriptions" that then is called on by the Producer's "fillDescription" function. But I don't see how it will resolve the problem of having a vector of ParameterSets with different fields.

I think that what we would need to do is change the way the ParameterSet is structured, but would prefer to not to do that here. I don't know how close we are to the next pre-release but we would like for this to get in here. Also, this exact same structure is used in GEMSegmentProducer and CSCSegmentProducer....I would hope that we would have a solution that matches across all three classes.

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@nickmccoll Fine with me: don't touch the ConfigurationDescription by now


wideParameters = stdParameters;
Expand Down Expand Up @@ -214,9 +216,18 @@ void ME0SegAlgoRU::lookForSegments( const SegmentParameters& params, const unsig
edm::LogVerbatim("ME0SegAlgoRU") << "[ME0SegAlgoRU::lookForSegments] # of hits in segment " << current_proto_segment.size() <<" min # "<< n_seg_min <<" => "<< (current_proto_segment.size() >= n_seg_min) << " chi2/ndof "<<current_fit->chi2()/current_fit->ndof()<<" => "<<(current_fit->chi2()/current_fit->ndof() < params.maxChi2GoodSeg ) <<std::endl;

if(current_proto_segment.size() < n_seg_min) continue;

const float current_metric = current_fit->chi2()/current_fit->ndof();
if(current_metric > params.maxChi2GoodSeg) continue;

if(params.requireCentralBX){
int nCentral = 0;
int nNonCentral =0;
for(const auto* rh : current_proto_segment ) {
if(std::fabs(rh->rh->tof()) < 2 ) nCentral++;
else nNonCentral++;
}
if(nNonCentral >= nCentral) continue;
}
// segok = true;

proto_segments.emplace_back( current_metric, current_proto_segment);
Expand Down Expand Up @@ -274,8 +285,9 @@ void ME0SegAlgoRU::addUniqueSegments(SegmentByMetricContainer& proto_segments, s

std::vector<const ME0RecHit*> bareRHs; bareRHs.reserve(currentProtoSegment.size());
for(const auto* rh : currentProtoSegment) bareRHs.push_back(rh->rh);
const float dPhi = ME0Segment::computeDeltaPhi(theChamber,current_fit->intercept(),current_fit->localdir());
ME0Segment temp(bareRHs, current_fit->intercept(),
current_fit->localdir(), current_fit->covarianceMatrix(), current_fit->chi2(), averageTime, timeUncrt);
current_fit->localdir(), current_fit->covarianceMatrix(), current_fit->chi2(), averageTime, timeUncrt, dPhi);
segments.push_back(temp);


Expand Down
1 change: 1 addition & 0 deletions RecoLocalMuon/GEMSegment/plugins/ME0SegAlgoRU.h
View file
Open in desktop
Expand Up @@ -40,6 +40,7 @@ class ME0SegAlgoRU : public ME0SegmentAlgorithmBase {
float maxChi2Prune ;
float maxChi2GoodSeg ;
float maxTOFDiff ;
bool requireCentralBX ;
unsigned int minNumberOfHits;
bool requireBeamConstr ;

Expand Down
10 changes: 6 additions & 4 deletions RecoLocalMuon/GEMSegment/plugins/ME0SegmentAlgorithm.cc
View file
Open in desktop
Expand Up @@ -81,7 +81,7 @@ std::vector<ME0Segment> ME0SegmentAlgorithm::run(const ME0Chamber * chamber, con
// clear the buffer for the subset of segments:
segments_temp.clear();
// build the subset of segments:
this->buildSegments((*sub_rechits), segments_temp );
this->buildSegments(chamber,(*sub_rechits), segments_temp );
// add the found subset of segments to the collection of all segments in this chamber:
segments.insert( segments.end(), segments_temp.begin(), segments_temp.end() );
}
Expand All @@ -90,7 +90,7 @@ std::vector<ME0Segment> ME0SegmentAlgorithm::run(const ME0Chamber * chamber, con
else {
HitAndPositionPtrContainer ptrRH; ptrRH.reserve(rechits.size());
for(const auto& rh : rechits) ptrRH.push_back(&rh);
this->buildSegments(ptrRH, segments);
this->buildSegments(chamber,ptrRH, segments);
return segments;
}
}
Expand Down Expand Up @@ -288,7 +288,7 @@ bool ME0SegmentAlgorithm::isGoodToMerge(const ME0Chamber * chamber, const HitAnd
return false;
}

void ME0SegmentAlgorithm::buildSegments(const HitAndPositionPtrContainer& rechits, std::vector<ME0Segment>& me0segs) {
void ME0SegmentAlgorithm::buildSegments(const ME0Chamber * chamber, const HitAndPositionPtrContainer& rechits, std::vector<ME0Segment>& me0segs) {
if (rechits.size() < minHitsPerSegment) return;

#ifdef EDM_ML_DEBUG // have lines below only compiled when in debug mode
Expand Down Expand Up @@ -344,9 +344,11 @@ void ME0SegmentAlgorithm::buildSegments(const HitAndPositionPtrContainer& rechit
if(rechits.size() > 1) timeUncrt=timeUncrt/(rechits.size()-1);
timeUncrt = sqrt(timeUncrt);

const float dPhi = ME0Segment::computeDeltaPhi(chamber,protoIntercept,protoDirection);

// save all information inside GEMCSCSegment
edm::LogVerbatim("ME0SegmentAlgorithm") << "[ME0SegmentAlgorithm::buildSegments] will now try to make ME0Segment from collection of "<<rechits.size()<<" ME0 RecHits";
ME0Segment tmp(bareRHs, protoIntercept, protoDirection, protoErrors, protoChi2, averageTime, timeUncrt);
ME0Segment tmp(bareRHs, protoIntercept, protoDirection, protoErrors, protoChi2, averageTime, timeUncrt, dPhi);

edm::LogVerbatim("ME0SegmentAlgorithm") << "[ME0SegmentAlgorithm::buildSegments] ME0Segment made";
edm::LogVerbatim("ME0SegmentAlgorithm") << "[ME0SegmentAlgorithm::buildSegments] "<<tmp;
Expand Down
2 changes: 1 addition & 1 deletion RecoLocalMuon/GEMSegment/plugins/ME0SegmentAlgorithm.h
View file
Open in desktop
Expand Up @@ -51,7 +51,7 @@ class ME0SegmentAlgorithm : public ME0SegmentAlgorithmBase {
bool isGoodToMerge(const ME0Chamber * chamber, const HitAndPositionPtrContainer& newChain, const HitAndPositionPtrContainer& oldChain);

// Build track segments in this chamber (this is where the actual segment-building algorithm hides.)
void buildSegments(const HitAndPositionPtrContainer& rechits, std::vector<ME0Segment>& me0segs);
void buildSegments(const ME0Chamber * chamber, const HitAndPositionPtrContainer& rechits, std::vector<ME0Segment>& me0segs);

// Member variables
const std::string myName;
Expand Down
7 changes: 4 additions & 3 deletions RecoLocalMuon/GEMSegment/python/ME0SegmentsRU_cfi.py
View file
Open in desktop
Expand Up @@ -6,10 +6,11 @@
maxChi2Additional = cms.double(100.0),
maxChi2Prune = cms.double(50),
maxChi2GoodSeg = cms.double(50),
maxPhiSeeds = cms.double(0.000547), #Assuming 768 strips
maxPhiAdditional = cms.double(0.000547), #Assuming 768 strips
maxETASeeds = cms.double(0.08), #Assuming 8 eta partitions
maxPhiSeeds = cms.double(0.001096605744), #Assuming 384 strips
maxPhiAdditional = cms.double(0.001096605744), #Assuming 384 strips
maxETASeeds = cms.double(0.1), #Assuming 8 eta partitions
maxTOFDiff = cms.double(25),
requireCentralBX = cms.bool(True), #require that a majority of hits come from central BX
minNumberOfHits = cms.uint32(4),
)

Expand Down