CMS 3D CMS Logo

AlignmentCorrelationsIORoot.cc
Go to the documentation of this file.
1 #include "TTree.h"
2 
6 
7 // this class's header
9 
10 // ----------------------------------------------------------------------------
11 // constructor
12 
14  treename = "AlignmentCorrelations";
15  treetxt = "Correlations";
16 }
17 
18 // ----------------------------------------------------------------------------
19 
21  tree->Branch("Ali1Id", &Ali1Id, "Ali1Id/i");
22  tree->Branch("Ali2Id", &Ali2Id, "Ali2Id/i");
23  tree->Branch("Ali1ObjId", &Ali1ObjId, "Ali1ObjId/I");
24  tree->Branch("Ali2ObjId", &Ali2ObjId, "Ali2ObjId/I");
25  tree->Branch("corSize", &corSize, "corSize/I");
26  tree->Branch("CorMatrix", &CorMatrix, "CorMatrix[corSize]/D");
27 }
28 
29 // ----------------------------------------------------------------------------
30 
32  tree->SetBranchAddress("corSize", &corSize);
33  tree->SetBranchAddress("Ali1Id", &Ali1Id);
34  tree->SetBranchAddress("Ali2Id", &Ali2Id);
35  tree->SetBranchAddress("Ali1ObjId", &Ali1ObjId);
36  tree->SetBranchAddress("Ali2ObjId", &Ali2ObjId);
37  tree->SetBranchAddress("CorMatrix", &CorMatrix);
38 }
39 
40 // ----------------------------------------------------------------------------
41 
42 int AlignmentCorrelationsIORoot::write(const align::Correlations& cor, bool validCheck) {
43  int icount = 0;
44 
45  for (align::Correlations::const_iterator it = cor.begin(); it != cor.end(); ++it) {
46  AlgebraicMatrix mat = (*it).second;
47  std::pair<Alignable*, Alignable*> Pair = (*it).first;
48  Alignable* ali1 = Pair.first;
49  Alignable* ali2 = Pair.second;
50  if ((ali1->alignmentParameters()->isValid() && ali2->alignmentParameters()->isValid()) || !(validCheck)) {
51  Ali1ObjId = ali1->alignableObjectId();
52  Ali2ObjId = ali2->alignableObjectId();
53  Ali1Id = ali1->id();
54  Ali2Id = ali2->id();
55  int maxColumn = mat.num_row();
56  corSize = maxColumn * maxColumn;
57  for (int row = 0; row < maxColumn; row++)
58  for (int col = 0; col < maxColumn; col++)
59  CorMatrix[row + col * maxColumn] = mat[row][col];
60  tree->Fill();
61  icount++;
62  }
63  }
64  edm::LogInfo("AlignmentCorrelationsIORoot") << "Writing correlations: all,written: " << cor.size() << "," << icount;
65  return 0;
66 }
67 
68 // ----------------------------------------------------------------------------
69 // read correlations for those alignables in vector given as argument
70 
72  align::Correlations theMap;
73 
74  // create ID map for all Alignables in alivec
75  align::Alignables::const_iterator it1;
76  std::map<std::pair<unsigned int, int>, Alignable*> idAlis;
77  for (it1 = alivec.begin(); it1 != alivec.end(); ++it1)
78  idAlis[std::make_pair((*it1)->id(), (*it1)->alignableObjectId())] = (*it1);
79 
80  std::map<std::pair<unsigned int, int>, Alignable*>::const_iterator aliSearch1;
81  std::map<std::pair<unsigned int, int>, Alignable*>::const_iterator aliSearch2;
82  int nfound = 0;
83  double maxEntry = tree->GetEntries();
84  for (int entry = 0; entry < maxEntry; entry++) {
85  tree->GetEntry(entry);
86  aliSearch1 = idAlis.find(std::make_pair(Ali1Id, Ali1ObjId));
87  aliSearch2 = idAlis.find(std::make_pair(Ali2Id, Ali2ObjId));
88  if (aliSearch1 != idAlis.end() && aliSearch2 != idAlis.end()) {
89  // Alignables for this pair found
90  nfound++;
91  Alignable* myAli1 = (*aliSearch1).second;
92  Alignable* myAli2 = (*aliSearch2).second;
93  // FIXME: instead of nParMax in the next few lines one should probably
94  // use something like sqrt(corSize) - but take care of rounding!
95  // I have no time to test... :-( GF
97  for (int row = 0; row < nParMax; row++)
98  for (int col = 0; col < nParMax; col++)
99  mat[row][col] = CorMatrix[row + col * nParMax];
100  theMap[std::make_pair(myAli1, myAli2)] = mat;
101  }
102  }
103 
104  edm::LogInfo("AlignmentCorrelationsIORoot") << "Read correlations: all,read: " << alivec.size() << "," << nfound;
105 
106  ierr = 0;
107  return theMap;
108 }
align::ID Ali1Id
correlation tree
align::ID id() const
Return the ID of Alignable, i.e. DetId of &#39;first&#39; component GeomDet(Unit).
Definition: Alignable.h:189
static const int nParMax
void setBranchAddresses(void) override
set root branches
std::map< std::pair< Alignable *, Alignable * >, AlgebraicMatrix > Correlations
Definition: Utilities.h:36
AlignmentParameters * alignmentParameters() const
Get the AlignmentParameters.
Definition: Alignable.h:61
virtual StructureType alignableObjectId() const =0
Return the alignable type identifier.
CLHEP::HepMatrix AlgebraicMatrix
int write(const align::Correlations &cor, bool validCheck) override
write correlations
void createBranches(void) override
create root branches
std::vector< Alignable * > Alignables
Definition: Utilities.h:32
align::Correlations read(const align::Alignables &alivec, int &ierr) override
read correlations
bool isValid(void) const
Get validity flag.
col
Definition: cuy.py:1010
Definition: tree.py:1