CMS 3D CMS Logo

/afs/cern.ch/work/a/aaltunda/public/www/CMSSW_5_3_14/src/Alignment/MuonAlignmentAlgorithms/plugins/CSCOverlapsAlignmentAlgorithm.cc

Go to the documentation of this file.
00001 #include "Alignment/MuonAlignmentAlgorithms/plugins/CSCOverlapsAlignmentAlgorithm.h"
00002 
00003 CSCOverlapsAlignmentAlgorithm::CSCOverlapsAlignmentAlgorithm(const edm::ParameterSet& iConfig)
00004    : AlignmentAlgorithmBase(iConfig)
00005    , m_minHitsPerChamber(iConfig.getParameter<int>("minHitsPerChamber"))
00006    , m_maxdrdz(iConfig.getParameter<double>("maxdrdz"))
00007    , m_fiducial(iConfig.getParameter<bool>("fiducial"))
00008    , m_useHitWeights(iConfig.getParameter<bool>("useHitWeights"))
00009    , m_slopeFromTrackRefit(iConfig.getParameter<bool>("slopeFromTrackRefit"))
00010    , m_minStationsInTrackRefits(iConfig.getParameter<int>("minStationsInTrackRefits"))
00011    , m_truncateSlopeResid(iConfig.getParameter<double>("truncateSlopeResid"))
00012    , m_truncateOffsetResid(iConfig.getParameter<double>("truncateOffsetResid"))
00013    , m_combineME11(iConfig.getParameter<bool>("combineME11"))
00014    , m_useTrackWeights(iConfig.getParameter<bool>("useTrackWeights"))
00015    , m_errorFromRMS(iConfig.getParameter<bool>("errorFromRMS"))
00016    , m_minTracksPerOverlap(iConfig.getParameter<int>("minTracksPerOverlap"))
00017    , m_makeHistograms(iConfig.getParameter<bool>("makeHistograms"))
00018    , m_mode_string(iConfig.getParameter<std::string>("mode"))
00019    , m_reportFileName(iConfig.getParameter<std::string>("reportFileName"))
00020    , m_minP(iConfig.getParameter<double>("minP"))
00021    , m_maxRedChi2(iConfig.getParameter<double>("maxRedChi2"))
00022    , m_writeTemporaryFile(iConfig.getParameter<std::string>("writeTemporaryFile"))
00023    , m_readTemporaryFiles(iConfig.getParameter<std::vector<std::string> >("readTemporaryFiles"))
00024    , m_doAlignment(iConfig.getParameter<bool>("doAlignment"))
00025 {
00026   if (m_mode_string == std::string("phiy")) m_mode = CSCPairResidualsConstraint::kModePhiy;
00027   else if (m_mode_string == std::string("phipos")) m_mode = CSCPairResidualsConstraint::kModePhiPos;
00028   else if (m_mode_string == std::string("phiz")) m_mode = CSCPairResidualsConstraint::kModePhiz;
00029   else if (m_mode_string == std::string("radius")) m_mode = CSCPairResidualsConstraint::kModeRadius;
00030   else throw cms::Exception("BadConfig") << "mode must be one of \"phiy\", \"phipos\", \"phiz\", \"radius\"" << std::endl;
00031 
00032   std::vector<edm::ParameterSet> fitters = iConfig.getParameter<std::vector<edm::ParameterSet> >("fitters");
00033   for (std::vector<edm::ParameterSet>::const_iterator fitter = fitters.begin();  fitter != fitters.end();  ++fitter) {
00034     m_fitters.push_back(CSCChamberFitter(*fitter, m_residualsConstraints));
00035   }
00036 
00037   for (std::vector<CSCPairResidualsConstraint*>::const_iterator residualsConstraint = m_residualsConstraints.begin();  residualsConstraint != m_residualsConstraints.end();  ++residualsConstraint) {
00038     (*residualsConstraint)->configure(this);
00039     m_quickChamberLookup[std::pair<CSCDetId,CSCDetId>((*residualsConstraint)->id_i(), (*residualsConstraint)->id_j())] = *residualsConstraint;
00040   }
00041 
00042   if (m_slopeFromTrackRefit) {
00043     m_trackTransformer = new TrackTransformer(iConfig.getParameter<edm::ParameterSet>("TrackTransformer"));
00044     m_propagatorName = iConfig.getParameter<edm::ParameterSet>("TrackTransformer").getParameter<std::string>("Propagator");
00045   }
00046   else {
00047     m_trackTransformer = NULL;
00048     m_propagatorName = std::string("");
00049   }
00050 
00051   m_propagatorPointer = NULL;
00052 
00053   if (m_makeHistograms) {
00054     edm::Service<TFileService> tFileService;
00055     m_histP10 = tFileService->make<TH1F>("P10", "", 100, 0, 10);
00056     m_histP100 = tFileService->make<TH1F>("P100", "", 100, 0, 100);
00057     m_histP1000 = tFileService->make<TH1F>("P1000", "", 100, 0, 1000);
00058 
00059     m_hitsPerChamber = tFileService->make<TH1F>("hitsPerChamber", "", 10, -0.5, 9.5);
00060 
00061     m_fiducial_ME11 = tFileService->make<TProfile>("fiducial_ME11", "", 100, 0.075, 0.100);
00062     m_fiducial_ME12 = tFileService->make<TProfile>("fiducial_ME12", "", 100, 0.080, 0.105);
00063     m_fiducial_MEx1 = tFileService->make<TProfile>("fiducial_MEx1", "", 100, 0.160, 0.210);
00064     m_fiducial_MEx2 = tFileService->make<TProfile>("fiducial_MEx2", "", 100, 0.080, 0.105);
00065 
00066     m_slope = tFileService->make<TH1F>("slope", "", 100, -0.5, 0.5);
00067     m_slope_MEp4 = tFileService->make<TH1F>("slope_MEp4", "", 100, -0.5, 0.5);
00068     m_slope_MEp3 = tFileService->make<TH1F>("slope_MEp3", "", 100, -0.5, 0.5);
00069     m_slope_MEp2 = tFileService->make<TH1F>("slope_MEp2", "", 100, -0.5, 0.5);
00070     m_slope_MEp1 = tFileService->make<TH1F>("slope_MEp1", "", 100, -0.5, 0.5);
00071     m_slope_MEm1 = tFileService->make<TH1F>("slope_MEm1", "", 100, -0.5, 0.5);
00072     m_slope_MEm2 = tFileService->make<TH1F>("slope_MEm2", "", 100, -0.5, 0.5);
00073     m_slope_MEm3 = tFileService->make<TH1F>("slope_MEm3", "", 100, -0.5, 0.5);
00074     m_slope_MEm4 = tFileService->make<TH1F>("slope_MEm4", "", 100, -0.5, 0.5);
00075 
00076     m_slopeResiduals = tFileService->make<TH1F>("slopeResiduals", "mrad", 300, -30., 30.);
00077     m_slopeResiduals_weighted = tFileService->make<TH1F>("slopeResiduals_weighted", "mrad", 300, -30., 30.);
00078     m_slopeResiduals_normalized = tFileService->make<TH1F>("slopeResiduals_normalized", "", 200, -20., 20.);
00079     m_offsetResiduals = tFileService->make<TH1F>("offsetResiduals", "mm", 300, -30., 30.);
00080     m_offsetResiduals_weighted = tFileService->make<TH1F>("offsetResiduals_weighted", "mm", 300, -30., 30.);
00081     m_offsetResiduals_normalized = tFileService->make<TH1F>("offsetResiduals_normalized", "", 200, -20., 20.);
00082 
00083     m_drdz = tFileService->make<TH1F>("drdz", "", 100, -0.5, 0.5);
00084 
00085     m_occupancy = tFileService->make<TH2F>("occupancy", "", 36, 1, 37, 20, 1, 21);
00086     for (int i = 1;  i <= 36;  i++) {
00087       std::stringstream pairname;
00088       pairname << i << "-";
00089       if (i+1 == 37) pairname << 1;
00090       else pairname << (i+1);
00091       m_occupancy->GetXaxis()->SetBinLabel(i, pairname.str().c_str());
00092     }
00093     m_occupancy->GetYaxis()->SetBinLabel(1, "ME-4/2");
00094     m_occupancy->GetYaxis()->SetBinLabel(2, "ME-4/1");
00095     m_occupancy->GetYaxis()->SetBinLabel(3, "ME-3/2");
00096     m_occupancy->GetYaxis()->SetBinLabel(4, "ME-3/1");
00097     m_occupancy->GetYaxis()->SetBinLabel(5, "ME-2/2");
00098     m_occupancy->GetYaxis()->SetBinLabel(6, "ME-2/1");
00099     m_occupancy->GetYaxis()->SetBinLabel(7, "ME-1/3");
00100     m_occupancy->GetYaxis()->SetBinLabel(8, "ME-1/2");
00101     if (!m_combineME11) {
00102       m_occupancy->GetYaxis()->SetBinLabel(9, "ME-1/1b");
00103       m_occupancy->GetYaxis()->SetBinLabel(10, "ME-1/1a");
00104       m_occupancy->GetYaxis()->SetBinLabel(11, "ME+1/1a");
00105       m_occupancy->GetYaxis()->SetBinLabel(12, "ME+1/1b");
00106     }
00107     else {
00108       m_occupancy->GetYaxis()->SetBinLabel(9, "ME-1/1");
00109       m_occupancy->GetYaxis()->SetBinLabel(10, "");
00110       m_occupancy->GetYaxis()->SetBinLabel(11, "");
00111       m_occupancy->GetYaxis()->SetBinLabel(12, "ME+1/1");
00112     }
00113     m_occupancy->GetYaxis()->SetBinLabel(13, "ME+1/2");
00114     m_occupancy->GetYaxis()->SetBinLabel(14, "ME+1/3");
00115     m_occupancy->GetYaxis()->SetBinLabel(15, "ME+2/1");
00116     m_occupancy->GetYaxis()->SetBinLabel(16, "ME+2/2");
00117     m_occupancy->GetYaxis()->SetBinLabel(17, "ME+3/1");
00118     m_occupancy->GetYaxis()->SetBinLabel(18, "ME+3/2");
00119     m_occupancy->GetYaxis()->SetBinLabel(19, "ME+4/1");
00120     m_occupancy->GetYaxis()->SetBinLabel(20, "ME+4/2");
00121 
00122     m_XYpos_mep1 = tFileService->make<TH2F>("XYpos_mep1", "Positions: ME+1", 140, -700., 700., 140, -700., 700.);
00123     m_XYpos_mep2 = tFileService->make<TH2F>("XYpos_mep2", "Positions: ME+2", 140, -700., 700., 140, -700., 700.);
00124     m_XYpos_mep3 = tFileService->make<TH2F>("XYpos_mep3", "Positions: ME+3", 140, -700., 700., 140, -700., 700.);
00125     m_XYpos_mep4 = tFileService->make<TH2F>("XYpos_mep4", "Positions: ME+4", 140, -700., 700., 140, -700., 700.);
00126     m_XYpos_mem1 = tFileService->make<TH2F>("XYpos_mem1", "Positions: ME-1", 140, -700., 700., 140, -700., 700.);
00127     m_XYpos_mem2 = tFileService->make<TH2F>("XYpos_mem2", "Positions: ME-2", 140, -700., 700., 140, -700., 700.);
00128     m_XYpos_mem3 = tFileService->make<TH2F>("XYpos_mem3", "Positions: ME-3", 140, -700., 700., 140, -700., 700.);
00129     m_XYpos_mem4 = tFileService->make<TH2F>("XYpos_mem4", "Positions: ME-4", 140, -700., 700., 140, -700., 700.);
00130     m_RPhipos_mep1 = tFileService->make<TH2F>("RPhipos_mep1", "Positions: ME+1", 144, -M_PI, M_PI, 21, 0., 700.);
00131     m_RPhipos_mep2 = tFileService->make<TH2F>("RPhipos_mep2", "Positions: ME+2", 144, -M_PI, M_PI, 21, 0., 700.);
00132     m_RPhipos_mep3 = tFileService->make<TH2F>("RPhipos_mep3", "Positions: ME+3", 144, -M_PI, M_PI, 21, 0., 700.);
00133     m_RPhipos_mep4 = tFileService->make<TH2F>("RPhipos_mep4", "Positions: ME+4", 144, -M_PI, M_PI, 21, 0., 700.);
00134     m_RPhipos_mem1 = tFileService->make<TH2F>("RPhipos_mem1", "Positions: ME-1", 144, -M_PI, M_PI, 21, 0., 700.);
00135     m_RPhipos_mem2 = tFileService->make<TH2F>("RPhipos_mem2", "Positions: ME-2", 144, -M_PI, M_PI, 21, 0., 700.);
00136     m_RPhipos_mem3 = tFileService->make<TH2F>("RPhipos_mem3", "Positions: ME-3", 144, -M_PI, M_PI, 21, 0., 700.);
00137     m_RPhipos_mem4 = tFileService->make<TH2F>("RPhipos_mem4", "Positions: ME-4", 144, -M_PI, M_PI, 21, 0., 700.);
00138   }
00139   else {
00140     m_histP10 = NULL;
00141     m_histP100 = NULL;
00142     m_histP1000 = NULL;
00143     m_hitsPerChamber = NULL;
00144     m_fiducial_ME11 = NULL;
00145     m_fiducial_ME12 = NULL;
00146     m_fiducial_MEx1 = NULL;
00147     m_fiducial_MEx2 = NULL;
00148     m_slope = NULL;
00149     m_slope_MEp4 = NULL;
00150     m_slope_MEp3 = NULL;
00151     m_slope_MEp2 = NULL;
00152     m_slope_MEp1 = NULL;
00153     m_slope_MEm1 = NULL;
00154     m_slope_MEm2 = NULL;
00155     m_slope_MEm3 = NULL;
00156     m_slope_MEm4 = NULL;
00157     m_slopeResiduals = NULL;
00158     m_slopeResiduals_weighted = NULL;
00159     m_slopeResiduals_normalized = NULL;
00160     m_offsetResiduals = NULL;
00161     m_offsetResiduals_weighted = NULL;
00162     m_offsetResiduals_normalized = NULL;
00163     m_drdz = NULL;
00164     m_occupancy = NULL;
00165     m_XYpos_mep1 = NULL;
00166     m_XYpos_mep2 = NULL;
00167     m_XYpos_mep3 = NULL;
00168     m_XYpos_mep4 = NULL;
00169     m_XYpos_mem1 = NULL;
00170     m_XYpos_mem2 = NULL;
00171     m_XYpos_mem3 = NULL;
00172     m_XYpos_mem4 = NULL;
00173     m_RPhipos_mep1 = NULL;
00174     m_RPhipos_mep2 = NULL;
00175     m_RPhipos_mep3 = NULL;
00176     m_RPhipos_mep4 = NULL;
00177     m_RPhipos_mem1 = NULL;
00178     m_RPhipos_mem2 = NULL;
00179     m_RPhipos_mem3 = NULL;
00180     m_RPhipos_mem4 = NULL;
00181   }
00182 }
00183 
00184 CSCOverlapsAlignmentAlgorithm::~CSCOverlapsAlignmentAlgorithm() {}
00185 
00186 void CSCOverlapsAlignmentAlgorithm::initialize(const edm::EventSetup& iSetup, AlignableTracker* alignableTracker, AlignableMuon* alignableMuon, AlignableExtras* alignableExtras, AlignmentParameterStore* alignmentParameterStore) {
00187   m_alignmentParameterStore = alignmentParameterStore;
00188   m_alignables = m_alignmentParameterStore->alignables();
00189 
00190   if (alignableTracker == NULL) m_alignableNavigator = new AlignableNavigator(alignableMuon);
00191   else m_alignableNavigator = new AlignableNavigator(alignableTracker, alignableMuon);
00192 
00193   for (std::vector<Alignable*>::const_iterator alignable = m_alignables.begin();  alignable != m_alignables.end();  ++alignable) {
00194     DetId id = (*alignable)->geomDetId();
00195     if (id.det() != DetId::Muon  ||  id.subdetId() != MuonSubdetId::CSC  ||  CSCDetId(id.rawId()).layer() != 0) {
00196       throw cms::Exception("BadConfig") << "Only CSC chambers may be alignable" << std::endl;
00197     }
00198 
00199     std::vector<bool> selector = (*alignable)->alignmentParameters()->selector();
00200     for (std::vector<bool>::const_iterator i = selector.begin();  i != selector.end();  ++i) {
00201       if (!(*i)) throw cms::Exception("BadConfig") << "All selector strings should be \"111111\"" << std::endl;
00202     }
00203   }
00204 
00205   edm::ESHandle<CSCGeometry> cscGeometry;
00206   iSetup.get<MuonGeometryRecord>().get(cscGeometry);
00207 
00208   for (std::vector<CSCPairResidualsConstraint*>::const_iterator residualsConstraint = m_residualsConstraints.begin();  residualsConstraint != m_residualsConstraints.end();  ++residualsConstraint) {
00209     (*residualsConstraint)->setZplane(&*cscGeometry);
00210   }
00211 
00212   if (m_readTemporaryFiles.size() != 0) {
00213     std::vector<std::ifstream*> input;
00214     for (std::vector<std::string>::const_iterator fileName = m_readTemporaryFiles.begin();  fileName != m_readTemporaryFiles.end();  ++fileName) {
00215       input.push_back(new std::ifstream(fileName->c_str()));
00216     }
00217 
00218     for (std::vector<CSCPairResidualsConstraint*>::const_iterator residualsConstraint = m_residualsConstraints.begin();  residualsConstraint != m_residualsConstraints.end();  ++residualsConstraint) {
00219       (*residualsConstraint)->read(input, m_readTemporaryFiles);
00220     }
00221 
00222     for (std::vector<std::ifstream*>::const_iterator file = input.begin();  file != input.end();  ++file) {
00223       delete (*file);
00224     }
00225   }
00226 }
00227 
00228 void CSCOverlapsAlignmentAlgorithm::run(const edm::EventSetup& iSetup, const EventInfo &eventInfo) {
00229   edm::ESHandle<Propagator> propagator;
00230   if (m_slopeFromTrackRefit) {
00231     iSetup.get<TrackingComponentsRecord>().get(m_propagatorName, propagator);
00232     if (m_propagatorPointer != &*propagator) {
00233       m_propagatorPointer = &*propagator;
00234 
00235       for (std::vector<CSCPairResidualsConstraint*>::const_iterator residualsConstraint = m_residualsConstraints.begin();  residualsConstraint != m_residualsConstraints.end();  ++residualsConstraint) {
00236         (*residualsConstraint)->setPropagator(m_propagatorPointer);
00237       }
00238     }
00239   }
00240 
00241   edm::ESHandle<TransientTrackBuilder> transientTrackBuilder;
00242   iSetup.get<TransientTrackRecord>().get("TransientTrackBuilder", transientTrackBuilder);
00243 
00244   if (m_trackTransformer != NULL) m_trackTransformer->setServices(iSetup);
00245 
00246   const ConstTrajTrackPairCollection &trajtracks = eventInfo.trajTrackPairs_;
00247   for (ConstTrajTrackPairCollection::const_iterator trajtrack = trajtracks.begin();  trajtrack != trajtracks.end();  ++trajtrack) {
00248     const Trajectory* traj = (*trajtrack).first;
00249     const reco::Track* track = (*trajtrack).second;
00250 
00251     if (m_makeHistograms) {
00252       m_histP10->Fill(track->p());
00253       m_histP100->Fill(track->p());
00254       m_histP1000->Fill(track->p());
00255     }
00256     if (track->p() >= m_minP) {
00257       std::vector<TrajectoryMeasurement> measurements = traj->measurements();
00258       reco::TransientTrack transientTrack = transientTrackBuilder->build(track);
00259 
00260       std::map<int,std::map<CSCDetId,bool> > stationsToChambers;
00261       for (std::vector<TrajectoryMeasurement>::const_iterator measurement = measurements.begin();  measurement != measurements.end();  ++measurement) {
00262         DetId id = measurement->recHit()->geographicalId();
00263         if (id.det() == DetId::Muon  &&  id.subdetId() == MuonSubdetId::CSC) {
00264           CSCDetId cscid(id.rawId());
00265           CSCDetId chamberId(cscid.endcap(), cscid.station(), cscid.ring(), cscid.chamber(), 0);
00266           if (m_combineME11  &&  cscid.station() == 1  &&  cscid.ring() == 4) chamberId = CSCDetId(cscid.endcap(), 1, 1, cscid.chamber(), 0);
00267           int station = (cscid.endcap() == 1 ? 1 : -1)*cscid.station();
00268 
00269           if (stationsToChambers.find(station) == stationsToChambers.end()) stationsToChambers[station] = std::map<CSCDetId,bool>();
00270           stationsToChambers[station][chamberId] = true;
00271 
00272           if (m_makeHistograms) {
00273             GlobalPoint pos = measurement->recHit()->globalPosition();
00274             if (cscid.endcap() == 1  &&  cscid.station() == 1) { m_XYpos_mep1->Fill(pos.x(), pos.y()); m_RPhipos_mep1->Fill(pos.phi(), pos.perp()); }
00275             if (cscid.endcap() == 1  &&  cscid.station() == 2) { m_XYpos_mep2->Fill(pos.x(), pos.y()); m_RPhipos_mep2->Fill(pos.phi(), pos.perp()); }
00276             if (cscid.endcap() == 1  &&  cscid.station() == 3) { m_XYpos_mep3->Fill(pos.x(), pos.y()); m_RPhipos_mep3->Fill(pos.phi(), pos.perp()); }
00277             if (cscid.endcap() == 1  &&  cscid.station() == 4) { m_XYpos_mep4->Fill(pos.x(), pos.y()); m_RPhipos_mep4->Fill(pos.phi(), pos.perp()); }
00278             if (cscid.endcap() == 2  &&  cscid.station() == 1) { m_XYpos_mem1->Fill(pos.x(), pos.y()); m_RPhipos_mem1->Fill(pos.phi(), pos.perp()); }
00279             if (cscid.endcap() == 2  &&  cscid.station() == 2) { m_XYpos_mem2->Fill(pos.x(), pos.y()); m_RPhipos_mem2->Fill(pos.phi(), pos.perp()); }
00280             if (cscid.endcap() == 2  &&  cscid.station() == 3) { m_XYpos_mem3->Fill(pos.x(), pos.y()); m_RPhipos_mem3->Fill(pos.phi(), pos.perp()); }
00281             if (cscid.endcap() == 2  &&  cscid.station() == 4) { m_XYpos_mem4->Fill(pos.x(), pos.y()); m_RPhipos_mem4->Fill(pos.phi(), pos.perp()); }
00282           }
00283         }
00284       }
00285       
00286       std::map<CSCPairResidualsConstraint*,bool> residualsConstraints;
00287       for (std::map<int,std::map<CSCDetId,bool> >::const_iterator iter = stationsToChambers.begin();  iter != stationsToChambers.end();  ++iter) {
00288         for (std::map<CSCDetId,bool>::const_iterator one = iter->second.begin();  one != iter->second.end();  ++one) {
00289           for (std::map<CSCDetId,bool>::const_iterator two = one;  two != iter->second.end();  ++two) {
00290             if (one != two) {
00291               std::map<std::pair<CSCDetId,CSCDetId>,CSCPairResidualsConstraint*>::const_iterator quick;
00292 
00293               quick = m_quickChamberLookup.find(std::pair<CSCDetId,CSCDetId>(one->first, two->first));
00294               if (quick != m_quickChamberLookup.end()) residualsConstraints[quick->second] = true;
00295 
00296               quick = m_quickChamberLookup.find(std::pair<CSCDetId,CSCDetId>(two->first, one->first));
00297               if (quick != m_quickChamberLookup.end()) residualsConstraints[quick->second] = true;
00298             }
00299           }
00300         }
00301       }
00302 
00303       for (std::map<CSCPairResidualsConstraint*,bool>::const_iterator residualsConstraint = residualsConstraints.begin();  residualsConstraint != residualsConstraints.end();  ++residualsConstraint) {
00304         residualsConstraint->first->addTrack(measurements, transientTrack, m_trackTransformer);
00305       }
00306     }
00307   }
00308 }
00309 
00310 void CSCOverlapsAlignmentAlgorithm::terminate() {
00311   // write residuals partial fits to temporary files for collection
00312   if (m_writeTemporaryFile != std::string("")) {
00313     std::ofstream output(m_writeTemporaryFile.c_str());
00314     for (std::vector<CSCPairResidualsConstraint*>::const_iterator residualsConstraint = m_residualsConstraints.begin();  residualsConstraint != m_residualsConstraints.end();  ++residualsConstraint) {
00315       (*residualsConstraint)->write(output);
00316     }
00317   }
00318 
00319   // write report for alignment results
00320   if (m_doAlignment) {
00321     std::ofstream report;
00322     bool writeReport = (m_reportFileName != std::string(""));
00323     if (writeReport) {
00324       report.open(m_reportFileName.c_str());
00325       report << "cscReports = []" << std::endl << std::endl
00326              << "class CSCChamberCorrection:" << std::endl
00327              << "    def __init__(self, name, detid, value):" << std::endl
00328              << "        self.name, self.detid, self.value = name, detid, value" << std::endl << std::endl
00329              << "class CSCErrorMode:" << std::endl
00330              << "    def __init__(self, error):" << std::endl
00331              << "        self.error = error" << std::endl
00332              << "        self.terms = {}" << std::endl
00333              << "        self.detids = {}" << std::endl
00334              << "    def addTerm(self, name, detid, coefficient):" << std::endl
00335              << "        self.terms[name] = coefficient" << std::endl
00336              << "        self.detids[name] = detid" << std::endl << std::endl
00337              << "class CSCConstraintResidual:" << std::endl
00338              << "    def __init__(self, i, j, before, uncert, residual, pull):" << std::endl
00339              << "        self.i, self.j, self.before, self.error, self.residual, self.pull = i, j, before, uncert, residual, pull" << std::endl << std::endl
00340              << "class CSCFitterReport:" << std::endl
00341              << "    def __init__(self, name, oldchi2, newchi2):" << std::endl
00342              << "        self.name, self.oldchi2, self.newchi2 = name, oldchi2, newchi2" << std::endl
00343              << "        self.chamberCorrections = []" << std::endl
00344              << "        self.errorModes = []" << std::endl
00345              << "        self.constraintResiduals = []" << std::endl << std::endl
00346              << "    def addChamberCorrection(self, name, detid, value):" << std::endl
00347              << "        self.chamberCorrections.append(CSCChamberCorrection(name, detid, value))" << std::endl << std::endl
00348              << "    def addErrorMode(self, error):" << std::endl
00349              << "        self.errorModes.append(CSCErrorMode(error))" << std::endl << std::endl
00350              << "    def addErrorModeTerm(self, name, detid, coefficient):" << std::endl
00351              << "        self.errorModes[-1].addTerm(name, detid, coefficient)" << std::endl << std::endl
00352              << "    def addCSCConstraintResidual(self, i, j, before, uncert, residual, pull):" << std::endl
00353              << "        self.constraintResiduals.append(CSCConstraintResidual(i, j, before, uncert, residual, pull))" << std::endl << std::endl
00354              << "import re" << std::endl
00355              << "def nameToKey(name):" << std::endl
00356              << "    match = re.match(\"ME([\\+\\-])([1-4])/([1-4])/([0-9]{2})\", name)" << std::endl
00357              << "    if match is None: return None" << std::endl
00358              << "    endcap, station, ring, chamber = match.groups()" << std::endl
00359              << "    if endcap == \"+\": endcap = 1" << std::endl
00360              << "    else: endcap = 2" << std::endl
00361              << "    station = int(station)" << std::endl
00362              << "    ring = int(ring)" << std::endl
00363              << "    chamber = int(chamber)" << std::endl
00364              << "    return endcap, station, ring, chamber" << std::endl << std::endl;
00365     }
00366 
00367     for (std::vector<CSCChamberFitter>::const_iterator fitter = m_fitters.begin();  fitter != m_fitters.end();  ++fitter) {
00368       if (m_mode == CSCPairResidualsConstraint::kModeRadius) {
00369         fitter->radiusCorrection(m_alignableNavigator, m_alignmentParameterStore, m_combineME11);
00370         
00371 
00372 
00373       }
00374       else {
00375         std::vector<CSCAlignmentCorrections*> corrections;
00376         fitter->fit(corrections);
00377          
00378         // corrections only exist if the fit was successful
00379         for (std::vector<CSCAlignmentCorrections*>::iterator correction = corrections.begin();  correction != corrections.end();  ++correction) {
00380            
00381            (*correction)->applyAlignment(m_alignableNavigator, m_alignmentParameterStore, m_mode, m_combineME11);
00382            if (m_makeHistograms) (*correction)->plot();
00383            if (writeReport) (*correction)->report(report);
00384         }
00385       }
00386     }
00387   }
00388 }
00389 
00390 #include "Alignment/CommonAlignmentAlgorithm/interface/AlignmentAlgorithmPluginFactory.h"
00391 DEFINE_EDM_PLUGIN(AlignmentAlgorithmPluginFactory, CSCOverlapsAlignmentAlgorithm, "CSCOverlapsAlignmentAlgorithm");