#include <HIPAlignmentAlgorithm.h>
Definition at line 24 of file HIPAlignmentAlgorithm.h.
HIPAlignmentAlgorithm::HIPAlignmentAlgorithm | ( | const edm::ParameterSet & | cfg | ) |
Constructor.
Definition at line 33 of file HIPAlignmentAlgorithm.cc.
References edm::ParameterSet::getParameter(), edm::ParameterSet::getUntrackedParameter(), isCollector, prof2calltree::l, AlignableObjectId::nameToType(), outfile, outfile2, outpath, salignedfile, siterationfile, smisalignedfile, sparameterfile, ssurveyfile, struefile, suvarfile, theAPEParameterSet, theApplyAPE, theCollectorNJobs, theCollectorPath, theCurrentPrescale, theEventPrescale, theFillTrackMonitoring, theLevels, theMaxAllowedHitPull, theMaxRelParameterError, and theMinimumNumberOfHits.
: AlignmentAlgorithmBase( cfg ) { // parse parameters verbose = cfg.getParameter<bool>("verbosity"); outpath = cfg.getParameter<string>("outpath"); outfile = cfg.getParameter<string>("outfile"); outfile2 = cfg.getParameter<string>("outfile2"); struefile = cfg.getParameter<string>("trueFile"); smisalignedfile = cfg.getParameter<string>("misalignedFile"); salignedfile = cfg.getParameter<string>("alignedFile"); siterationfile = cfg.getParameter<string>("iterationFile"); suvarfile = cfg.getParameter<string>("uvarFile"); sparameterfile = cfg.getParameter<string>("parameterFile"); ssurveyfile = cfg.getParameter<string>("surveyFile"); outfile =outpath+outfile;//Eventwise tree outfile2 =outpath+outfile2;//Alignablewise tree struefile =outpath+struefile; smisalignedfile=outpath+smisalignedfile; salignedfile =outpath+salignedfile; siterationfile =outpath+siterationfile; suvarfile =outpath+suvarfile; sparameterfile =outpath+sparameterfile; ssurveyfile =outpath+ssurveyfile; // parameters for APE theApplyAPE = cfg.getParameter<bool>("applyAPE"); theAPEParameterSet = cfg.getParameter<std::vector<edm::ParameterSet> >("apeParam"); theMaxAllowedHitPull = cfg.getParameter<double>("maxAllowedHitPull"); theMinimumNumberOfHits = cfg.getParameter<int>("minimumNumberOfHits"); theMaxRelParameterError = cfg.getParameter<double>("maxRelParameterError"); // for collector mode (parallel processing) isCollector=cfg.getParameter<bool>("collectorActive"); theCollectorNJobs=cfg.getParameter<int>("collectorNJobs"); theCollectorPath=cfg.getParameter<string>("collectorPath"); theFillTrackMonitoring=cfg.getUntrackedParameter<bool>("fillTrackMonitoring"); if (isCollector) edm::LogWarning("Alignment") << "[HIPAlignmentAlgorithm] Collector mode"; theEventPrescale = cfg.getParameter<int>("eventPrescale"); theCurrentPrescale = theEventPrescale; AlignableObjectId dummy; const std::vector<std::string>& levels = cfg.getUntrackedParameter<std::vector<std::string> >("surveyResiduals"); for (unsigned int l = 0; l < levels.size(); ++l) { theLevels.push_back( dummy.nameToType(levels[l]) ); } edm::LogWarning("Alignment") << "[HIPAlignmentAlgorithm] constructed."; }
HIPAlignmentAlgorithm::~HIPAlignmentAlgorithm | ( | ) | [inline] |
void HIPAlignmentAlgorithm::bookRoot | ( | void | ) | [private] |
Definition at line 1029 of file HIPAlignmentAlgorithm.cc.
References m2_Eta, m2_Id, m2_Layer, m2_Nhit, m2_ObjId, m2_Phi, m2_Type, m2_Xpos, m2_Ypos, m2_Zpos, m3_Id, m3_ObjId, m3_par, m_Chi2n, m_d0, m_dz, m_Eta, m_Nhits, m_nhPXB, m_nhPXF, m_Ntracks, m_P, m_Phi, m_Pt, outfile, outfile2, ssurveyfile, theFile, theFile2, theFile3, theIteration, theLevels, theTree, theTree2, and theTree3.
Referenced by startNewLoop().
{ // create ROOT files theFile = new TFile(outfile.c_str(),"update"); theFile->cd(); // book event-wise ROOT Tree TString tname="T1"; char iterString[5]; snprintf(iterString, sizeof(iterString), "%i",theIteration); tname.Append("_"); tname.Append(iterString); theTree = new TTree(tname,"Eventwise tree"); //theTree->Branch("Run", &m_Run, "Run/I"); //theTree->Branch("Event", &m_Event, "Event/I"); theTree->Branch("Ntracks", &m_Ntracks, "Ntracks/I"); theTree->Branch("Nhits", m_Nhits, "Nhits[Ntracks]/I"); theTree->Branch("nhPXB", m_nhPXB, "nhPXB[Ntracks]/I"); theTree->Branch("nhPXF", m_nhPXF, "nhPXF[Ntracks]/I"); theTree->Branch("Pt", m_Pt, "Pt[Ntracks]/F"); theTree->Branch("P", m_P, "P[Ntracks]/F"); theTree->Branch("Eta", m_Eta, "Eta[Ntracks]/F"); theTree->Branch("Phi", m_Phi, "Phi[Ntracks]/F"); theTree->Branch("Chi2n", m_Chi2n, "Chi2n[Ntracks]/F"); theTree->Branch("d0", m_d0, "d0[Ntracks]/F"); theTree->Branch("dz", m_dz, "dz[Ntracks]/F"); // book Alignable-wise ROOT Tree theFile2 = new TFile(outfile2.c_str(),"update"); theFile2->cd(); theTree2 = new TTree("T2","Alignablewise tree"); theTree2->Branch("Nhit", &m2_Nhit, "Nhit/I"); theTree2->Branch("Type", &m2_Type, "Type/I"); theTree2->Branch("Layer", &m2_Layer, "Layer/I"); theTree2->Branch("Xpos", &m2_Xpos, "Xpos/F"); theTree2->Branch("Ypos", &m2_Ypos, "Ypos/F"); theTree2->Branch("Zpos", &m2_Zpos, "Zpos/F"); theTree2->Branch("Eta", &m2_Eta, "Eta/F"); theTree2->Branch("Phi", &m2_Phi, "Phi/F"); theTree2->Branch("Id", &m2_Id, "Id/i"); theTree2->Branch("ObjId", &m2_ObjId, "ObjId/I"); // book survey-wise ROOT Tree only if survey is enabled if (theLevels.size() > 0){ edm::LogWarning("Alignment") << "[HIPAlignmentAlgorithm::bookRoot] Survey trees booked."; theFile3 = new TFile(ssurveyfile.c_str(),"update"); theFile3->cd(); theTree3 = new TTree(tname, "Survey Tree"); theTree3->Branch("Id", &m3_Id, "Id/i"); theTree3->Branch("ObjId", &m3_ObjId, "ObjId/I"); theTree3->Branch("Par", &m3_par, "Par[6]/F"); } edm::LogWarning("Alignment") << "[HIPAlignmentAlgorithm::bookRoot] Root trees booked."; }
double HIPAlignmentAlgorithm::calcAPE | ( | double * | par, |
int | iter, | ||
double | function | ||
) | [private] |
Definition at line 1001 of file HIPAlignmentAlgorithm.cc.
References create_public_lumi_plots::exp, max(), funct::pow(), and launcher::step.
Referenced by setAlignmentPositionError().
{ double diter=(double)iter; // The following floating-point equality check is safe because this // "0." and this "1." are generated by the compiler, in the very // same file. Whatever approximization scheme it uses to turn "1." // into 0.9999999999998 in the HIPAlignmentAlgorithm::initialize is // also used here. If I'm wrong, you'll get an assertion. if (function == 0.) { return max(par[1],par[0]+((par[1]-par[0])/par[2])*diter); } else if (function == 1.) { return max(0.,par[0]*(exp(-pow(diter,par[1])/par[2]))); } else if (function == 2.) { int ipar2 = (int) par[2]; int step = iter/ipar2; double dstep = (double) step; return max(0.0, par[0] - par[1]*dstep); } else assert(false); // should have been caught in the constructor }
bool HIPAlignmentAlgorithm::calcParameters | ( | Alignable * | ali | ) | [private] |
Definition at line 1171 of file HIPAlignmentAlgorithm.cc.
References abs, Alignable::alignmentParameters(), AlignmentParameters::cloneFromSelected(), i, HIPUserVariables::jtve, HIPUserVariables::jtvj, HIPUserVariables::nhit, Alignable::setAlignmentParameters(), AlignmentParameters::setValid(), mathSSE::sqrt(), theMaxRelParameterError, theMinimumNumberOfHits, and AlignmentParameters::userVariables().
Referenced by terminate().
{ // Alignment parameters AlignmentParameters* par = ali->alignmentParameters(); // access user variables HIPUserVariables* uservar = dynamic_cast<HIPUserVariables*>(par->userVariables()); int nhit = uservar->nhit; // The following variable is needed for the extended 1D/2D hit fix using // matrix shrinkage and expansion // int hitdim = uservar->hitdim; if (nhit < theMinimumNumberOfHits) { par->setValid(false); return false; } AlgebraicSymMatrix jtvj = uservar->jtvj; AlgebraicVector jtve = uservar->jtve; // Shrink input in case of 1D hits and 'v' selected // in alignment parameters // if (hitdim==1 && selector[1]==true) { // int iremove = 1; // if (selector[0]==false) iremove--; // AlgebraicSymMatrix tempjtvj(jtvj.num_row()-1); // int nr = 0, nc = 0; // for (int r=0;r<jtvj.num_row();r++) { // if (r==iremove) continue; // nc = 0; // for (int c=0;c<jtvj.num_col();c++) { // if (c==iremove) continue; // tempjtvj[nr][nc] = jtvj[r][c]; // nc++; // } // nr++; // } // jtvj = tempjtvj; // AlgebraicVector tempjtve(jtve.num_row()-1); // nr = 0; // for (int r=0;r<jtve.num_row();r++) { // if (r==iremove) continue; // tempjtve[nr] = jtve[r]; // nr++; // } // jtve = tempjtve; // } int ierr; AlgebraicSymMatrix jtvjinv = jtvj.inverse(ierr); if (ierr !=0) { edm::LogError("Alignment") << "Matrix inversion failed!"; return false; } // these are the alignment corrections+covariance (for selected params) AlgebraicVector params = - (jtvjinv * jtve); AlgebraicSymMatrix cov = jtvjinv; edm::LogInfo("Alignment") << "parameters " << params; // errors of parameters int npar = params.num_row(); AlgebraicVector paramerr(npar); AlgebraicVector relerr(npar); for (int i=0;i<npar;i++) { if (abs(cov[i][i])>0) paramerr[i] = sqrt(abs(cov[i][i])); else paramerr[i] = params[i]; relerr[i] = abs(paramerr[i]/params[i]); if (relerr[i] >= theMaxRelParameterError) { params[i] = 0; paramerr[i]=0; } } // expand output in case of 1D hits and 'v' selected // in alignment parameters // if (hitdim==1 && selector[1]==true) { // int iremove = 1; // if (selector[0]==false) iremove--; // AlgebraicSymMatrix tempcov(cov.num_row()+1,0); // int nr = 0, nc = 0; // for (int r=0;r<cov.num_row();r++) { // if (r==iremove) nr++; // nc = 0; // for (int c=0;c<cov.num_col();c++) { // if (c==iremove) nc++; // tempcov[nr][nc] = cov[r][c]; // nc++; // } // nr++; // } // cov = tempcov; // AlgebraicVector tempparams(params.num_row()+1,0); // nr = 0; // for (int r=0;r<params.num_row();r++) { // if (r==iremove) nr++; // tempparams[nr] = params[r]; // nr++; // } // params = tempparams; // } // store alignment parameters AlignmentParameters* parnew = par->cloneFromSelected(params,cov); ali->setAlignmentParameters(parnew); parnew->setValid(true); return true; }
void HIPAlignmentAlgorithm::collector | ( | void | ) | [private] |
Definition at line 1288 of file HIPAlignmentAlgorithm.cc.
References HIPUserVariables::alichi2, Alignable::alignmentParameters(), HIPUserVariables::alindof, AlignmentParameterStore::attachUserVariables(), HIPUserVariables::clone(), fillEventwiseTree(), ioerr, HIPUserVariables::jtve, HIPUserVariables::jtvj, HIPUserVariables::nhit, HIPUserVariablesIORoot::readHIPUserVariables(), theAlignables, theAlignmentParameterStore, theCollectorNJobs, theCollectorPath, theFillTrackMonitoring, theIteration, and AlignmentParameters::userVariables().
Referenced by startNewLoop().
{ edm::LogWarning("Alignment") << "[HIPAlignmentAlgorithm::collector] called for iteration " << theIteration << std::endl; HIPUserVariablesIORoot HIPIO; for (int ijob=1;ijob<=theCollectorNJobs;ijob++) { edm::LogWarning("Alignment") << "reading uservar for job " << ijob; stringstream ss; string str; ss << ijob; ss >> str; string uvfile = theCollectorPath+"/job"+str+"/IOUserVariables.root"; vector<AlignmentUserVariables*> uvarvec = HIPIO.readHIPUserVariables(theAlignables, (char*)uvfile.c_str(), theIteration, ioerr); if (ioerr!=0) { edm::LogWarning("Alignment") << "[HIPAlignmentAlgorithm::collector] could not read user variable files for job " << ijob; continue; } // add vector<AlignmentUserVariables*> uvarvecadd; vector<AlignmentUserVariables*>::const_iterator iuvarnew=uvarvec.begin(); for (vector<Alignable*>::const_iterator it=theAlignables.begin(); it!=theAlignables.end(); ++it) { Alignable* ali = *it; AlignmentParameters* ap = ali->alignmentParameters(); HIPUserVariables* uvarold = dynamic_cast<HIPUserVariables*>(ap->userVariables()); HIPUserVariables* uvarnew = dynamic_cast<HIPUserVariables*>(*iuvarnew); HIPUserVariables* uvar = uvarold->clone(); if (uvarnew!=0) { uvar->nhit = (uvarold->nhit)+(uvarnew->nhit); uvar->jtvj = (uvarold->jtvj)+(uvarnew->jtvj); uvar->jtve = (uvarold->jtve)+(uvarnew->jtve); uvar->alichi2 = (uvarold->alichi2)+(uvarnew->alichi2); uvar->alindof = (uvarold->alindof)+(uvarnew->alindof); delete uvarnew; } uvarvecadd.push_back(uvar); iuvarnew++; } theAlignmentParameterStore->attachUserVariables(theAlignables, uvarvecadd, ioerr); //fill Eventwise Tree if (theFillTrackMonitoring) { uvfile = theCollectorPath+"/job"+str+"/HIPAlignmentEvents.root"; edm::LogWarning("Alignment") << "Added to the tree " << fillEventwiseTree(uvfile.c_str(), theIteration, ioerr) << "tracks"; } }//end loop on jobs }
int HIPAlignmentAlgorithm::fillEventwiseTree | ( | const char * | filename, |
int | iter, | ||
int | ierr | ||
) | [private] |
Definition at line 1356 of file HIPAlignmentAlgorithm.cc.
References m_Chi2n, m_d0, m_dz, m_Eta, m_Nhits, m_nhPXB, m_nhPXF, m_Ntracks, m_P, m_Phi, m_Pt, MAXREC, and theTree.
Referenced by collector().
{ int totntrk = 0; char treeName[64]; snprintf(treeName, sizeof(treeName), "T1_%d", iter); //open the file "HIPAlignmentEvents.root" in the job directory TFile *jobfile = new TFile(filename, "READ"); //grab the tree corresponding to this iteration TTree *jobtree = (TTree*)jobfile->Get(treeName); //address and read the variables static const int nmaxtrackperevent = 1000; int jobNtracks, jobNhitspertrack[nmaxtrackperevent], jobnhPXB[nmaxtrackperevent], jobnhPXF[nmaxtrackperevent]; float jobP[nmaxtrackperevent], jobPt[nmaxtrackperevent], jobEta[nmaxtrackperevent] , jobPhi[nmaxtrackperevent]; float jobd0[nmaxtrackperevent], jobdz[nmaxtrackperevent] , jobChi2n[nmaxtrackperevent]; jobtree->SetBranchAddress("Ntracks", &jobNtracks); jobtree->SetBranchAddress("Nhits", jobNhitspertrack); jobtree->SetBranchAddress("nhPXB", jobnhPXB); jobtree->SetBranchAddress("nhPXF", jobnhPXF); jobtree->SetBranchAddress("Pt", jobPt); jobtree->SetBranchAddress("P", jobP); jobtree->SetBranchAddress("d0", jobd0); jobtree->SetBranchAddress("dz", jobdz); jobtree->SetBranchAddress("Eta", jobEta); jobtree->SetBranchAddress("Phi", jobPhi); jobtree->SetBranchAddress("Chi2n", jobChi2n); int ievent = 0; for (ievent=0;ievent<jobtree->GetEntries();++ievent) { jobtree->GetEntry(ievent); //fill the collector tree with them // TO BE IMPLEMENTED: a prescale factor like in run() m_Ntracks = jobNtracks; int ntrk = 0; while (ntrk<m_Ntracks) { if (ntrk<MAXREC) { totntrk = ntrk+1; m_Nhits[ntrk] = jobNhitspertrack[ntrk]; m_Pt[ntrk] = jobPt[ntrk]; m_P[ntrk] = jobP[ntrk]; m_nhPXB[ntrk] = jobnhPXB[ntrk]; m_nhPXF[ntrk] = jobnhPXF[ntrk]; m_Eta[ntrk] = jobEta[ntrk]; m_Phi[ntrk] = jobPhi[ntrk]; m_Chi2n[ntrk] = jobChi2n[ntrk]; m_d0[ntrk] = jobd0[ntrk]; m_dz[ntrk] = jobdz[ntrk]; }//end if j<MAXREC else{ edm::LogWarning("Alignment") << "[HIPAlignmentAlgorithm::fillEventwiseTree] Number of tracks in Eventwise tree exceeds MAXREC: " << m_Ntracks << " Skipping exceeding tracks."; ntrk = m_Ntracks+1; } ++ntrk; }//end while loop theTree->Fill(); }//end loop on i - entries in the job tree //clean up delete jobtree; delete jobfile; return totntrk; }//end fillEventwiseTree
void HIPAlignmentAlgorithm::fillRoot | ( | void | ) | [private] |
Definition at line 1096 of file HIPAlignmentAlgorithm.cc.
References Alignable::alignableObjectId(), Alignable::alignmentParameters(), PV3DBase< T, PVType, FrameType >::eta(), Alignable::id(), AlignmentParameters::isValid(), m2_Eta, m2_Id, m2_Layer, m2_Nhit, m2_ObjId, m2_Phi, m2_Type, m2_Xpos, m2_Ypos, m2_Zpos, HIPUserVariables::nhit, AlignmentParameters::parameters(), PV3DBase< T, PVType, FrameType >::phi(), pos, GloballyPositioned< T >::position(), Alignable::surface(), theAlignables, theAlignmentParameterStore, theFile2, theTree2, AlignmentParameterStore::typeAndLayer(), AlignmentParameters::userVariables(), PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().
Referenced by terminate().
{ theFile2->cd(); int naligned=0; for (vector<Alignable*>::const_iterator it=theAlignables.begin(); it!=theAlignables.end(); ++it) { Alignable* ali = (*it); AlignmentParameters* dap = ali->alignmentParameters(); // consider only those parameters classified as 'valid' if (dap->isValid()) { // get number of hits from user variable HIPUserVariables* uservar = dynamic_cast<HIPUserVariables*>(dap->userVariables()); m2_Nhit = uservar->nhit; // get type/layer std::pair<int,int> tl = theAlignmentParameterStore->typeAndLayer(ali); m2_Type = tl.first; m2_Layer = tl.second; // get identifier (as for IO) m2_Id = ali->id(); m2_ObjId = ali->alignableObjectId(); // get position GlobalPoint pos = ali->surface().position(); m2_Xpos = pos.x(); m2_Ypos = pos.y(); m2_Zpos = pos.z(); m2_Eta = pos.eta(); m2_Phi = pos.phi(); AlgebraicVector pars = dap->parameters(); if (verbose) { edm::LogVerbatim("Alignment") << "------------------------------------------------------------------------\n" << " ALIGNABLE: " << setw(6) << naligned << '\n' << "hits: " << setw(4) << m2_Nhit << " type: " << setw(4) << m2_Type << " layer: " << setw(4) << m2_Layer << " id: " << setw(4) << m2_Id << " objId: " << setw(4) << m2_ObjId << '\n' << fixed << setprecision(5) << "x,y,z: " << setw(12) << m2_Xpos << setw(12) << m2_Ypos << setw(12) << m2_Zpos << " eta,phi: " << setw(12) << m2_Eta << setw(12) << m2_Phi << '\n' << "params: " << setw(12) << pars[0] << setw(12) << pars[1] << setw(12) << pars[2] << setw(12) << pars[3] << setw(12) << pars[4] << setw(12) << pars[5]; } naligned++; theTree2->Fill(); } } }
void HIPAlignmentAlgorithm::initialize | ( | const edm::EventSetup & | setup, |
AlignableTracker * | tracker, | ||
AlignableMuon * | muon, | ||
AlignableExtras * | extras, | ||
AlignmentParameterStore * | store | ||
) | [virtual] |
Call at beginning of job.
Implements AlignmentAlgorithmBase.
Definition at line 96 of file HIPAlignmentAlgorithm.cc.
References AlignmentParameterSelector::addSelections(), AlignmentParameterStore::alignables(), MuonAlignmentFromReference_cff::alignParams, AlignmentParameterSelector::clear(), Exception, edm::ParameterSet::getParameter(), i, AlignmentParameterSelector::selectedAlignables(), theAlignableDetAccessor, theAlignables, theAlignmentParameterStore, theAPEParameters, theAPEParameterSet, and theApplyAPE.
{ edm::LogWarning("Alignment") << "[HIPAlignmentAlgorithm] Initializing..."; // accessor Det->AlignableDet if ( !muon ) theAlignableDetAccessor = new AlignableNavigator(tracker); else if ( !tracker ) theAlignableDetAccessor = new AlignableNavigator(muon); else theAlignableDetAccessor = new AlignableNavigator(tracker, muon); // set alignmentParameterStore theAlignmentParameterStore=store; // get alignables theAlignables = theAlignmentParameterStore->alignables(); // clear theAPEParameters, if necessary theAPEParameters.clear(); // get APE parameters if(theApplyAPE){ AlignmentParameterSelector selector(tracker, muon); for (std::vector<edm::ParameterSet>::const_iterator setiter = theAPEParameterSet.begin(); setiter != theAPEParameterSet.end(); ++setiter) { std::vector<Alignable*> alignables; selector.clear(); edm::ParameterSet selectorPSet = setiter->getParameter<edm::ParameterSet>("Selector"); std::vector<std::string> alignParams = selectorPSet.getParameter<std::vector<std::string> >("alignParams"); if (alignParams.size() == 1 && alignParams[0] == std::string("selected")) { alignables = theAlignables; } else { selector.addSelections(selectorPSet); alignables = selector.selectedAlignables(); } std::vector<double> apeSPar = setiter->getParameter<std::vector<double> >("apeSPar"); std::vector<double> apeRPar = setiter->getParameter<std::vector<double> >("apeRPar"); std::string function = setiter->getParameter<std::string>("function"); if (apeSPar.size() != 3 || apeRPar.size() != 3) throw cms::Exception("BadConfig") << "apeSPar and apeRPar must have 3 values each" << std::endl; for (std::vector<double>::const_iterator i = apeRPar.begin(); i != apeRPar.end(); ++i) { apeSPar.push_back(*i); } if (function == std::string("linear")) { apeSPar.push_back(0.); // c.f. note in calcAPE } else if (function == std::string("exponential")) { apeSPar.push_back(1.); // c.f. note in calcAPE } else if (function == std::string("step")) { apeSPar.push_back(2.); // c.f. note in calcAPE } else { throw cms::Exception("BadConfig") << "APE function must be \"linear\" or \"exponential\"." << std::endl; } theAPEParameters.push_back(std::pair<std::vector<Alignable*>, std::vector<double> >(alignables, apeSPar)); } } }
bool HIPAlignmentAlgorithm::processHit1D | ( | const AlignableDetOrUnitPtr | alidet, |
const Alignable * | ali, | ||
const TrajectoryStateOnSurface & | tsos, | ||
const TransientTrackingRecHit * | hit | ||
) | [private] |
Definition at line 434 of file HIPAlignmentAlgorithm.cc.
References HIPUserVariables::alichi2, Alignable::alignmentParameters(), HIPUserVariables::alindof, HIPUserVariables::jtve, HIPUserVariables::jtvj, TrajectoryStateOnSurface::localError(), TrackingRecHit::localPosition(), TrajectoryStateOnSurface::localPosition(), TrackingRecHit::localPositionError(), HIPUserVariables::nhit, pos, LocalTrajectoryError::positionError(), AlignmentParameters::selectedDerivatives(), mathSSE::sqrt(), theMaxAllowedHitPull, AlignmentParameters::userVariables(), PV3DBase< T, PVType, FrameType >::x(), and LocalError::xx().
Referenced by run().
{ static const unsigned int hitDim = 1; // get trajectory impact point LocalPoint alvec = tsos.localPosition(); AlgebraicVector pos(hitDim); pos[0] = alvec.x(); // get impact point covariance AlgebraicSymMatrix ipcovmat(hitDim); ipcovmat[0][0] = tsos.localError().positionError().xx(); // get hit local position and covariance AlgebraicVector coor(hitDim); coor[0] = hit->localPosition().x(); AlgebraicSymMatrix covmat(hitDim); covmat[0][0] = hit->localPositionError().xx(); // add hit and impact point covariance matrices covmat = covmat + ipcovmat; // calculate the x pull of this hit double xpull = 0.; if (covmat[0][0] != 0.) xpull = (pos[0] - coor[0])/sqrt(fabs(covmat[0][0])); // get Alignment Parameters AlignmentParameters* params = ali->alignmentParameters(); // get derivatives AlgebraicMatrix derivs2D = params->selectedDerivatives(tsos, alidet); // calculate user parameters int npar = derivs2D.num_row(); // std::cout << "Dumping derivsSEL: \n" << derivs2D <<std::endl; AlgebraicMatrix derivs(npar, hitDim, 0); for (int ipar=0;ipar<npar;ipar++) { for (unsigned int idim=0;idim<hitDim;idim++) { derivs[ipar][idim] = derivs2D[ipar][idim]; } } // std::cout << "Dumping derivs: \n" << derivs << std::endl; // invert covariance matrix int ierr; // if (covmat[1][1]>4.0) nhitDim = hitDim; covmat.invert(ierr); if (ierr != 0) { edm::LogError("Alignment") << "Matrix inversion failed!"; return false; } bool useThisHit = (theMaxAllowedHitPull <= 0.); useThisHit |= (fabs(xpull) < theMaxAllowedHitPull); // bailing out if (!useThisHit) return false; // std::cout << "We use this hit in " << subDet << std::endl; // std::cout << "Npars= " << npar << " NhitDim=1 Size of derivs: " // << derivs.num_row() << " x " << derivs.num_col() << std::endl; // AlgebraicMatrix jtvjtmp(); AlgebraicMatrix covtmp(covmat); // std::cout << "Preapring JTVJTMP -> " << std::flush; AlgebraicMatrix jtvjtmp(derivs * covtmp *derivs.T()); // std::cout << "TMP JTVJ= \n" << jtvjtmp << std::endl; AlgebraicSymMatrix thisjtvj(npar); AlgebraicVector thisjtve(npar); // std::cout << "Preparing ThisJtVJ" << std::flush; // thisjtvj = covmat.similarity(derivs); thisjtvj.assign(jtvjtmp); // cout<<" ThisJtVE"<<std::endl; thisjtve=derivs * covmat * (pos-coor); AlgebraicVector hitresidual(hitDim); hitresidual[0] = (pos[0] - coor[0]); AlgebraicMatrix hitresidualT; hitresidualT = hitresidual.T(); // std::cout << "HitResidualT = \n" << hitresidualT << std::endl; // access user variables (via AlignmentParameters) HIPUserVariables* uservar = dynamic_cast<HIPUserVariables*>(params->userVariables()); uservar->jtvj += thisjtvj; uservar->jtve += thisjtve; uservar->nhit++; // The following variable is needed for the extended 1D/2D hit fix using // matrix shrinkage and expansion // uservar->hitdim = hitDim; //for alignable chi squared float thischi2; thischi2 = (hitresidualT *covmat *hitresidual)[0]; if ( verbose && ((thischi2/ static_cast<float>(uservar->nhit)) >10.0) ) { edm::LogWarning("Alignment") << "Added to Chi2 the number " << thischi2 <<" having " << uservar->nhit << " dof " << std::endl << "X-resid " << hitresidual[0] << " Y-resid " << hitresidual[1] << std::endl << " Cov^-1 matr (covmat): [0][0]= " << covmat[0][0] << " [0][1]= " << covmat[0][1] << " [1][0]= " << covmat[1][0] << " [1][1]= " << covmat[1][1] << std::endl; } uservar->alichi2 += thischi2; // a bit weird, but vector.transposed * matrix * vector doesn't give a double in CMSSW's opinion uservar->alindof += hitDim; // 2D hits contribute twice to the ndofs return true; }
bool HIPAlignmentAlgorithm::processHit2D | ( | const AlignableDetOrUnitPtr | alidet, |
const Alignable * | ali, | ||
const TrajectoryStateOnSurface & | tsos, | ||
const TransientTrackingRecHit * | hit | ||
) | [private] |
Definition at line 550 of file HIPAlignmentAlgorithm.cc.
References HIPUserVariables::alichi2, Alignable::alignmentParameters(), HIPUserVariables::alindof, HIPUserVariables::jtve, HIPUserVariables::jtvj, TrajectoryStateOnSurface::localError(), TrackingRecHit::localPosition(), TrajectoryStateOnSurface::localPosition(), TrackingRecHit::localPositionError(), HIPUserVariables::nhit, pos, LocalTrajectoryError::positionError(), AlignmentParameters::selectedDerivatives(), mathSSE::sqrt(), theMaxAllowedHitPull, AlignmentParameters::userVariables(), PV3DBase< T, PVType, FrameType >::x(), LocalError::xx(), LocalError::xy(), PV3DBase< T, PVType, FrameType >::y(), and LocalError::yy().
Referenced by run().
{ static const unsigned int hitDim = 2; // get trajectory impact point LocalPoint alvec = tsos.localPosition(); AlgebraicVector pos(hitDim); pos[0] = alvec.x(); pos[1] = alvec.y(); // get impact point covariance AlgebraicSymMatrix ipcovmat(hitDim); ipcovmat[0][0] = tsos.localError().positionError().xx(); ipcovmat[1][1] = tsos.localError().positionError().yy(); ipcovmat[0][1] = tsos.localError().positionError().xy(); // get hit local position and covariance AlgebraicVector coor(hitDim); coor[0] = hit->localPosition().x(); coor[1] = hit->localPosition().y(); AlgebraicSymMatrix covmat(hitDim); covmat[0][0] = hit->localPositionError().xx(); covmat[1][1] = hit->localPositionError().yy(); covmat[0][1] = hit->localPositionError().xy(); // add hit and impact point covariance matrices covmat = covmat + ipcovmat; // calculate the x pull and y pull of this hit double xpull = 0.; double ypull = 0.; if (covmat[0][0] != 0.) xpull = (pos[0] - coor[0])/sqrt(fabs(covmat[0][0])); if (covmat[1][1] != 0.) ypull = (pos[1] - coor[1])/sqrt(fabs(covmat[1][1])); // get Alignment Parameters AlignmentParameters* params = ali->alignmentParameters(); // get derivatives AlgebraicMatrix derivs2D = params->selectedDerivatives(tsos, alidet); // calculate user parameters int npar = derivs2D.num_row(); // std::cout << "Dumping derivsSEL: \n"<< derivs2D <<std::endl; AlgebraicMatrix derivs(npar, hitDim, 0); for (int ipar=0;ipar<npar;ipar++) { for (unsigned int idim=0;idim<hitDim;idim++) { derivs[ipar][idim] = derivs2D[ipar][idim]; } } // std::cout << "Dumping derivs: \n" << derivs << std::endl; // invert covariance matrix int ierr; // if (covmat[1][1]>4.0) nhitDim = hitDim; covmat.invert(ierr); if (ierr != 0) { edm::LogError("Alignment") << "Matrix inversion failed!"; return false; } bool useThisHit = (theMaxAllowedHitPull <= 0.); useThisHit |= (fabs(xpull) < theMaxAllowedHitPull && fabs(ypull) < theMaxAllowedHitPull); // bailing out if (!useThisHit) return false; // std::cout << "We use this hit in " << subDet << std::endl; // std::cout << "Npars= " << npar << " NhitDim=2 Size of derivs: " // << derivs.num_row() << " x " << derivs.num_col() << std::endl; // AlgebraicMatrix jtvjtmp(); AlgebraicMatrix covtmp(covmat); // std::cout << "Preapring JTVJTMP -> " << std::flush; AlgebraicMatrix jtvjtmp(derivs * covtmp *derivs.T()); // std::cout << "TMP JTVJ= \n" << jtvjtmp << std::endl; AlgebraicSymMatrix thisjtvj(npar); AlgebraicVector thisjtve(npar); // std::cout << "Preparing ThisJtVJ" << std::flush; // thisjtvj = covmat.similarity(derivs); thisjtvj.assign(jtvjtmp); // cout<<" ThisJtVE"<<std::endl; thisjtve=derivs * covmat * (pos-coor); AlgebraicVector hitresidual(hitDim); hitresidual[0] = (pos[0] - coor[0]); hitresidual[1] = (pos[1] - coor[1]); // if(nhitDim>1) { // hitresidual[1] =0.0; // nhitDim=1; // } AlgebraicMatrix hitresidualT; hitresidualT = hitresidual.T(); // std::cout << "HitResidualT = \n" << hitresidualT << std::endl; // access user variables (via AlignmentParameters) HIPUserVariables* uservar = dynamic_cast<HIPUserVariables*>(params->userVariables()); uservar->jtvj += thisjtvj; uservar->jtve += thisjtve; uservar->nhit++; // The following variable is needed for the extended 1D/2D hit fix using // matrix shrinkage and expansion // uservar->hitdim = hitDim; //for alignable chi squared float thischi2; thischi2 = (hitresidualT *covmat *hitresidual)[0]; if ( verbose && ((thischi2/ static_cast<float>(uservar->nhit)) >10.0) ) { edm::LogWarning("Alignment") << "Added to Chi2 the number " << thischi2 <<" having " << uservar->nhit << " dof " << std::endl << "X-resid " << hitresidual[0] << " Y-resid " << hitresidual[1] << std::endl << " Cov^-1 matr (covmat): [0][0]= " << covmat[0][0] << " [0][1]= " << covmat[0][1] << " [1][0]= " << covmat[1][0] << " [1][1]= " << covmat[1][1] << std::endl; } uservar->alichi2 += thischi2; // a bit weird, but vector.transposed * matrix * vector doesn't give a double in CMSSW's opinion uservar->alindof += hitDim; // 2D hits contribute twice to the ndofs return true; }
int HIPAlignmentAlgorithm::readIterationFile | ( | std::string | filename | ) | [private] |
Definition at line 909 of file HIPAlignmentAlgorithm.cc.
References recoMuon::in, and query::result.
Referenced by startNewLoop().
{ int result; ifstream inIterFile((char*)filename.c_str(), ios::in); if (!inIterFile) { edm::LogError("Alignment") << "[HIPAlignmentAlgorithm::readIterationFile] ERROR! " << "Unable to open Iteration file"; result = -1; } else { inIterFile >> result; edm::LogWarning("Alignment") << "[HIPAlignmentAlgorithm::readIterationFile] " << "Read last iteration number from file: " << result; } inIterFile.close(); return result; }
void HIPAlignmentAlgorithm::run | ( | const edm::EventSetup & | setup, |
const EventInfo & | eventInfo | ||
) | [virtual] |
Run the algorithm.
Implements AlignmentAlgorithmBase.
Definition at line 679 of file HIPAlignmentAlgorithm.cc.
References CompositeAlignmentParameters::alignableFromAlignableDet(), AlignableNavigator::alignableFromGeomDet(), AlignableNavigator::alignablesFromHits(), TrajectoryMeasurement::backwardPredictedState(), className(), AlignmentAlgorithmBase::EventInfo::clusterValueMap_, TrajectoryStateCombiner::combine(), reco::TrackBase::d0(), AlignableNavigator::detAndSubdetInMap(), reco::TrackBase::dz(), eta(), reco::TrackBase::eta(), Exception, TrajectoryMeasurement::forwardPredictedState(), TrackingRecHit::geographicalId(), TransientTrackingRecHit::hit(), reco::TrackBase::hitPattern(), isCollector, AlignmentClusterFlag::isTaken(), TrackingRecHit::isValid(), TrajectoryStateOnSurface::isValid(), m_Chi2n, m_d0, m_dz, m_Eta, m_Nhits, m_nhPXB, m_nhPXF, m_Ntracks, m_P, m_Phi, m_Pt, MAXREC, Trajectory::measurements(), reco::TrackBase::normalizedChi2(), reco::TrackBase::numberOfValidHits(), reco::TrackBase::p(), AlCaHLTBitMon_ParallelJobs::p, phi, reco::TrackBase::phi(), processHit1D(), processHit2D(), reco::TrackBase::pt(), TrajectoryMeasurement::recHit(), AlignmentParameterStore::selectParameters(), DetId::subdetId(), theAlignableDetAccessor, theAlignmentParameterStore, theCurrentPrescale, theEventPrescale, theFile, theFillTrackMonitoring, theTree, testEve_cfg::tracks, and AlignmentAlgorithmBase::EventInfo::trajTrackPairs_.
{ if (isCollector) return; TrajectoryStateCombiner tsoscomb; // AM: not really needed // AM: m_Ntracks = 0 should be sufficient int itr=0; m_Ntracks=0; for(itr=0;itr<MAXREC;++itr){ m_Nhits[itr]=0; m_Pt[itr]=-5.0; m_P[itr]=-5.0; m_nhPXB[itr]=0; m_nhPXF[itr]=0; m_Eta[itr]=-99.0; m_Phi[itr]=-4.0; m_Chi2n[itr]=-11.0; m_d0[itr]=-999; m_dz[itr]=-999; } itr=0; // AM: what is this needed for? theFile->cd(); // loop over tracks const ConstTrajTrackPairCollection &tracks = eventInfo.trajTrackPairs_; for (ConstTrajTrackPairCollection::const_iterator it=tracks.begin(); it!=tracks.end(); ++it) { const Trajectory* traj = (*it).first; const reco::Track* track = (*it).second; float pt = track->pt(); float eta = track->eta(); float phi = track->phi(); float p = track->p(); float chi2n = track->normalizedChi2(); int nhit = track->numberOfValidHits(); float d0 = track->d0(); float dz = track->dz(); int nhpxb = track->hitPattern().numberOfValidPixelBarrelHits(); int nhpxf = track->hitPattern().numberOfValidPixelEndcapHits(); if (verbose) edm::LogInfo("Alignment") << "New track pt,eta,phi,chi2n,hits: " << pt << "," << eta << "," << phi << "," << chi2n << "," << nhit; // edm::LogWarning("Alignment") << "New track pt,eta,phi,chi2n,hits: " // << pt << "," // << eta << "," // << phi << "," // << chi2n << "," // << nhit; // fill track parameters in root tree if (itr<MAXREC) { m_Nhits[itr]=nhit; m_Pt[itr]=pt; m_P[itr]=p; m_Eta[itr]=eta; m_Phi[itr]=phi; m_Chi2n[itr]=chi2n; m_nhPXB[itr]=nhpxb; m_nhPXF[itr]=nhpxf; m_d0[itr]=d0; m_dz[itr]=dz; itr++; m_Ntracks=itr; } // AM: Can be simplified vector<const TransientTrackingRecHit*> hitvec; vector<TrajectoryStateOnSurface> tsosvec; // loop over measurements vector<TrajectoryMeasurement> measurements = traj->measurements(); for (vector<TrajectoryMeasurement>::iterator im=measurements.begin(); im!=measurements.end(); ++im) { TrajectoryMeasurement meas = *im; // const TransientTrackingRecHit* ttrhit = &(*meas.recHit()); // const TrackingRecHit *hit = ttrhit->hit(); const TransientTrackingRecHit* hit = &(*meas.recHit()); if (hit->isValid() && theAlignableDetAccessor->detAndSubdetInMap( hit->geographicalId() )) { // this is the updated state (including the current hit) // TrajectoryStateOnSurface tsos=meas.updatedState(); // combine fwd and bwd predicted state to get state // which excludes current hit bool skiphit = false; if (eventInfo.clusterValueMap_) { // check from the PrescalingMap if the hit was taken. // If not skip to the next TM // bool hitTaken=false; AlignmentClusterFlag myflag; int subDet = hit->geographicalId().subdetId(); //take the actual RecHit out of the Transient one const TrackingRecHit *rechit=hit->hit(); if (subDet>2) { // AM: if possible use enum instead of hard-coded value const std::type_info &type = typeid(*rechit); if (type == typeid(SiStripRecHit1D)) { const SiStripRecHit1D* stripHit1D = dynamic_cast<const SiStripRecHit1D*>(rechit); if (stripHit1D) { SiStripRecHit1D::ClusterRef stripclust(stripHit1D->cluster()); // myflag=PrescMap[stripclust]; myflag = (*eventInfo.clusterValueMap_)[stripclust]; } else { edm::LogError("HIPAlignmentAlgorithm") << "ERROR in <HIPAlignmentAlgorithm::run>: Dynamic cast of Strip RecHit failed! " << "TypeId of the RecHit: " << className(*hit) <<endl; } }//end if type = SiStripRecHit1D else if(type == typeid(SiStripRecHit2D)){ const SiStripRecHit2D* stripHit2D = dynamic_cast<const SiStripRecHit2D*>(rechit); if (stripHit2D) { SiStripRecHit2D::ClusterRef stripclust(stripHit2D->cluster()); // myflag=PrescMap[stripclust]; myflag = (*eventInfo.clusterValueMap_)[stripclust]; } else { edm::LogError("HIPAlignmentAlgorithm") << "ERROR in <HIPAlignmentAlgorithm::run>: Dynamic cast of Strip RecHit failed! " // << "TypeId of the TTRH: " << className(*ttrhit) << endl; << "TypeId of the TTRH: " << className(*hit) << endl; } } //end if type == SiStripRecHit2D } //end if hit from strips else { const SiPixelRecHit* pixelhit= dynamic_cast<const SiPixelRecHit*>(rechit); if (pixelhit) { SiPixelClusterRefNew pixelclust(pixelhit->cluster()); // myflag=PrescMap[pixelclust]; myflag = (*eventInfo.clusterValueMap_)[pixelclust]; } else { edm::LogError("HIPAlignmentAlgorithm") << "ERROR in <HIPAlignmentAlgorithm::run>: Dynamic cast of Pixel RecHit failed! " // << "TypeId of the TTRH: " << className(*ttrhit) << endl; << "TypeId of the TTRH: " << className(*hit) << endl; } } //end 'else' it is a pixel hit // bool hitTaken=myflag.isTaken(); if (!myflag.isTaken()) { skiphit=true; //cout<<"Hit from subdet "<<rechit->geographicalId().subdetId()<<" prescaled out."<<endl; continue; } }//end if Prescaled Hits if (skiphit) { throw cms::Exception("LogicError") << "ERROR in <HIPAlignmentAlgorithm::run>: this hit should have been skipped!" << endl; } TrajectoryStateOnSurface tsos = tsoscomb.combine(meas.forwardPredictedState(), meas.backwardPredictedState()); if(tsos.isValid()){ // hitvec.push_back(ttrhit); hitvec.push_back(hit); tsosvec.push_back(tsos); } } //hit valid } // transform RecHit vector to AlignableDet vector vector <AlignableDetOrUnitPtr> alidetvec = theAlignableDetAccessor->alignablesFromHits(hitvec); // get concatenated alignment parameters for list of alignables CompositeAlignmentParameters aap = theAlignmentParameterStore->selectParameters(alidetvec); vector<TrajectoryStateOnSurface>::const_iterator itsos=tsosvec.begin(); vector<const TransientTrackingRecHit*>::const_iterator ihit=hitvec.begin(); // loop over vectors(hit,tsos) while (itsos != tsosvec.end()) { // get AlignableDet for this hit const GeomDet* det = (*ihit)->det(); // int subDet= (*ihit)->geographicalId().subdetId(); uint32_t nhitDim = (*ihit)->dimension(); AlignableDetOrUnitPtr alidet = theAlignableDetAccessor->alignableFromGeomDet(det); // get relevant Alignable Alignable* ali = aap.alignableFromAlignableDet(alidet); if (ali!=0) { if (nhitDim==1) { processHit1D(alidet, ali, *itsos, *ihit); } else if (nhitDim==2) { processHit2D(alidet, ali, *itsos, *ihit); } } itsos++; ihit++; } } // end of track loop // fill eventwise root tree (with prescale defined in pset) if (theFillTrackMonitoring) { theCurrentPrescale--; if (theCurrentPrescale<=0) { theTree->Fill(); theCurrentPrescale = theEventPrescale; } } }
void HIPAlignmentAlgorithm::setAlignmentPositionError | ( | void | ) | [private] |
Definition at line 948 of file HIPAlignmentAlgorithm.cc.
References calcAPE(), i, AlignmentParameterStore::setAlignmentPositionError(), theAlignmentParameterStore, theAPEParameters, theApplyAPE, and theIteration.
Referenced by startNewLoop().
{ // Check if user wants to override APE if ( !theApplyAPE ) { edm::LogWarning("Alignment") <<"[HIPAlignmentAlgorithm::setAlignmentPositionError] No APE applied"; return; // NO APE APPLIED } edm::LogWarning("Alignment") <<"[HIPAlignmentAlgorithm::setAlignmentPositionError] Apply APE!"; double apeSPar[3], apeRPar[3]; for (std::vector<std::pair<std::vector<Alignable*>, std::vector<double> > >::const_iterator alipars = theAPEParameters.begin(); alipars != theAPEParameters.end(); ++alipars) { const std::vector<Alignable*> &alignables = alipars->first; const std::vector<double> &pars = alipars->second; apeSPar[0] = pars[0]; apeSPar[1] = pars[1]; apeSPar[2] = pars[2]; apeRPar[0] = pars[3]; apeRPar[1] = pars[4]; apeRPar[2] = pars[5]; double function = pars[6]; // Printout for debug printf("APE: %u alignables\n", (unsigned int)alignables.size()); for ( int i=0; i<21; ++i ) { double apelinstest=calcAPE(apeSPar,i,0.); double apeexpstest=calcAPE(apeSPar,i,1.); double apelinrtest=calcAPE(apeRPar,i,0.); double apeexprtest=calcAPE(apeRPar,i,1.); printf("APE: iter slin sexp rlin rexp: %5d %12.5f %12.5f %12.5f %12.5f\n", i,apelinstest,apeexpstest,apelinrtest,apeexprtest); } // set APE double apeshift=calcAPE(apeSPar,theIteration,function); double aperot =calcAPE(apeRPar,theIteration,function); theAlignmentParameterStore->setAlignmentPositionError( alignables, apeshift, aperot ); } }
void HIPAlignmentAlgorithm::startNewLoop | ( | void | ) | [virtual] |
Called at start of new loop.
Reimplemented from AlignmentAlgorithmBase.
Definition at line 168 of file HIPAlignmentAlgorithm.cc.
References AlignmentParameterStore::applyAlignableAbsolutePositions(), bookRoot(), collector(), ioerr, isCollector, AlignmentParameters::numSelected(), AlignmentIORoot::readAlignableAbsolutePositions(), readIterationFile(), salignedfile, setAlignmentPositionError(), AlignmentParameters::setUserVariables(), siterationfile, smisalignedfile, struefile, theAlignables, theAlignmentParameterStore, theIO, theIteration, AlignmentIORoot::writeAlignableAbsolutePositions(), and AlignmentIORoot::writeAlignableOriginalPositions().
{ // iterate over all alignables and attach user variables for( vector<Alignable*>::const_iterator it=theAlignables.begin(); it!=theAlignables.end(); it++ ) { AlignmentParameters* ap = (*it)->alignmentParameters(); int npar=ap->numSelected(); HIPUserVariables* userpar = new HIPUserVariables(npar); ap->setUserVariables(userpar); } // try to read in alignment parameters from a previous iteration AlignablePositions theAlignablePositionsFromFile = theIO.readAlignableAbsolutePositions(theAlignables, (char*)salignedfile.c_str(),-1,ioerr); int numAlignablesFromFile = theAlignablePositionsFromFile.size(); if (numAlignablesFromFile==0) { // file not there: first iteration // set iteration number to 1 if (isCollector) theIteration=0; else theIteration=1; edm::LogWarning("Alignment") << "[HIPAlignmentAlgorithm] File not found => iteration "<<theIteration; // get true (de-misaligned positions) and write to root file // hardcoded iteration=1 theIO.writeAlignableOriginalPositions(theAlignables, (char*)struefile.c_str(),1,false,ioerr); // get misaligned positions and write to root file // hardcoded iteration=1 theIO.writeAlignableAbsolutePositions(theAlignables, (char*)smisalignedfile.c_str(),1,false,ioerr); } else { // there have been previous iterations edm::LogWarning("Alignment") << "[HIPAlignmentAlgorithm] Alignables Read " << numAlignablesFromFile; // get iteration number from file theIteration = readIterationFile(siterationfile); // increase iteration theIteration++; edm::LogWarning("Alignment") <<"[HIPAlignmentAlgorithm] Iteration increased by one!"; // now apply psotions of file from prev iteration edm::LogWarning("Alignment") <<"[HIPAlignmentAlgorithm] Apply positions from file ..."; theAlignmentParameterStore->applyAlignableAbsolutePositions(theAlignables, theAlignablePositionsFromFile,ioerr); } edm::LogWarning("Alignment") <<"[HIPAlignmentAlgorithm] Current Iteration number: " << theIteration; // book root trees bookRoot(); /*---------------------moved to terminate------------------------------ if (theLevels.size() > 0) { edm::LogWarning("Alignment") << "[HIPAlignmentAlgorithm] Using survey constraint"; unsigned int nAlignable = theAlignables.size(); for (unsigned int i = 0; i < nAlignable; ++i) { const Alignable* ali = theAlignables[i]; AlignmentParameters* ap = ali->alignmentParameters(); HIPUserVariables* uservar = dynamic_cast<HIPUserVariables*>(ap->userVariables()); for (unsigned int l = 0; l < theLevels.size(); ++l) { SurveyResidual res(*ali, theLevels[l], true); if ( res.valid() ) { AlgebraicSymMatrix invCov = res.inverseCovariance(); // variable for tree AlgebraicVector sensResid = res.sensorResidual(); m3_Id = ali->id(); m3_ObjId = theLevels[l]; m3_par[0] = sensResid[0]; m3_par[1] = sensResid[1]; m3_par[2] = sensResid[2]; m3_par[3] = sensResid[3]; m3_par[4] = sensResid[4]; m3_par[5] = sensResid[5]; uservar->jtvj += invCov; uservar->jtve += invCov * sensResid; theTree3->Fill(); } } // align::LocalVectors residuals = res1.pointsResidual(); // unsigned int nPoints = residuals.size(); // for (unsigned int k = 0; k < nPoints; ++k) // { // AlgebraicMatrix J = term->survey()->derivatives(k); // AlgebraicVector e(3); // local residual // const align::LocalVector& lr = residuals[k]; // e(1) = lr.x(); e(2) = lr.y(); e(3) = lr.z(); // uservar->jtvj += invCov1.similarity(J); // uservar->jtve += J * (invCov1 * e); // } } } //------------------------------------------------------*/ // set alignment position error setAlignmentPositionError(); // run collector job if we are in parallel mode if (isCollector) collector(); }
void HIPAlignmentAlgorithm::terminate | ( | void | ) | [virtual] |
Call at end of job.
Implements AlignmentAlgorithmBase.
Definition at line 303 of file HIPAlignmentAlgorithm.cc.
References Alignable::alignmentParameters(), AlignmentParameterStore::applyParameters(), calcParameters(), fillRoot(), i, Alignable::id(), SurveyResidual::inverseCovariance(), ioerr, HIPUserVariables::jtve, HIPUserVariables::jtvj, prof2calltree::l, m3_Id, m3_ObjId, m3_par, salignedfile, SurveyResidual::sensorResidual(), AlignmentParameters::setValid(), siterationfile, sparameterfile, suvarfile, run_regression::test, theAlignables, theAlignmentParameterStore, theFile, theFile2, theFile3, theIO, theIteration, theLevels, theTree, theTree2, theTree3, AlignmentParameters::userVariables(), SurveyResidual::valid(), AlignmentIORoot::writeAlignableAbsolutePositions(), AlignmentIORoot::writeAlignmentParameters(), HIPUserVariablesIORoot::writeHIPUserVariables(), and writeIterationFile().
{ edm::LogWarning("Alignment") << "[HIPAlignmentAlgorithm] Terminating"; // calculating survey residuals if (theLevels.size() > 0) { edm::LogWarning("Alignment") << "[HIPAlignmentAlgorithm] Using survey constraint"; unsigned int nAlignable = theAlignables.size(); for (unsigned int i = 0; i < nAlignable; ++i) { const Alignable* ali = theAlignables[i]; AlignmentParameters* ap = ali->alignmentParameters(); HIPUserVariables* uservar = dynamic_cast<HIPUserVariables*>(ap->userVariables()); for (unsigned int l = 0; l < theLevels.size(); ++l) { SurveyResidual res(*ali, theLevels[l], true); if ( res.valid() ) { AlgebraicSymMatrix invCov = res.inverseCovariance(); // variable for tree AlgebraicVector sensResid = res.sensorResidual(); m3_Id = ali->id(); m3_ObjId = theLevels[l]; m3_par[0] = sensResid[0]; m3_par[1] = sensResid[1]; m3_par[2] = sensResid[2]; m3_par[3] = sensResid[3]; m3_par[4] = sensResid[4]; m3_par[5] = sensResid[5]; uservar->jtvj += invCov; uservar->jtve += invCov * sensResid; theTree3->Fill(); } } // align::LocalVectors residuals = res1.pointsResidual(); // unsigned int nPoints = residuals.size(); // for (unsigned int k = 0; k < nPoints; ++k) // { // AlgebraicMatrix J = term->survey()->derivatives(k); // AlgebraicVector e(3); // local residual // const align::LocalVector& lr = residuals[k]; // e(1) = lr.x(); e(2) = lr.y(); e(3) = lr.z(); // uservar->jtvj += invCov1.similarity(J); // uservar->jtve += J * (invCov1 * e); // } } } // write user variables HIPUserVariablesIORoot HIPIO; HIPIO.writeHIPUserVariables (theAlignables,(char*)suvarfile.c_str(), theIteration,false,ioerr); // now calculate alignment corrections ... int ialigned=0; // iterate over alignment parameters for(vector<Alignable*>::const_iterator it=theAlignables.begin(); it!=theAlignables.end(); it++) { Alignable* ali=(*it); // Alignment parameters AlignmentParameters* par = ali->alignmentParameters(); // try to calculate parameters bool test = calcParameters(ali); // if successful, apply parameters if (test) { edm::LogInfo("Alignment") << "now apply params"; theAlignmentParameterStore->applyParameters(ali); // set these parameters 'valid' ali->alignmentParameters()->setValid(true); // increase counter ialigned++; } else par->setValid(false); } edm::LogWarning("Alignment") << "[HIPAlignmentAlgorithm::terminate] Aligned units: " << ialigned; // fill alignable wise root tree fillRoot(); edm::LogWarning("Alignment") << "[HIPAlignmentAlgorithm] Writing aligned parameters to file: " << theAlignables.size(); // write new absolute positions to disk theIO.writeAlignableAbsolutePositions(theAlignables, (char*)salignedfile.c_str(),theIteration,false,ioerr); // write alignment parameters to disk theIO.writeAlignmentParameters(theAlignables, (char*)sparameterfile.c_str(),theIteration,false,ioerr); // write iteration number to file writeIterationFile(siterationfile,theIteration); // write out trees and close root file // eventwise tree theFile->cd(); theTree->Write(); delete theFile; if (theLevels.size() > 0){ theFile3->cd(); theTree3->Write(); delete theFile3; } // alignable-wise tree is only filled once if (theIteration==1) { // only for 1st iteration theFile2->cd(); theTree2->Write(); delete theFile2; } }
void HIPAlignmentAlgorithm::writeIterationFile | ( | std::string | filename, |
int | iter | ||
) | [private] |
Definition at line 931 of file HIPAlignmentAlgorithm.cc.
References dbtoconf::out.
Referenced by terminate().
{ ofstream outIterFile((char*)(filename.c_str()), ios::out); if (!outIterFile) { edm::LogError("Alignment") << "[HIPAlignmentAlgorithm::writeIterationFile] ERROR: Unable to write Iteration file"; } else { outIterFile << iter; edm::LogWarning("Alignment") <<"[HIPAlignmentAlgorithm::writeIterationFile] writing iteration number to file: " << iter; } outIterFile.close(); }
int HIPAlignmentAlgorithm::ioerr [private] |
Definition at line 79 of file HIPAlignmentAlgorithm.h.
Referenced by collector(), startNewLoop(), and terminate().
bool HIPAlignmentAlgorithm::isCollector [private] |
Definition at line 101 of file HIPAlignmentAlgorithm.h.
Referenced by HIPAlignmentAlgorithm(), run(), and startNewLoop().
float HIPAlignmentAlgorithm::m2_Eta [private] |
Definition at line 125 of file HIPAlignmentAlgorithm.h.
Referenced by bookRoot(), and fillRoot().
align::ID HIPAlignmentAlgorithm::m2_Id [private] |
Definition at line 126 of file HIPAlignmentAlgorithm.h.
Referenced by bookRoot(), and fillRoot().
int HIPAlignmentAlgorithm::m2_Layer [private] |
Definition at line 124 of file HIPAlignmentAlgorithm.h.
Referenced by bookRoot(), and fillRoot().
int HIPAlignmentAlgorithm::m2_Nhit [private] |
Definition at line 124 of file HIPAlignmentAlgorithm.h.
Referenced by bookRoot(), and fillRoot().
Definition at line 127 of file HIPAlignmentAlgorithm.h.
Referenced by bookRoot(), and fillRoot().
float HIPAlignmentAlgorithm::m2_Phi [private] |
Definition at line 125 of file HIPAlignmentAlgorithm.h.
Referenced by bookRoot(), and fillRoot().
int HIPAlignmentAlgorithm::m2_Type [private] |
Definition at line 124 of file HIPAlignmentAlgorithm.h.
Referenced by bookRoot(), and fillRoot().
float HIPAlignmentAlgorithm::m2_Xpos [private] |
Definition at line 125 of file HIPAlignmentAlgorithm.h.
Referenced by bookRoot(), and fillRoot().
float HIPAlignmentAlgorithm::m2_Ypos [private] |
Definition at line 125 of file HIPAlignmentAlgorithm.h.
Referenced by bookRoot(), and fillRoot().
float HIPAlignmentAlgorithm::m2_Zpos [private] |
Definition at line 125 of file HIPAlignmentAlgorithm.h.
Referenced by bookRoot(), and fillRoot().
align::ID HIPAlignmentAlgorithm::m3_Id [private] |
Definition at line 130 of file HIPAlignmentAlgorithm.h.
Referenced by bookRoot(), and terminate().
Definition at line 131 of file HIPAlignmentAlgorithm.h.
Referenced by bookRoot(), and terminate().
float HIPAlignmentAlgorithm::m3_par[6] [private] |
Definition at line 132 of file HIPAlignmentAlgorithm.h.
Referenced by bookRoot(), and terminate().
float HIPAlignmentAlgorithm::m_Chi2n[MAXREC] [private] |
Definition at line 121 of file HIPAlignmentAlgorithm.h.
Referenced by bookRoot(), fillEventwiseTree(), and run().
float HIPAlignmentAlgorithm::m_d0[MAXREC] [private] |
Definition at line 121 of file HIPAlignmentAlgorithm.h.
Referenced by bookRoot(), fillEventwiseTree(), and run().
float HIPAlignmentAlgorithm::m_dz[MAXREC] [private] |
Definition at line 121 of file HIPAlignmentAlgorithm.h.
Referenced by bookRoot(), fillEventwiseTree(), and run().
float HIPAlignmentAlgorithm::m_Eta[MAXREC] [private] |
Definition at line 121 of file HIPAlignmentAlgorithm.h.
Referenced by bookRoot(), fillEventwiseTree(), and run().
int HIPAlignmentAlgorithm::m_Nhits[MAXREC] [private] |
Definition at line 120 of file HIPAlignmentAlgorithm.h.
Referenced by bookRoot(), fillEventwiseTree(), and run().
int HIPAlignmentAlgorithm::m_nhPXB[MAXREC] [private] |
Definition at line 120 of file HIPAlignmentAlgorithm.h.
Referenced by bookRoot(), fillEventwiseTree(), and run().
int HIPAlignmentAlgorithm::m_nhPXF[MAXREC] [private] |
Definition at line 120 of file HIPAlignmentAlgorithm.h.
Referenced by bookRoot(), fillEventwiseTree(), and run().
int HIPAlignmentAlgorithm::m_Ntracks [private] |
Definition at line 120 of file HIPAlignmentAlgorithm.h.
Referenced by bookRoot(), fillEventwiseTree(), and run().
float HIPAlignmentAlgorithm::m_P[MAXREC] [private] |
Definition at line 121 of file HIPAlignmentAlgorithm.h.
Referenced by bookRoot(), fillEventwiseTree(), and run().
float HIPAlignmentAlgorithm::m_Phi[MAXREC] [private] |
Definition at line 121 of file HIPAlignmentAlgorithm.h.
Referenced by bookRoot(), fillEventwiseTree(), and run().
float HIPAlignmentAlgorithm::m_Pt[MAXREC] [private] |
Definition at line 121 of file HIPAlignmentAlgorithm.h.
Referenced by bookRoot(), fillEventwiseTree(), and run().
const int HIPAlignmentAlgorithm::MAXREC = 99 [static, private] |
Definition at line 118 of file HIPAlignmentAlgorithm.h.
Referenced by fillEventwiseTree(), and run().
std::string HIPAlignmentAlgorithm::outfile [private] |
Definition at line 87 of file HIPAlignmentAlgorithm.h.
Referenced by bookRoot(), and HIPAlignmentAlgorithm().
std::string HIPAlignmentAlgorithm::outfile2 [private] |
Definition at line 87 of file HIPAlignmentAlgorithm.h.
Referenced by bookRoot(), and HIPAlignmentAlgorithm().
std::string HIPAlignmentAlgorithm::outpath [private] |
Definition at line 87 of file HIPAlignmentAlgorithm.h.
Referenced by HIPAlignmentAlgorithm().
std::string HIPAlignmentAlgorithm::salignedfile [private] |
Definition at line 88 of file HIPAlignmentAlgorithm.h.
Referenced by HIPAlignmentAlgorithm(), startNewLoop(), and terminate().
std::string HIPAlignmentAlgorithm::siterationfile [private] |
Definition at line 88 of file HIPAlignmentAlgorithm.h.
Referenced by HIPAlignmentAlgorithm(), startNewLoop(), and terminate().
std::string HIPAlignmentAlgorithm::smisalignedfile [private] |
Definition at line 88 of file HIPAlignmentAlgorithm.h.
Referenced by HIPAlignmentAlgorithm(), and startNewLoop().
std::string HIPAlignmentAlgorithm::sparameterfile [private] |
Definition at line 87 of file HIPAlignmentAlgorithm.h.
Referenced by HIPAlignmentAlgorithm(), and terminate().
std::string HIPAlignmentAlgorithm::ssurveyfile [private] |
Definition at line 88 of file HIPAlignmentAlgorithm.h.
Referenced by bookRoot(), and HIPAlignmentAlgorithm().
std::string HIPAlignmentAlgorithm::struefile [private] |
Definition at line 88 of file HIPAlignmentAlgorithm.h.
Referenced by HIPAlignmentAlgorithm(), and startNewLoop().
std::string HIPAlignmentAlgorithm::suvarfile [private] |
Definition at line 87 of file HIPAlignmentAlgorithm.h.
Referenced by HIPAlignmentAlgorithm(), and terminate().
Definition at line 76 of file HIPAlignmentAlgorithm.h.
Referenced by initialize(), and run().
std::vector<Alignable*> HIPAlignmentAlgorithm::theAlignables [private] |
Definition at line 75 of file HIPAlignmentAlgorithm.h.
Referenced by collector(), fillRoot(), initialize(), startNewLoop(), and terminate().
Definition at line 74 of file HIPAlignmentAlgorithm.h.
Referenced by collector(), fillRoot(), initialize(), run(), setAlignmentPositionError(), startNewLoop(), and terminate().
std::vector<std::pair<std::vector<Alignable*>, std::vector<double> > > HIPAlignmentAlgorithm::theAPEParameters [private] |
Definition at line 93 of file HIPAlignmentAlgorithm.h.
Referenced by initialize(), and setAlignmentPositionError().
std::vector<edm::ParameterSet> HIPAlignmentAlgorithm::theAPEParameterSet [private] |
Definition at line 92 of file HIPAlignmentAlgorithm.h.
Referenced by HIPAlignmentAlgorithm(), and initialize().
bool HIPAlignmentAlgorithm::theApplyAPE [private] |
Definition at line 91 of file HIPAlignmentAlgorithm.h.
Referenced by HIPAlignmentAlgorithm(), initialize(), and setAlignmentPositionError().
int HIPAlignmentAlgorithm::theCollectorNJobs [private] |
Definition at line 102 of file HIPAlignmentAlgorithm.h.
Referenced by collector(), and HIPAlignmentAlgorithm().
std::string HIPAlignmentAlgorithm::theCollectorPath [private] |
Definition at line 103 of file HIPAlignmentAlgorithm.h.
Referenced by collector(), and HIPAlignmentAlgorithm().
int HIPAlignmentAlgorithm::theCurrentPrescale [private] |
Definition at line 104 of file HIPAlignmentAlgorithm.h.
Referenced by HIPAlignmentAlgorithm(), and run().
int HIPAlignmentAlgorithm::theEventPrescale [private] |
Definition at line 104 of file HIPAlignmentAlgorithm.h.
Referenced by HIPAlignmentAlgorithm(), and run().
TFile* HIPAlignmentAlgorithm::theFile [private] |
Definition at line 110 of file HIPAlignmentAlgorithm.h.
Referenced by bookRoot(), run(), and terminate().
TFile* HIPAlignmentAlgorithm::theFile2 [private] |
Definition at line 112 of file HIPAlignmentAlgorithm.h.
Referenced by bookRoot(), fillRoot(), and terminate().
TFile* HIPAlignmentAlgorithm::theFile3 [private] |
Definition at line 114 of file HIPAlignmentAlgorithm.h.
Referenced by bookRoot(), and terminate().
bool HIPAlignmentAlgorithm::theFillTrackMonitoring [private] |
Definition at line 105 of file HIPAlignmentAlgorithm.h.
Referenced by collector(), HIPAlignmentAlgorithm(), and run().
AlignmentIORoot HIPAlignmentAlgorithm::theIO [private] |
Definition at line 78 of file HIPAlignmentAlgorithm.h.
Referenced by startNewLoop(), and terminate().
int HIPAlignmentAlgorithm::theIteration [private] |
Definition at line 80 of file HIPAlignmentAlgorithm.h.
Referenced by bookRoot(), collector(), setAlignmentPositionError(), startNewLoop(), and terminate().
std::vector<align::StructureType> HIPAlignmentAlgorithm::theLevels [private] |
Definition at line 107 of file HIPAlignmentAlgorithm.h.
Referenced by bookRoot(), HIPAlignmentAlgorithm(), and terminate().
double HIPAlignmentAlgorithm::theMaxAllowedHitPull [private] |
Definition at line 95 of file HIPAlignmentAlgorithm.h.
Referenced by HIPAlignmentAlgorithm(), processHit1D(), and processHit2D().
double HIPAlignmentAlgorithm::theMaxRelParameterError [private] |
Definition at line 99 of file HIPAlignmentAlgorithm.h.
Referenced by calcParameters(), and HIPAlignmentAlgorithm().
int HIPAlignmentAlgorithm::theMinimumNumberOfHits [private] |
Definition at line 97 of file HIPAlignmentAlgorithm.h.
Referenced by calcParameters(), and HIPAlignmentAlgorithm().
TTree* HIPAlignmentAlgorithm::theTree [private] |
Definition at line 111 of file HIPAlignmentAlgorithm.h.
Referenced by bookRoot(), fillEventwiseTree(), run(), and terminate().
TTree* HIPAlignmentAlgorithm::theTree2 [private] |
Definition at line 113 of file HIPAlignmentAlgorithm.h.
Referenced by bookRoot(), fillRoot(), and terminate().
TTree* HIPAlignmentAlgorithm::theTree3 [private] |
Definition at line 115 of file HIPAlignmentAlgorithm.h.
Referenced by bookRoot(), and terminate().
bool HIPAlignmentAlgorithm::verbose [private] |
Definition at line 85 of file HIPAlignmentAlgorithm.h.