File indexing completed on 2024-04-06 11:56:01
0001 #include "Alignment/CocoaModel/interface/CocoaDaqReaderRoot.h"
0002 #include "TFile.h"
0003 #include "Alignment/CocoaDaq/interface/CocoaDaqRootEvent.h"
0004 #include "Alignment/CocoaModel/interface/Measurement.h"
0005 #include "Alignment/CocoaModel/interface/Model.h"
0006 #include "Alignment/CocoaUtilities/interface/ALIUtils.h"
0007
0008 #include "CondFormats/OptAlignObjects/interface/OpticalAlignMeasurements.h"
0009
0010 #include <iostream>
0011 #include <string>
0012
0013 #include "TClonesArray.h"
0014
0015
0016 CocoaDaqReaderRoot::CocoaDaqReaderRoot(const std::string& m_inFileName) {
0017 if (ALIUtils::debug >= 3)
0018 std::cout << " CocoaDaqReaderRoot opening file: " << m_inFileName << std::endl;
0019
0020 theFile = new TFile(m_inFileName.c_str());
0021 if (!theTree) {
0022 std::cerr << " CocoaDaqReaderRoot TTree file not found " << m_inFileName << std::endl;
0023 throw std::exception();
0024 }
0025
0026
0027 theTree = (TTree*)theFile->Get("CocoaDaq");
0028
0029
0030 if (!theTree) {
0031 std::cerr << " CocoaDaqReaderRoot TTree in file " << m_inFileName << " should be called 'CocoaDaq' " << std::endl;
0032 throw std::exception();
0033 }
0034 TBranch* branch = theTree->GetBranch("Alignment_Cocoa");
0035
0036 nev = branch->GetEntries();
0037
0038
0039 nextEvent = 0;
0040
0041
0042 theEvent = new CocoaDaqRootEvent();
0043
0044
0045 theTree->SetBranchAddress("Alignment_Cocoa", &theEvent);
0046
0047
0048 CocoaDaqReader::SetDaqReader(this);
0049 }
0050
0051
0052 CocoaDaqReaderRoot::~CocoaDaqReaderRoot() { theFile->Close(); }
0053
0054
0055 bool CocoaDaqReaderRoot::ReadNextEvent() { return ReadEvent(nextEvent); }
0056
0057
0058 bool CocoaDaqReaderRoot::ReadEvent(int nev) {
0059 std::vector<OpticalAlignMeasurementInfo> measList;
0060
0061 int nb = 0;
0062
0063 nb = theTree->GetEntry(nev);
0064
0065 if (ALIUtils::debug >= 3)
0066 std::cout << "CocoaDaqReaderRoot reading event " << nev << " " << nb << std::endl;
0067 if (nb == 0)
0068 return false;
0069
0070
0071 int n = 1;
0072 if (nev % n == 0 && ALIUtils::debug >= 3)
0073 theEvent->DumpIt();
0074
0075
0076
0077 if (ALIUtils::debug >= 3)
0078 std::cout << " CocoaDaqReaderRoot::ReadEvent npos2D " << theEvent->GetNumPos2D() << " nCOPS "
0079 << theEvent->GetNumPosCOPS() << std::endl;
0080
0081 for (int ii = 0; ii < theEvent->GetNumPos2D(); ii++) {
0082 AliDaqPosition2D* pos2D = (AliDaqPosition2D*)theEvent->GetArray_Position2D()->At(ii);
0083 if (ALIUtils::debug >= 4)
0084 std::cout << "2D sensor " << ii << " has ID = " << pos2D->GetID() << std::endl;
0085 pos2D->DumpIt("2DSENSOR");
0086 measList.push_back(GetMeasFromPosition2D(pos2D));
0087 }
0088 for (int ii = 0; ii < theEvent->GetNumPosCOPS(); ii++) {
0089 AliDaqPositionCOPS* posCOPS = (AliDaqPositionCOPS*)theEvent->GetArray_PositionCOPS()->At(ii);
0090 measList.push_back(GetMeasFromPositionCOPS(posCOPS));
0091 if (ALIUtils::debug >= 4) {
0092 std::cout << "COPS sensor " << ii << " has ID = " << posCOPS->GetID() << std::endl;
0093 posCOPS->DumpIt("COPS");
0094 }
0095 }
0096 for (int ii = 0; ii < theEvent->GetNumTilt(); ii++) {
0097 AliDaqTilt* tilt = (AliDaqTilt*)theEvent->GetArray_Tilt()->At(ii);
0098 measList.push_back(GetMeasFromTilt(tilt));
0099 if (ALIUtils::debug >= 4) {
0100 std::cout << "TILT sensor " << ii << " has ID = " << tilt->GetID() << std::endl;
0101 tilt->DumpIt("TILT");
0102 }
0103 }
0104 for (int ii = 0; ii < theEvent->GetNumDist(); ii++) {
0105 AliDaqDistance* dist = (AliDaqDistance*)theEvent->GetArray_Dist()->At(ii);
0106 measList.push_back(GetMeasFromDist(dist));
0107 if (ALIUtils::debug >= 4) {
0108 std::cout << "DIST sensor " << ii << " has ID = " << dist->GetID() << std::endl;
0109 dist->DumpIt("DIST");
0110 }
0111 }
0112
0113 nextEvent = nev + 1;
0114
0115 BuildMeasurementsFromOptAlign(measList);
0116
0117 return true;
0118 }
0119
0120
0121 OpticalAlignMeasurementInfo CocoaDaqReaderRoot::GetMeasFromPosition2D(AliDaqPosition2D* pos2D) {
0122 OpticalAlignMeasurementInfo meas;
0123
0124 meas.type_ = "SENSOR2D";
0125 meas.name_ = std::string(pos2D->GetID().Data());
0126
0127 std::vector<bool> isSimu;
0128 for (size_t jj = 0; jj < 2; jj++) {
0129 isSimu.push_back(false);
0130 }
0131 meas.isSimulatedValue_ = isSimu;
0132 std::vector<OpticalAlignParam> paramList;
0133 OpticalAlignParam oaParam1;
0134 oaParam1.name_ = "H:";
0135 oaParam1.value_ = pos2D->GetX() / 100.;
0136 oaParam1.error_ = pos2D->GetXerror() / 100.;
0137 paramList.push_back(oaParam1);
0138
0139 OpticalAlignParam oaParam2;
0140 oaParam2.name_ = "V:";
0141 oaParam2.value_ = pos2D->GetY() / 100.;
0142 oaParam2.error_ = pos2D->GetYerror() / 100.;
0143 paramList.push_back(oaParam2);
0144
0145 meas.values_ = paramList;
0146
0147 return meas;
0148 }
0149
0150
0151 OpticalAlignMeasurementInfo CocoaDaqReaderRoot::GetMeasFromPositionCOPS(AliDaqPositionCOPS* posCOPS) {
0152 OpticalAlignMeasurementInfo meas;
0153
0154 meas.type_ = "COPS";
0155 meas.name_ = std::string(posCOPS->GetID().Data());
0156
0157 std::vector<bool> isSimu;
0158 for (size_t jj = 0; jj < 4; jj++) {
0159 isSimu.push_back(false);
0160 }
0161 meas.isSimulatedValue_ = isSimu;
0162
0163 std::vector<OpticalAlignParam> paramList;
0164 OpticalAlignParam oaParam1;
0165 oaParam1.name_ = "U:";
0166 oaParam1.value_ = posCOPS->GetUp() / 100.;
0167 oaParam1.error_ = posCOPS->GetUpError() / 100.;
0168 paramList.push_back(oaParam1);
0169
0170 OpticalAlignParam oaParam2;
0171 oaParam2.name_ = "U:";
0172 oaParam2.value_ = posCOPS->GetDown() / 100.;
0173 oaParam2.error_ = posCOPS->GetDownError() / 100.;
0174 paramList.push_back(oaParam2);
0175
0176 OpticalAlignParam oaParam3;
0177 oaParam3.name_ = "U:";
0178 oaParam3.value_ = posCOPS->GetRight() / 100.;
0179 oaParam3.error_ = posCOPS->GetRightError() / 100.;
0180 paramList.push_back(oaParam3);
0181
0182 OpticalAlignParam oaParam4;
0183 oaParam4.name_ = "U:";
0184 oaParam4.value_ = posCOPS->GetLeft() / 100.;
0185 oaParam4.error_ = posCOPS->GetLeftError() / 100.;
0186 paramList.push_back(oaParam4);
0187
0188 meas.values_ = paramList;
0189
0190 return meas;
0191 }
0192
0193
0194 OpticalAlignMeasurementInfo CocoaDaqReaderRoot::GetMeasFromTilt(AliDaqTilt* tilt) {
0195 OpticalAlignMeasurementInfo meas;
0196
0197 meas.type_ = "TILTMETER";
0198 meas.name_ = std::string(tilt->GetID().Data());
0199
0200 std::vector<bool> isSimu;
0201 for (size_t jj = 0; jj < 2; jj++) {
0202 isSimu.push_back(false);
0203 }
0204 meas.isSimulatedValue_ = isSimu;
0205 std::vector<OpticalAlignParam> paramList;
0206 OpticalAlignParam oaParam;
0207 oaParam.name_ = "T:";
0208 oaParam.value_ = tilt->GetTilt();
0209 oaParam.error_ = tilt->GetTiltError();
0210 paramList.push_back(oaParam);
0211
0212 meas.values_ = paramList;
0213
0214 return meas;
0215 }
0216
0217
0218 OpticalAlignMeasurementInfo CocoaDaqReaderRoot::GetMeasFromDist(AliDaqDistance* dist) {
0219 OpticalAlignMeasurementInfo meas;
0220
0221 meas.type_ = "DISTANCEMETER";
0222 meas.name_ = std::string(dist->GetID().Data());
0223
0224 std::vector<bool> isSimu;
0225 for (size_t jj = 0; jj < 2; jj++) {
0226 isSimu.push_back(false);
0227 }
0228 meas.isSimulatedValue_ = isSimu;
0229 std::vector<OpticalAlignParam> paramList;
0230 OpticalAlignParam oaParam;
0231 oaParam.name_ = "D:";
0232 oaParam.value_ = dist->GetDistance() / 100.;
0233 oaParam.error_ = dist->GetDistanceError() / 100.;
0234 paramList.push_back(oaParam);
0235
0236 meas.values_ = paramList;
0237
0238 return meas;
0239 }
0240
0241
0242 void CocoaDaqReaderRoot::BuildMeasurementsFromOptAlign(std::vector<OpticalAlignMeasurementInfo>& measList) {
0243 if (ALIUtils::debug >= 3)
0244 std::cout << "@@@ CocoaDaqReaderRoot::BuildMeasurementsFromOptAlign " << std::endl;
0245
0246
0247
0248
0249
0250
0251
0252
0253 ALIint nMeasRoot = measList.size();
0254 if (ALIUtils::debug >= 4) {
0255 std::cout << " Building " << nMeasRoot << " measurements from ROOT file " << std::endl;
0256 }
0257
0258
0259 std::vector<Measurement*>::const_iterator vmcite;
0260 for (vmcite = Model::MeasurementList().begin(); vmcite != Model::MeasurementList().end(); ++vmcite) {
0261 ALIint fcolon = (*vmcite)->name().find(':');
0262 ALIstring oname = (*vmcite)->name();
0263 oname = oname.substr(fcolon + 1, oname.length());
0264
0265
0266 ALIint ii;
0267 for (ii = 0; ii < nMeasRoot; ii++) {
0268 OpticalAlignMeasurementInfo measInfo = measList[ii];
0269 std::cout << " measurement name ROOT " << measInfo.name_ << " Model= " << (*vmcite)->name() << " short " << oname
0270 << std::endl;
0271
0272 if (oname == measInfo.name_) {
0273
0274
0275 if ((*vmcite)->type() != measInfo.type_) {
0276 std::cerr << "!!! Measurement from ROOT file: type in file is " << measInfo.type_ << " and should be "
0277 << (*vmcite)->type() << std::endl;
0278 exit(1);
0279 }
0280
0281 std::cout << " NOBJECTS IN MEAS " << (*vmcite)->OptOList().size() << " NMEAS "
0282 << Model::MeasurementList().size() << std::endl;
0283
0284 std::vector<OpticalAlignParam> measValues = measInfo.values_;
0285
0286 for (size_t jj = 0; jj < measValues.size(); jj++) {
0287 (*vmcite)->fillData(jj, &(measValues[jj]));
0288 }
0289
0290 std::cout << " NOBJECTS IN MEAS after " << (*vmcite)->OptOList().size() << " NMEAS "
0291 << Model::MeasurementList().size() << std::endl;
0292
0293 break;
0294 }
0295 }
0296 if (ii == nMeasRoot) {
0297 std::cerr << "!!! Reading measurement from file: measurement not found! Type in list is " << oname << std::endl;
0298 exit(1);
0299 }
0300 }
0301 }