Merge pull request #36257 from CTPPS/pps-alignment-global

PPS global alignment - reformat

CalibPPS/AlignmentGlobal/BuildFile.xml

@@ -0,0 +1,5 @@
+<export>
+  <lib name="1"/>
+</export>
+<use name="rootcore"/>
+<use name="rootgraphics"/>

CalibPPS/AlignmentGlobal/interface/utils.h

@@ -0,0 +1,25 @@
+/****************************************************************************
+* Authors: 
+*  Jan Kašpar ([email protected]) 
+*  Mateusz Kocot ([email protected])
+****************************************************************************/
+
+#ifndef CalibPPS_AlignmentGlobal_utils_h
+#define CalibPPS_AlignmentGlobal_utils_h
+
+#include "TProfile.h"
+
+namespace alig_utils {
+
+  // Fits a linear function to a TProfile.
+  int fitProfile(TProfile* p,
+                 double x_mean,
+                 double x_rms,
+                 unsigned int minBinEntries,
+                 unsigned int minNBinsReasonable,
+                 double& sl,
+                 double& sl_unc);
+
+}  // namespace alig_utils
+
+#endif

CalibPPS/AlignmentGlobal/plugins/BuildFile.xml

@@ -10,6 +10,7 @@
   <use name="DataFormats/CTPPSReco"/>
   <use name="CondFormats/PPSObjects"/>
   <use name="CondFormats/DataRecord"/>
+  <use name="CalibPPS/AlignmentGlobal"/>
   <use name="rootcore"/>
   <use name="rootgraphics"/>
 </library>

CalibPPS/AlignmentGlobal/plugins/PPSAlignmentWorker.cc

