File indexing completed on 2024-04-06 11:57:24
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012 #include "TFile.h"
0013 #include "TTree.h"
0014 #include "Math/EulerAngles.h"
0015
0016 #include "CondFormats/Alignment/interface/Alignments.h"
0017 #include "CondFormats/Alignment/interface/SurveyErrors.h"
0018 #include "CondFormats/AlignmentRecord/interface/TrackerSurveyRcd.h"
0019 #include "CondFormats/AlignmentRecord/interface/TrackerSurveyErrorExtendedRcd.h"
0020 #include "FWCore/Framework/interface/ESHandle.h"
0021 #include "FWCore/Framework/interface/EventSetup.h"
0022 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0023 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0024 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0025
0026 typedef AlignTransform SurveyValue;
0027 typedef Alignments SurveyValues;
0028
0029 class SurveyDBReader : public edm::one::EDAnalyzer<> {
0030 public:
0031
0032 SurveyDBReader(const edm::ParameterSet&);
0033
0034
0035 virtual void beginJob() { theFirstEvent = true; }
0036
0037 virtual void analyze(const edm::Event&, const edm::EventSetup&);
0038
0039 private:
0040 const edm::ESGetToken<SurveyValues, TrackerSurveyRcd> tkSurveyToken_;
0041 const edm::ESGetToken<SurveyErrors, TrackerSurveyErrorExtendedRcd> tkSurveyErrToken_;
0042 const std::string theFileName;
0043
0044 bool theFirstEvent;
0045 };
0046
0047 SurveyDBReader::SurveyDBReader(const edm::ParameterSet& cfg)
0048 : tkSurveyToken_(esConsumes()),
0049 tkSurveyErrToken_(esConsumes()),
0050 theFileName(cfg.getParameter<std::string>("fileName")),
0051 theFirstEvent(true) {}
0052
0053 void SurveyDBReader::analyze(const edm::Event&, const edm::EventSetup& setup) {
0054 if (theFirstEvent) {
0055 const SurveyValues* valuesHandle = &setup.getData(tkSurveyToken_);
0056 const SurveyErrors* errorsHandle = &setup.getData(tkSurveyErrToken_);
0057
0058 const std::vector<SurveyValue>& values = valuesHandle->m_align;
0059 const std::vector<SurveyError>& errors = errorsHandle->m_surveyErrors;
0060
0061 unsigned int size = values.size();
0062
0063 if (errors.size() != size) {
0064 edm::LogError("SurveyDBReader") << "Value and error records have different sizes. Abort!";
0065
0066 return;
0067 }
0068
0069 uint8_t type = 0;
0070 uint32_t id = 0;
0071
0072 ROOT::Math::Cartesian3D<double>* pos(0);
0073 ROOT::Math::EulerAngles* rot(0);
0074 const align::ErrorMatrix* cov(0);
0075
0076 TFile fout(theFileName.c_str(), "RECREATE");
0077 TTree tree("survey", "");
0078
0079 tree.Branch("type", &type, "type/b");
0080 tree.Branch("id", &id, "id/i");
0081 tree.Branch("pos", "ROOT::Math::Cartesian3D<double>", &pos);
0082 tree.Branch("rot", "ROOT::Math::EulerAngles", &rot);
0083 tree.Branch("cov", "align::ErrorMatrix", &cov);
0084
0085 for (unsigned int i = 0; i < size; ++i) {
0086 const SurveyValue& value = values[i];
0087 const SurveyError& error = errors[i];
0088
0089 const CLHEP::Hep3Vector& transl = value.translation();
0090 const CLHEP::HepRotation& orient = value.rotation();
0091 const align::ErrorMatrix& matrix = error.matrix();
0092
0093 ROOT::Math::Cartesian3D<double> coords(transl.x(), transl.y(), transl.z());
0094 ROOT::Math::EulerAngles angles(orient.phi(), orient.theta(), orient.psi());
0095
0096 type = error.structureType();
0097 id = error.rawId();
0098 pos = &coords;
0099 rot = &angles;
0100 cov = &matrix;
0101
0102 tree.Fill();
0103 }
0104
0105 fout.Write();
0106
0107 edm::LogInfo("SurveyDBReader") << "Number of alignables read " << size << std::endl;
0108
0109 theFirstEvent = false;
0110 }
0111 }
0112
0113 #include "FWCore/Framework/interface/MakerMacros.h"
0114 DEFINE_FWK_MODULE(SurveyDBReader);