File indexing completed on 2024-04-06 11:56:43
0001 #include "Alignment/MuonAlignmentAlgorithms/plugins/CSCOverlapsAlignmentAlgorithm.h"
0002
0003 CSCOverlapsAlignmentAlgorithm::CSCOverlapsAlignmentAlgorithm(const edm::ParameterSet& iConfig,
0004 edm::ConsumesCollector& iC)
0005 : AlignmentAlgorithmBase(iConfig, iC),
0006 m_minHitsPerChamber(iConfig.getParameter<int>("minHitsPerChamber")),
0007 m_maxdrdz(iConfig.getParameter<double>("maxdrdz")),
0008 m_fiducial(iConfig.getParameter<bool>("fiducial")),
0009 m_useHitWeights(iConfig.getParameter<bool>("useHitWeights")),
0010 m_slopeFromTrackRefit(iConfig.getParameter<bool>("slopeFromTrackRefit")),
0011 m_minStationsInTrackRefits(iConfig.getParameter<int>("minStationsInTrackRefits")),
0012 m_truncateSlopeResid(iConfig.getParameter<double>("truncateSlopeResid")),
0013 m_truncateOffsetResid(iConfig.getParameter<double>("truncateOffsetResid")),
0014 m_combineME11(iConfig.getParameter<bool>("combineME11")),
0015 m_useTrackWeights(iConfig.getParameter<bool>("useTrackWeights")),
0016 m_errorFromRMS(iConfig.getParameter<bool>("errorFromRMS")),
0017 m_minTracksPerOverlap(iConfig.getParameter<int>("minTracksPerOverlap")),
0018 m_makeHistograms(iConfig.getParameter<bool>("makeHistograms")),
0019 m_cscGeometryToken(iC.esConsumes<edm::Transition::BeginRun>()),
0020 m_propToken(iC.esConsumes(edm::ESInputTag(
0021 "",
0022 m_slopeFromTrackRefit
0023 ? iConfig.getParameter<edm::ParameterSet>("TrackTransformer").getParameter<std::string>("Propagator")
0024 : std::string("")))),
0025 m_tthbToken(iC.esConsumes(edm::ESInputTag("", "TransientTrackBuilder"))),
0026 m_mode_string(iConfig.getParameter<std::string>("mode")),
0027 m_reportFileName(iConfig.getParameter<std::string>("reportFileName")),
0028 m_minP(iConfig.getParameter<double>("minP")),
0029 m_maxRedChi2(iConfig.getParameter<double>("maxRedChi2")),
0030 m_writeTemporaryFile(iConfig.getParameter<std::string>("writeTemporaryFile")),
0031 m_readTemporaryFiles(iConfig.getParameter<std::vector<std::string> >("readTemporaryFiles")),
0032 m_doAlignment(iConfig.getParameter<bool>("doAlignment")) {
0033 if (m_mode_string == std::string("phiy"))
0034 m_mode = CSCPairResidualsConstraint::kModePhiy;
0035 else if (m_mode_string == std::string("phipos"))
0036 m_mode = CSCPairResidualsConstraint::kModePhiPos;
0037 else if (m_mode_string == std::string("phiz"))
0038 m_mode = CSCPairResidualsConstraint::kModePhiz;
0039 else if (m_mode_string == std::string("radius"))
0040 m_mode = CSCPairResidualsConstraint::kModeRadius;
0041 else
0042 throw cms::Exception("BadConfig") << "mode must be one of \"phiy\", \"phipos\", \"phiz\", \"radius\"" << std::endl;
0043
0044 std::vector<edm::ParameterSet> fitters = iConfig.getParameter<std::vector<edm::ParameterSet> >("fitters");
0045 for (std::vector<edm::ParameterSet>::const_iterator fitter = fitters.begin(); fitter != fitters.end(); ++fitter) {
0046 m_fitters.push_back(CSCChamberFitter(*fitter, m_residualsConstraints));
0047 }
0048
0049 for (std::vector<CSCPairResidualsConstraint*>::const_iterator residualsConstraint = m_residualsConstraints.begin();
0050 residualsConstraint != m_residualsConstraints.end();
0051 ++residualsConstraint) {
0052 (*residualsConstraint)->configure(this);
0053 m_quickChamberLookup[std::pair<CSCDetId, CSCDetId>((*residualsConstraint)->id_i(), (*residualsConstraint)->id_j())] =
0054 *residualsConstraint;
0055 }
0056
0057 if (m_slopeFromTrackRefit) {
0058 m_trackTransformer = new TrackTransformer(iConfig.getParameter<edm::ParameterSet>("TrackTransformer"), iC);
0059 } else {
0060 m_trackTransformer = nullptr;
0061 }
0062
0063 m_propagatorPointer = nullptr;
0064
0065 if (m_makeHistograms) {
0066 edm::Service<TFileService> tFileService;
0067 m_histP10 = tFileService->make<TH1F>("P10", "", 100, 0, 10);
0068 m_histP100 = tFileService->make<TH1F>("P100", "", 100, 0, 100);
0069 m_histP1000 = tFileService->make<TH1F>("P1000", "", 100, 0, 1000);
0070
0071 m_hitsPerChamber = tFileService->make<TH1F>("hitsPerChamber", "", 10, -0.5, 9.5);
0072
0073 m_fiducial_ME11 = tFileService->make<TProfile>("fiducial_ME11", "", 100, 0.075, 0.100);
0074 m_fiducial_ME12 = tFileService->make<TProfile>("fiducial_ME12", "", 100, 0.080, 0.105);
0075 m_fiducial_MEx1 = tFileService->make<TProfile>("fiducial_MEx1", "", 100, 0.160, 0.210);
0076 m_fiducial_MEx2 = tFileService->make<TProfile>("fiducial_MEx2", "", 100, 0.080, 0.105);
0077
0078 m_slope = tFileService->make<TH1F>("slope", "", 100, -0.5, 0.5);
0079 m_slope_MEp4 = tFileService->make<TH1F>("slope_MEp4", "", 100, -0.5, 0.5);
0080 m_slope_MEp3 = tFileService->make<TH1F>("slope_MEp3", "", 100, -0.5, 0.5);
0081 m_slope_MEp2 = tFileService->make<TH1F>("slope_MEp2", "", 100, -0.5, 0.5);
0082 m_slope_MEp1 = tFileService->make<TH1F>("slope_MEp1", "", 100, -0.5, 0.5);
0083 m_slope_MEm1 = tFileService->make<TH1F>("slope_MEm1", "", 100, -0.5, 0.5);
0084 m_slope_MEm2 = tFileService->make<TH1F>("slope_MEm2", "", 100, -0.5, 0.5);
0085 m_slope_MEm3 = tFileService->make<TH1F>("slope_MEm3", "", 100, -0.5, 0.5);
0086 m_slope_MEm4 = tFileService->make<TH1F>("slope_MEm4", "", 100, -0.5, 0.5);
0087
0088 m_slopeResiduals = tFileService->make<TH1F>("slopeResiduals", "mrad", 300, -30., 30.);
0089 m_slopeResiduals_weighted = tFileService->make<TH1F>("slopeResiduals_weighted", "mrad", 300, -30., 30.);
0090 m_slopeResiduals_normalized = tFileService->make<TH1F>("slopeResiduals_normalized", "", 200, -20., 20.);
0091 m_offsetResiduals = tFileService->make<TH1F>("offsetResiduals", "mm", 300, -30., 30.);
0092 m_offsetResiduals_weighted = tFileService->make<TH1F>("offsetResiduals_weighted", "mm", 300, -30., 30.);
0093 m_offsetResiduals_normalized = tFileService->make<TH1F>("offsetResiduals_normalized", "", 200, -20., 20.);
0094
0095 m_drdz = tFileService->make<TH1F>("drdz", "", 100, -0.5, 0.5);
0096
0097 m_occupancy = tFileService->make<TH2F>("occupancy", "", 36, 1, 37, 20, 1, 21);
0098 for (int i = 1; i <= 36; i++) {
0099 std::stringstream pairname;
0100 pairname << i << "-";
0101 if (i + 1 == 37)
0102 pairname << 1;
0103 else
0104 pairname << (i + 1);
0105 m_occupancy->GetXaxis()->SetBinLabel(i, pairname.str().c_str());
0106 }
0107 m_occupancy->GetYaxis()->SetBinLabel(1, "ME-4/2");
0108 m_occupancy->GetYaxis()->SetBinLabel(2, "ME-4/1");
0109 m_occupancy->GetYaxis()->SetBinLabel(3, "ME-3/2");
0110 m_occupancy->GetYaxis()->SetBinLabel(4, "ME-3/1");
0111 m_occupancy->GetYaxis()->SetBinLabel(5, "ME-2/2");
0112 m_occupancy->GetYaxis()->SetBinLabel(6, "ME-2/1");
0113 m_occupancy->GetYaxis()->SetBinLabel(7, "ME-1/3");
0114 m_occupancy->GetYaxis()->SetBinLabel(8, "ME-1/2");
0115 if (!m_combineME11) {
0116 m_occupancy->GetYaxis()->SetBinLabel(9, "ME-1/1b");
0117 m_occupancy->GetYaxis()->SetBinLabel(10, "ME-1/1a");
0118 m_occupancy->GetYaxis()->SetBinLabel(11, "ME+1/1a");
0119 m_occupancy->GetYaxis()->SetBinLabel(12, "ME+1/1b");
0120 } else {
0121 m_occupancy->GetYaxis()->SetBinLabel(9, "ME-1/1");
0122 m_occupancy->GetYaxis()->SetBinLabel(10, "");
0123 m_occupancy->GetYaxis()->SetBinLabel(11, "");
0124 m_occupancy->GetYaxis()->SetBinLabel(12, "ME+1/1");
0125 }
0126 m_occupancy->GetYaxis()->SetBinLabel(13, "ME+1/2");
0127 m_occupancy->GetYaxis()->SetBinLabel(14, "ME+1/3");
0128 m_occupancy->GetYaxis()->SetBinLabel(15, "ME+2/1");
0129 m_occupancy->GetYaxis()->SetBinLabel(16, "ME+2/2");
0130 m_occupancy->GetYaxis()->SetBinLabel(17, "ME+3/1");
0131 m_occupancy->GetYaxis()->SetBinLabel(18, "ME+3/2");
0132 m_occupancy->GetYaxis()->SetBinLabel(19, "ME+4/1");
0133 m_occupancy->GetYaxis()->SetBinLabel(20, "ME+4/2");
0134
0135 m_XYpos_mep1 = tFileService->make<TH2F>("XYpos_mep1", "Positions: ME+1", 140, -700., 700., 140, -700., 700.);
0136 m_XYpos_mep2 = tFileService->make<TH2F>("XYpos_mep2", "Positions: ME+2", 140, -700., 700., 140, -700., 700.);
0137 m_XYpos_mep3 = tFileService->make<TH2F>("XYpos_mep3", "Positions: ME+3", 140, -700., 700., 140, -700., 700.);
0138 m_XYpos_mep4 = tFileService->make<TH2F>("XYpos_mep4", "Positions: ME+4", 140, -700., 700., 140, -700., 700.);
0139 m_XYpos_mem1 = tFileService->make<TH2F>("XYpos_mem1", "Positions: ME-1", 140, -700., 700., 140, -700., 700.);
0140 m_XYpos_mem2 = tFileService->make<TH2F>("XYpos_mem2", "Positions: ME-2", 140, -700., 700., 140, -700., 700.);
0141 m_XYpos_mem3 = tFileService->make<TH2F>("XYpos_mem3", "Positions: ME-3", 140, -700., 700., 140, -700., 700.);
0142 m_XYpos_mem4 = tFileService->make<TH2F>("XYpos_mem4", "Positions: ME-4", 140, -700., 700., 140, -700., 700.);
0143 m_RPhipos_mep1 = tFileService->make<TH2F>("RPhipos_mep1", "Positions: ME+1", 144, -M_PI, M_PI, 21, 0., 700.);
0144 m_RPhipos_mep2 = tFileService->make<TH2F>("RPhipos_mep2", "Positions: ME+2", 144, -M_PI, M_PI, 21, 0., 700.);
0145 m_RPhipos_mep3 = tFileService->make<TH2F>("RPhipos_mep3", "Positions: ME+3", 144, -M_PI, M_PI, 21, 0., 700.);
0146 m_RPhipos_mep4 = tFileService->make<TH2F>("RPhipos_mep4", "Positions: ME+4", 144, -M_PI, M_PI, 21, 0., 700.);
0147 m_RPhipos_mem1 = tFileService->make<TH2F>("RPhipos_mem1", "Positions: ME-1", 144, -M_PI, M_PI, 21, 0., 700.);
0148 m_RPhipos_mem2 = tFileService->make<TH2F>("RPhipos_mem2", "Positions: ME-2", 144, -M_PI, M_PI, 21, 0., 700.);
0149 m_RPhipos_mem3 = tFileService->make<TH2F>("RPhipos_mem3", "Positions: ME-3", 144, -M_PI, M_PI, 21, 0., 700.);
0150 m_RPhipos_mem4 = tFileService->make<TH2F>("RPhipos_mem4", "Positions: ME-4", 144, -M_PI, M_PI, 21, 0., 700.);
0151 } else {
0152 m_histP10 = nullptr;
0153 m_histP100 = nullptr;
0154 m_histP1000 = nullptr;
0155 m_hitsPerChamber = nullptr;
0156 m_fiducial_ME11 = nullptr;
0157 m_fiducial_ME12 = nullptr;
0158 m_fiducial_MEx1 = nullptr;
0159 m_fiducial_MEx2 = nullptr;
0160 m_slope = nullptr;
0161 m_slope_MEp4 = nullptr;
0162 m_slope_MEp3 = nullptr;
0163 m_slope_MEp2 = nullptr;
0164 m_slope_MEp1 = nullptr;
0165 m_slope_MEm1 = nullptr;
0166 m_slope_MEm2 = nullptr;
0167 m_slope_MEm3 = nullptr;
0168 m_slope_MEm4 = nullptr;
0169 m_slopeResiduals = nullptr;
0170 m_slopeResiduals_weighted = nullptr;
0171 m_slopeResiduals_normalized = nullptr;
0172 m_offsetResiduals = nullptr;
0173 m_offsetResiduals_weighted = nullptr;
0174 m_offsetResiduals_normalized = nullptr;
0175 m_drdz = nullptr;
0176 m_occupancy = nullptr;
0177 m_XYpos_mep1 = nullptr;
0178 m_XYpos_mep2 = nullptr;
0179 m_XYpos_mep3 = nullptr;
0180 m_XYpos_mep4 = nullptr;
0181 m_XYpos_mem1 = nullptr;
0182 m_XYpos_mem2 = nullptr;
0183 m_XYpos_mem3 = nullptr;
0184 m_XYpos_mem4 = nullptr;
0185 m_RPhipos_mep1 = nullptr;
0186 m_RPhipos_mep2 = nullptr;
0187 m_RPhipos_mep3 = nullptr;
0188 m_RPhipos_mep4 = nullptr;
0189 m_RPhipos_mem1 = nullptr;
0190 m_RPhipos_mem2 = nullptr;
0191 m_RPhipos_mem3 = nullptr;
0192 m_RPhipos_mem4 = nullptr;
0193 }
0194 }
0195
0196 CSCOverlapsAlignmentAlgorithm::~CSCOverlapsAlignmentAlgorithm() {}
0197
0198 void CSCOverlapsAlignmentAlgorithm::initialize(const edm::EventSetup& iSetup,
0199 AlignableTracker* alignableTracker,
0200 AlignableMuon* alignableMuon,
0201 AlignableExtras* alignableExtras,
0202 AlignmentParameterStore* alignmentParameterStore) {
0203 m_alignmentParameterStore = alignmentParameterStore;
0204 m_alignables = m_alignmentParameterStore->alignables();
0205
0206 if (alignableTracker == nullptr)
0207 m_alignableNavigator = new AlignableNavigator(alignableMuon);
0208 else
0209 m_alignableNavigator = new AlignableNavigator(alignableTracker, alignableMuon);
0210
0211 for (const auto& alignable : m_alignables) {
0212 DetId id = alignable->geomDetId();
0213 if (id.det() != DetId::Muon || id.subdetId() != MuonSubdetId::CSC || CSCDetId(id.rawId()).layer() != 0) {
0214 throw cms::Exception("BadConfig") << "Only CSC chambers may be alignable" << std::endl;
0215 }
0216
0217 std::vector<bool> selector = alignable->alignmentParameters()->selector();
0218 for (std::vector<bool>::const_iterator i = selector.begin(); i != selector.end(); ++i) {
0219 if (!(*i))
0220 throw cms::Exception("BadConfig") << "All selector strings should be \"111111\"" << std::endl;
0221 }
0222 }
0223
0224 const CSCGeometry* cscGeometry = &iSetup.getData(m_cscGeometryToken);
0225
0226 for (std::vector<CSCPairResidualsConstraint*>::const_iterator residualsConstraint = m_residualsConstraints.begin();
0227 residualsConstraint != m_residualsConstraints.end();
0228 ++residualsConstraint) {
0229 (*residualsConstraint)->setZplane(cscGeometry);
0230 }
0231
0232 if (!m_readTemporaryFiles.empty()) {
0233 std::vector<std::ifstream*> input;
0234 for (std::vector<std::string>::const_iterator fileName = m_readTemporaryFiles.begin();
0235 fileName != m_readTemporaryFiles.end();
0236 ++fileName) {
0237 input.push_back(new std::ifstream(fileName->c_str()));
0238 }
0239
0240 for (std::vector<CSCPairResidualsConstraint*>::const_iterator residualsConstraint = m_residualsConstraints.begin();
0241 residualsConstraint != m_residualsConstraints.end();
0242 ++residualsConstraint) {
0243 (*residualsConstraint)->read(input, m_readTemporaryFiles);
0244 }
0245
0246 for (std::vector<std::ifstream*>::const_iterator file = input.begin(); file != input.end(); ++file) {
0247 delete (*file);
0248 }
0249 }
0250 }
0251
0252 void CSCOverlapsAlignmentAlgorithm::run(const edm::EventSetup& iSetup, const EventInfo& eventInfo) {
0253 edm::ESHandle<Propagator> propagator;
0254 if (m_slopeFromTrackRefit) {
0255 iSetup.getHandle(m_propToken);
0256 if (m_propagatorPointer != &*propagator) {
0257 m_propagatorPointer = &*propagator;
0258
0259 for (std::vector<CSCPairResidualsConstraint*>::const_iterator residualsConstraint =
0260 m_residualsConstraints.begin();
0261 residualsConstraint != m_residualsConstraints.end();
0262 ++residualsConstraint) {
0263 (*residualsConstraint)->setPropagator(m_propagatorPointer);
0264 }
0265 }
0266 }
0267
0268 const TransientTrackBuilder* transientTrackBuilder = &iSetup.getData(m_tthbToken);
0269
0270 if (m_trackTransformer != nullptr)
0271 m_trackTransformer->setServices(iSetup);
0272
0273 const ConstTrajTrackPairCollection& trajtracks = eventInfo.trajTrackPairs();
0274 for (ConstTrajTrackPairCollection::const_iterator trajtrack = trajtracks.begin(); trajtrack != trajtracks.end();
0275 ++trajtrack) {
0276 const Trajectory* traj = (*trajtrack).first;
0277 const reco::Track* track = (*trajtrack).second;
0278
0279 if (m_makeHistograms) {
0280 m_histP10->Fill(track->p());
0281 m_histP100->Fill(track->p());
0282 m_histP1000->Fill(track->p());
0283 }
0284 if (track->p() >= m_minP) {
0285 std::vector<TrajectoryMeasurement> measurements = traj->measurements();
0286 reco::TransientTrack transientTrack = transientTrackBuilder->build(track);
0287
0288 std::map<int, std::map<CSCDetId, bool> > stationsToChambers;
0289 for (std::vector<TrajectoryMeasurement>::const_iterator measurement = measurements.begin();
0290 measurement != measurements.end();
0291 ++measurement) {
0292 DetId id = measurement->recHit()->geographicalId();
0293 if (id.det() == DetId::Muon && id.subdetId() == MuonSubdetId::CSC) {
0294 CSCDetId cscid(id.rawId());
0295 CSCDetId chamberId(cscid.endcap(), cscid.station(), cscid.ring(), cscid.chamber(), 0);
0296 if (m_combineME11 && cscid.station() == 1 && cscid.ring() == 4)
0297 chamberId = CSCDetId(cscid.endcap(), 1, 1, cscid.chamber(), 0);
0298 int station = (cscid.endcap() == 1 ? 1 : -1) * cscid.station();
0299
0300 if (stationsToChambers.find(station) == stationsToChambers.end())
0301 stationsToChambers[station] = std::map<CSCDetId, bool>();
0302 stationsToChambers[station][chamberId] = true;
0303
0304 if (m_makeHistograms) {
0305 GlobalPoint pos = measurement->recHit()->globalPosition();
0306 if (cscid.endcap() == 1 && cscid.station() == 1) {
0307 m_XYpos_mep1->Fill(pos.x(), pos.y());
0308 m_RPhipos_mep1->Fill(pos.phi(), pos.perp());
0309 }
0310 if (cscid.endcap() == 1 && cscid.station() == 2) {
0311 m_XYpos_mep2->Fill(pos.x(), pos.y());
0312 m_RPhipos_mep2->Fill(pos.phi(), pos.perp());
0313 }
0314 if (cscid.endcap() == 1 && cscid.station() == 3) {
0315 m_XYpos_mep3->Fill(pos.x(), pos.y());
0316 m_RPhipos_mep3->Fill(pos.phi(), pos.perp());
0317 }
0318 if (cscid.endcap() == 1 && cscid.station() == 4) {
0319 m_XYpos_mep4->Fill(pos.x(), pos.y());
0320 m_RPhipos_mep4->Fill(pos.phi(), pos.perp());
0321 }
0322 if (cscid.endcap() == 2 && cscid.station() == 1) {
0323 m_XYpos_mem1->Fill(pos.x(), pos.y());
0324 m_RPhipos_mem1->Fill(pos.phi(), pos.perp());
0325 }
0326 if (cscid.endcap() == 2 && cscid.station() == 2) {
0327 m_XYpos_mem2->Fill(pos.x(), pos.y());
0328 m_RPhipos_mem2->Fill(pos.phi(), pos.perp());
0329 }
0330 if (cscid.endcap() == 2 && cscid.station() == 3) {
0331 m_XYpos_mem3->Fill(pos.x(), pos.y());
0332 m_RPhipos_mem3->Fill(pos.phi(), pos.perp());
0333 }
0334 if (cscid.endcap() == 2 && cscid.station() == 4) {
0335 m_XYpos_mem4->Fill(pos.x(), pos.y());
0336 m_RPhipos_mem4->Fill(pos.phi(), pos.perp());
0337 }
0338 }
0339 }
0340 }
0341
0342 std::map<CSCPairResidualsConstraint*, bool> residualsConstraints;
0343 for (std::map<int, std::map<CSCDetId, bool> >::const_iterator iter = stationsToChambers.begin();
0344 iter != stationsToChambers.end();
0345 ++iter) {
0346 for (std::map<CSCDetId, bool>::const_iterator one = iter->second.begin(); one != iter->second.end(); ++one) {
0347 for (std::map<CSCDetId, bool>::const_iterator two = one; two != iter->second.end(); ++two) {
0348 if (one != two) {
0349 std::map<std::pair<CSCDetId, CSCDetId>, CSCPairResidualsConstraint*>::const_iterator quick;
0350
0351 quick = m_quickChamberLookup.find(std::pair<CSCDetId, CSCDetId>(one->first, two->first));
0352 if (quick != m_quickChamberLookup.end())
0353 residualsConstraints[quick->second] = true;
0354
0355 quick = m_quickChamberLookup.find(std::pair<CSCDetId, CSCDetId>(two->first, one->first));
0356 if (quick != m_quickChamberLookup.end())
0357 residualsConstraints[quick->second] = true;
0358 }
0359 }
0360 }
0361 }
0362
0363 for (std::map<CSCPairResidualsConstraint*, bool>::const_iterator residualsConstraint =
0364 residualsConstraints.begin();
0365 residualsConstraint != residualsConstraints.end();
0366 ++residualsConstraint) {
0367 residualsConstraint->first->addTrack(measurements, transientTrack, m_trackTransformer);
0368 }
0369 }
0370 }
0371 }
0372
0373 void CSCOverlapsAlignmentAlgorithm::terminate(const edm::EventSetup& iSetup) {
0374
0375 if (m_writeTemporaryFile != std::string("")) {
0376 std::ofstream output(m_writeTemporaryFile.c_str());
0377 for (std::vector<CSCPairResidualsConstraint*>::const_iterator residualsConstraint = m_residualsConstraints.begin();
0378 residualsConstraint != m_residualsConstraints.end();
0379 ++residualsConstraint) {
0380 (*residualsConstraint)->write(output);
0381 }
0382 }
0383
0384
0385 if (m_doAlignment) {
0386 std::ofstream report;
0387 bool writeReport = (m_reportFileName != std::string(""));
0388 if (writeReport) {
0389 report.open(m_reportFileName.c_str());
0390 report << "cscReports = []" << std::endl
0391 << std::endl
0392 << "class CSCChamberCorrection:" << std::endl
0393 << " def __init__(self, name, detid, value):" << std::endl
0394 << " self.name, self.detid, self.value = name, detid, value" << std::endl
0395 << std::endl
0396 << "class CSCErrorMode:" << std::endl
0397 << " def __init__(self, error):" << std::endl
0398 << " self.error = error" << std::endl
0399 << " self.terms = {}" << std::endl
0400 << " self.detids = {}" << std::endl
0401 << " def addTerm(self, name, detid, coefficient):" << std::endl
0402 << " self.terms[name] = coefficient" << std::endl
0403 << " self.detids[name] = detid" << std::endl
0404 << std::endl
0405 << "class CSCConstraintResidual:" << std::endl
0406 << " def __init__(self, i, j, before, uncert, residual, pull):" << std::endl
0407 << " self.i, self.j, self.before, self.error, self.residual, self.pull = i, j, before, uncert, "
0408 "residual, pull"
0409 << std::endl
0410 << std::endl
0411 << "class CSCFitterReport:" << std::endl
0412 << " def __init__(self, name, oldchi2, newchi2):" << std::endl
0413 << " self.name, self.oldchi2, self.newchi2 = name, oldchi2, newchi2" << std::endl
0414 << " self.chamberCorrections = []" << std::endl
0415 << " self.errorModes = []" << std::endl
0416 << " self.constraintResiduals = []" << std::endl
0417 << std::endl
0418 << " def addChamberCorrection(self, name, detid, value):" << std::endl
0419 << " self.chamberCorrections.append(CSCChamberCorrection(name, detid, value))" << std::endl
0420 << std::endl
0421 << " def addErrorMode(self, error):" << std::endl
0422 << " self.errorModes.append(CSCErrorMode(error))" << std::endl
0423 << std::endl
0424 << " def addErrorModeTerm(self, name, detid, coefficient):" << std::endl
0425 << " self.errorModes[-1].addTerm(name, detid, coefficient)" << std::endl
0426 << std::endl
0427 << " def addCSCConstraintResidual(self, i, j, before, uncert, residual, pull):" << std::endl
0428 << " self.constraintResiduals.append(CSCConstraintResidual(i, j, before, uncert, residual, pull))"
0429 << std::endl
0430 << std::endl
0431 << "import re" << std::endl
0432 << "def nameToKey(name):" << std::endl
0433 << " match = re.match(\"ME([\\+\\-])([1-4])/([1-4])/([0-9]{2})\", name)" << std::endl
0434 << " if match is None: return None" << std::endl
0435 << " endcap, station, ring, chamber = match.groups()" << std::endl
0436 << " if endcap == \"+\": endcap = 1" << std::endl
0437 << " else: endcap = 2" << std::endl
0438 << " station = int(station)" << std::endl
0439 << " ring = int(ring)" << std::endl
0440 << " chamber = int(chamber)" << std::endl
0441 << " return endcap, station, ring, chamber" << std::endl
0442 << std::endl;
0443 }
0444
0445 for (std::vector<CSCChamberFitter>::const_iterator fitter = m_fitters.begin(); fitter != m_fitters.end();
0446 ++fitter) {
0447 if (m_mode == CSCPairResidualsConstraint::kModeRadius) {
0448 fitter->radiusCorrection(m_alignableNavigator, m_alignmentParameterStore, m_combineME11);
0449
0450 } else {
0451 std::vector<CSCAlignmentCorrections*> corrections;
0452 fitter->fit(corrections);
0453
0454
0455 for (std::vector<CSCAlignmentCorrections*>::iterator correction = corrections.begin();
0456 correction != corrections.end();
0457 ++correction) {
0458 (*correction)->applyAlignment(m_alignableNavigator, m_alignmentParameterStore, m_mode, m_combineME11);
0459 if (m_makeHistograms)
0460 (*correction)->plot();
0461 if (writeReport)
0462 (*correction)->report(report);
0463 }
0464 }
0465 }
0466 }
0467 }
0468
0469 #include "Alignment/CommonAlignmentAlgorithm/interface/AlignmentAlgorithmPluginFactory.h"
0470 DEFINE_EDM_PLUGIN(AlignmentAlgorithmPluginFactory, CSCOverlapsAlignmentAlgorithm, "CSCOverlapsAlignmentAlgorithm");