@@ -48,9 +48,7 @@ class PPSAlignmentWorker : public DQMEDAnalyzer {
 
     // hit distributions
     std::map<unsigned int, MonitorElement*> m_h2_y_vs_x_bef_sel;
-
     std::map<unsigned int, MonitorElement*> m_h2_y_vs_x_mlt_sel;
-
     std::map<unsigned int, MonitorElement*> m_h2_y_vs_x_aft_sel;
 
     // cut plots
@@ -75,14 +73,15 @@ class PPSAlignmentWorker : public DQMEDAnalyzer {
 
       SlicePlots();
       SlicePlots(DQMStore::IBooker& iBooker, const PPSAlignmentConfiguration& cfg, bool debug);
+      void fill(const double y, const double yDiff, const bool debug);
     };
 
     std::map<unsigned int, SlicePlots> x_slice_plots_N, x_slice_plots_F;
 
     void init(DQMStore::IBooker& iBooker,
               const PPSAlignmentConfiguration& cfg,
               const PPSAlignmentConfiguration::SectorConfig& scfg,
-              const std::string& folder,
+              const std::string& rootDir,
               bool debug);
 
     unsigned int process(const CTPPSLocalTrackLiteCollection& tracks, const PPSAlignmentConfiguration& cfg, bool debug);
@@ -97,7 +96,7 @@ class PPSAlignmentWorker : public DQMEDAnalyzer {
   SectorData sectorData45_;
   SectorData sectorData56_;
 
-  std::string folder_;
+  std::string dqmDir_;
   bool debug_;
 };
 
@@ -110,22 +109,22 @@ PPSAlignmentWorker::PPSAlignmentWorker(const edm::ParameterSet& iConfig)
       esTokenAnalyze_(esConsumes<PPSAlignmentConfiguration, PPSAlignmentConfigurationRcd>(
           edm::ESInputTag("", iConfig.getParameter<std::string>("label")))),
       tracksToken_(consumes<CTPPSLocalTrackLiteCollection>(iConfig.getParameter<edm::InputTag>("tagTracks"))),
-      folder_(iConfig.getParameter<std::string>("folder")),
+      dqmDir_(iConfig.getParameter<std::string>("dqm_dir")),
       debug_(iConfig.getParameter<bool>("debug")) {
-  edm::LogInfo("PPS").log([&](auto& li) {
-    li << "[worker] parameters:\n";
+  edm::LogInfo("PPSAlignmentWorker").log([&](auto& li) {
+    li << "parameters:\n";
     li << "* label: " << iConfig.getParameter<std::string>("label") << "\n";
     li << "* tagTracks: " << iConfig.getParameter<edm::InputTag>("tagTracks") << "\n";
-    li << "* folder: " << folder_ << "\n";
+    li << "* dqm_dir: " << dqmDir_ << "\n";
     li << "* debug: " << std::boolalpha << debug_;
   });
 }
 
 void PPSAlignmentWorker::bookHistograms(DQMStore::IBooker& iBooker, edm::Run const&, edm::EventSetup const& iSetup) {
   const auto& cfg = iSetup.getData(esTokenBookHistograms_);
 
-  sectorData45_.init(iBooker, cfg, cfg.sectorConfig45(), folder_ + "/worker", debug_);
-  sectorData56_.init(iBooker, cfg, cfg.sectorConfig56(), folder_ + "/worker", debug_);
+  sectorData45_.init(iBooker, cfg, cfg.sectorConfig45(), dqmDir_ + "/worker", debug_);
+  sectorData56_.init(iBooker, cfg, cfg.sectorConfig56(), dqmDir_ + "/worker", debug_);
 }
 
 void PPSAlignmentWorker::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
@@ -142,7 +141,7 @@ void PPSAlignmentWorker::fillDescriptions(edm::ConfigurationDescriptions& descri
 
   desc.add<std::string>("label", "");
   desc.add<edm::InputTag>("tagTracks", edm::InputTag("ctppsLocalTrackLiteProducer"));
-  desc.add<std::string>("folder", "AlCaReco/PPSAlignment");
+  desc.add<std::string>("dqm_dir", "AlCaReco/PPSAlignment");
   desc.add<bool>("debug", false);
 
   descriptions.addWithDefaultLabel(desc);
@@ -173,10 +172,17 @@ PPSAlignmentWorker::SectorData::SlicePlots::SlicePlots(DQMStore::IBooker& iBooke
                                        cfg.binning().slice_y_max_);
 }
 
+void PPSAlignmentWorker::SectorData::SlicePlots::fill(const double y, const double yDiff, const bool debug) {
+  h_y->Fill(y);
+  p_y_diffFN_vs_y->Fill(y, yDiff);
+  if (debug)
+    h2_y_diffFN_vs_y->Fill(y, yDiff);
+}
+
 void PPSAlignmentWorker::SectorData::init(DQMStore::IBooker& iBooker,
                                           const PPSAlignmentConfiguration& cfg,
                                           const PPSAlignmentConfiguration::SectorConfig& scfg,
-                                          const std::string& folder,
+                                          const std::string& rootDir,
                                           bool debug) {
   scfg_ = scfg;
 
@@ -193,44 +199,44 @@ void PPSAlignmentWorker::SectorData::init(DQMStore::IBooker& iBooker,
   const double y_min = cfg.binning().y_min_, y_max = cfg.binning().y_max_;
 
   // hit distributions
-  iBooker.setCurrentFolder(folder + "/" + scfg_.name_ + "/before selection/" + scfg_.rp_N_.name_);
+  iBooker.setCurrentFolder(rootDir + "/" + scfg_.name_ + "/before selection/" + scfg_.rp_N_.name_);
   m_h2_y_vs_x_bef_sel[scfg_.rp_N_.id_] =
       iBooker.book2DD("h2_y_vs_x", ";x;y", n_bins_x, x_min_str, x_max_str, n_bins_y, y_min, y_max);
-  iBooker.setCurrentFolder(folder + "/" + scfg_.name_ + "/before selection/" + scfg_.rp_F_.name_);
+  iBooker.setCurrentFolder(rootDir + "/" + scfg_.name_ + "/before selection/" + scfg_.rp_F_.name_);
   m_h2_y_vs_x_bef_sel[scfg_.rp_F_.id_] =
       iBooker.book2DD("h2_y_vs_x", ";x;y", n_bins_x, x_min_pix, x_max_pix, n_bins_y, y_min, y_max);
 
-  iBooker.setCurrentFolder(folder + "/" + scfg_.name_ + "/multiplicity selection/" + scfg_.rp_N_.name_);
+  iBooker.setCurrentFolder(rootDir + "/" + scfg_.name_ + "/multiplicity selection/" + scfg_.rp_N_.name_);
   m_h2_y_vs_x_mlt_sel[scfg_.rp_N_.id_] =
       iBooker.book2DD("h2_y_vs_x", ";x;y", n_bins_x, x_min_str, x_max_str, n_bins_y, y_min, y_max);
-  iBooker.setCurrentFolder(folder + "/" + scfg_.name_ + "/multiplicity selection/" + scfg_.rp_F_.name_);
+  iBooker.setCurrentFolder(rootDir + "/" + scfg_.name_ + "/multiplicity selection/" + scfg_.rp_F_.name_);
   m_h2_y_vs_x_mlt_sel[scfg_.rp_F_.id_] =
       iBooker.book2DD("h2_y_vs_x", ";x;y", n_bins_x, x_min_pix, x_max_pix, n_bins_y, y_min, y_max);
 
-  iBooker.setCurrentFolder(folder + "/" + scfg_.name_ + "/after selection/" + scfg_.rp_N_.name_);
+  iBooker.setCurrentFolder(rootDir + "/" + scfg_.name_ + "/after selection/" + scfg_.rp_N_.name_);
   m_h2_y_vs_x_aft_sel[scfg_.rp_N_.id_] =
       iBooker.book2DD("h2_y_vs_x", ";x;y", n_bins_x, x_min_str, x_max_str, n_bins_y, y_min, y_max);
-  iBooker.setCurrentFolder(folder + "/" + scfg_.name_ + "/after selection/" + scfg_.rp_F_.name_);
+  iBooker.setCurrentFolder(rootDir + "/" + scfg_.name_ + "/after selection/" + scfg_.rp_F_.name_);
   m_h2_y_vs_x_aft_sel[scfg_.rp_F_.id_] =
       iBooker.book2DD("h2_y_vs_x", ";x;y", n_bins_x, x_min_pix, x_max_pix, n_bins_y, y_min, y_max);
 
   // cut plots
-  iBooker.setCurrentFolder(folder + "/" + scfg_.name_ + "/cuts/cut_h");
+  iBooker.setCurrentFolder(rootDir + "/" + scfg_.name_ + "/cuts/cut_h");
   h_q_cut_h_bef = iBooker.book1DD("h_q_cut_h_bef", ";cq_h", 400, -2., 2.);
   h_q_cut_h_aft = iBooker.book1DD("h_q_cut_h_aft", ";cq_h", 400, -2., 2.);
   h2_cut_h_bef =
       iBooker.book2DD("h2_cut_h_bef", ";x_up;x_dw", n_bins_x, x_min_str, x_max_str, n_bins_x, x_min_pix, x_max_pix);
   h2_cut_h_aft =
       iBooker.book2DD("h2_cut_h_aft", ";x_up;x_dw", n_bins_x, x_min_str, x_max_str, n_bins_x, x_min_pix, x_max_pix);
 
-  iBooker.setCurrentFolder(folder + "/" + scfg_.name_ + "/cuts/cut_v");
+  iBooker.setCurrentFolder(rootDir + "/" + scfg_.name_ + "/cuts/cut_v");
   h_q_cut_v_bef = iBooker.book1DD("h_q_cut_v_bef", ";cq_v", 400, -2., 2.);
   h_q_cut_v_aft = iBooker.book1DD("h_q_cut_v_aft", ";cq_v", 400, -2., 2.);
   h2_cut_v_bef = iBooker.book2DD("h2_cut_v_bef", ";y_up;y_dw", n_bins_y, y_min, y_max, n_bins_y, y_min, y_max);
   h2_cut_v_aft = iBooker.book2DD("h2_cut_v_aft", ";y_up;y_dw", n_bins_y, y_min, y_max, n_bins_y, y_min, y_max);
 
   // near-far plots
-  iBooker.setCurrentFolder(folder + "/" + scfg_.name_ + "/near_far");
+  iBooker.setCurrentFolder(rootDir + "/" + scfg_.name_ + "/near_far");
 
   auto profilePtr = std::make_unique<TProfile>("",
                                                ";x_{N};x_{F} - x_{N}",
@@ -239,14 +245,15 @@ void PPSAlignmentWorker::SectorData::init(DQMStore::IBooker& iBooker,
                                                cfg.binning().diffFN_x_max_);
   p_x_diffFN_vs_x_N = iBooker.bookProfile("p_x_diffFN_vs_x_N", profilePtr.get());
 
+  // slice plots
   for (int i = 0; i < scfg_.rp_N_.x_slice_n_; i++) {
     const double x_min = scfg_.rp_N_.x_slice_min_ + i * scfg_.rp_N_.x_slice_w_;
     const double x_max = scfg_.rp_N_.x_slice_min_ + (i + 1) * scfg_.rp_N_.x_slice_w_;
 
     char buf[100];
     sprintf(buf, "%.1f-%.1f", x_min, x_max);
 
-    iBooker.setCurrentFolder(folder + "/" + scfg_.name_ + "/near_far/x slices, N/" + buf);
+    iBooker.setCurrentFolder(rootDir + "/" + scfg_.name_ + "/near_far/x slices, N/" + buf);
     x_slice_plots_N.insert({i, SlicePlots(iBooker, cfg, debug)});
   }
 
@@ -257,7 +264,7 @@ void PPSAlignmentWorker::SectorData::init(DQMStore::IBooker& iBooker,
     char buf[100];
     sprintf(buf, "%.1f-%.1f", x_min, x_max);
 
-    iBooker.setCurrentFolder(folder + "/" + scfg_.name_ + "/near_far/x slices, F/" + buf);
+    iBooker.setCurrentFolder(rootDir + "/" + scfg_.name_ + "/near_far/x slices, F/" + buf);
     x_slice_plots_F.insert({i, SlicePlots(iBooker, cfg, debug)});
   }
 }
@@ -304,16 +311,17 @@ unsigned int PPSAlignmentWorker::SectorData::process(const CTPPSLocalTrackLiteCo
 
   // do the selection
   unsigned int pairsSelected = 0;
-
   for (const auto& trUp : tracksUp) {
     for (const auto& trDw : tracksDw) {
       h2_cut_h_bef->Fill(trUp.x(), trDw.x());
       h2_cut_v_bef->Fill(trUp.y(), trDw.y());
 
+      // horizontal cut
       const double cq_h = trDw.x() + scfg_.cut_h_a_ * trUp.x() + scfg_.cut_h_c_;
       h_q_cut_h_bef->Fill(cq_h);
       const bool cv_h = (std::fabs(cq_h) < cfg.n_si() * scfg_.cut_h_si_);
 
+      // vertical cut
       const double cq_v = trDw.y() + scfg_.cut_v_a_ * trUp.y() + scfg_.cut_v_c_;
       h_q_cut_v_bef->Fill(cq_v);
       const bool cv_v = (std::fabs(cq_v) < cfg.n_si() * scfg_.cut_v_si_);
@@ -340,18 +348,12 @@ unsigned int PPSAlignmentWorker::SectorData::process(const CTPPSLocalTrackLiteCo
 
         int idx = (trUp.x() - scfg_.rp_N_.x_slice_min_) / scfg_.rp_N_.x_slice_w_;
         if (idx >= 0 && idx < scfg_.rp_N_.x_slice_n_) {
-          x_slice_plots_N[idx].h_y->Fill(trUp.y());
-          x_slice_plots_N[idx].p_y_diffFN_vs_y->Fill(trUp.y(), trDw.y() - trUp.y());
-          if (debug)
-            x_slice_plots_N[idx].h2_y_diffFN_vs_y->Fill(trUp.y(), trDw.y() - trUp.y());
+          x_slice_plots_N[idx].fill(trUp.y(), trDw.y() - trUp.y(), debug);
         }
 
         idx = (trDw.x() - scfg_.rp_F_.x_slice_min_) / scfg_.rp_F_.x_slice_w_;
         if (idx >= 0 && idx < scfg_.rp_F_.x_slice_n_) {
-          x_slice_plots_F[idx].h_y->Fill(trDw.y());
-          x_slice_plots_F[idx].p_y_diffFN_vs_y->Fill(trDw.y(), trDw.y() - trUp.y());
-          if (debug)
-            x_slice_plots_F[idx].h2_y_diffFN_vs_y->Fill(trDw.y(), trDw.y() - trUp.y());
+          x_slice_plots_F[idx].fill(trDw.y(), trDw.y() - trUp.y(), debug);
         }
       }
     }

CalibPPS/AlignmentGlobal/src/utils.cc

@@ -0,0 +1,45 @@
+/****************************************************************************
+* Authors: 
+*  Jan Kašpar ([email protected]) 
+*  Mateusz Kocot ([email protected])
+****************************************************************************/
+
+#include "CalibPPS/AlignmentGlobal/interface/utils.h"
+
+#include "TF1.h"
+
+#include <memory>
+
+// Fits a linear function to a TProfile.
+int alig_utils::fitProfile(TProfile* p,
+                           double x_mean,
+                           double x_rms,
+                           unsigned int minBinEntries,
+                           unsigned int minNBinsReasonable,
+                           double& sl,
+                           double& sl_unc) {
+  unsigned int n_reasonable = 0;
+  for (int bi = 1; bi <= p->GetNbinsX(); bi++) {
+    if (p->GetBinEntries(bi) < minBinEntries) {
+      p->SetBinContent(bi, 0.);
+      p->SetBinError(bi, 0.);
+    } else {
+      n_reasonable++;
+    }
+  }
+
+  if (n_reasonable < minNBinsReasonable)
+    return 1;
+
+  double x_min = x_mean - x_rms, x_max = x_mean + x_rms;
+
+  auto ff_pol1 = std::make_unique<TF1>("ff_pol1", "[0] + [1]*x");
+
+  ff_pol1->SetParameter(0., 0.);
+  p->Fit(ff_pol1.get(), "Q", "", x_min, x_max);
+
+  sl = ff_pol1->GetParameter(1);
+  sl_unc = ff_pol1->GetParError(1);
+
+  return 0;
+}

CalibPPS/AlignmentGlobal/test/README.md

@@ -4,12 +4,12 @@
 ### Reference dataset
  - `config_reference_cff.py` - configuration (Event Setup) for the reference dataset
  - `input_files_reference_cff.py` - file with vstring of ROOT input files for the reference dataset 
- - `run_distributions_reference_cfg.py` - process configuration for PPSAlignmentWorker. Produces standard ROOT file with reference histograms for x alignment of the test dataset.
+ - `run_distributions_reference_cfg.py` - process configuration for `PPSAlignmentWorker`. Produces a standard ROOT file with reference histograms for the horizontal alignment.
 ### Test dataset
  - `config_cff.py` - configuration (Event Setup) for the test dataset
  - `input_files_cff.py` - file with vstring of ROOT input files for the test dataset
- - `run_distributions_cfg.py` - process configuration for PPSAlignmentWorker. Produces DQMIO ROOT file with histograms for the harvester.
- - `run_analysis_manual_cfg.py` - process configuration for PPSAlignmentHarvester. Produces alignment results.
+ - `run_distributions_cfg.py` - process configuration for `PPSAlignmentWorker`. Produces a DQMIO ROOT file with histograms for the harvester.
+ - `run_analysis_cfg.py` - process configuration for `PPSAlignmentHarvester`. Runs the horizontal, the horizontal relative and the vertical alignment. Produces a text file with results and two ROOT files - one with DQM plots, and the other one with debug plots. Can be easily configured to produce an SQLite file with the results as well - check out the variables at the top of the file.
 
 ## Running instructions
 ```

CalibPPS/AlignmentGlobal/test/run_analysis_cfg.py

@@ -19,7 +19,7 @@
 output_conditions = 'sqlite_file:alignment_results.db'
 
 # Database tag. Used only if write_sqlite_results is set to True.
-output_db_tag = 'CTPPSRPAlignmentCorrectionsData_test'
+output_db_tag = 'CTPPSRPAlignment_real_pcl'
 ###################################
 
 import sys 

CalibPPS/ESProducers/plugins/BuildFile.xml

@@ -2,6 +2,7 @@
   <flags EDM_PLUGIN="1"/>
 
   <use name="CalibPPS/ESProducers"/>
+  <use name="CalibPPS/AlignmentGlobal"/>
 
   <use name="FWCore/Framework"/>
   <use name="FWCore/Utilities"/>