Merge branch 'EK' into 'master'

Modify BorderFinder in 1D mode, add gitignore

See merge request pwg-c2f/centrality!2

.gitignore

+# Created by .ignore support plugin (hsz.mobi)
+### JetBrains template
+# Covers JetBrains IDEs: IntelliJ, RubyMine, PhpStorm, AppCode, PyCharm, CLion, Android Studio and WebStorm
+# Reference: https://intellij-support.jetbrains.com/hc/en-us/articles/206544839
+
+# User-specific stuff
+.idea/**/workspace.xml
+.idea/**/tasks.xml
+.idea/**/usage.statistics.xml
+.idea/**/dictionaries
+.idea/**/shelf
+
+# Sensitive or high-churn files
+.idea/**/dataSources/
+.idea/**/dataSources.ids
+.idea/**/dataSources.local.xml
+.idea/**/sqlDataSources.xml
+.idea/**/dynamic.xml
+.idea/**/uiDesigner.xml
+.idea/**/dbnavigator.xml
+
+# Gradle
+.idea/**/gradle.xml
+.idea/**/libraries
+
+# Gradle and Maven with auto-import
+# When using Gradle or Maven with auto-import, you should exclude module files,
+# since they will be recreated, and may cause churn.  Uncomment if using
+# auto-import.
+# .idea/modules.xml
+# .idea/*.iml
+# .idea/modules
+
+# CMake
+cmake-build-*/
+
+# Mongo Explorer plugin
+.idea/**/mongoSettings.xml
+
+# File-based project format
+*.iws
+
+# IntelliJ
+out/
+
+# mpeltonen/sbt-idea plugin
+.idea_modules/
+
+# JIRA plugin
+atlassian-ide-plugin.xml
+
+# Cursive Clojure plugin
+.idea/replstate.xml
+
+# Crashlytics plugin (for Android Studio and IntelliJ)
+com_crashlytics_export_strings.xml
+crashlytics.properties
+crashlytics-build.properties
+fabric.properties
+
+# Editor-based Rest Client
+.idea/httpRequests
+### C++ template
+# Prerequisites
+*.d
+
+# Compiled Object files
+*.slo
+*.lo
+*.o
+*.obj
+
+# Precompiled Headers
+*.gch
+*.pch
+
+# Compiled Dynamic libraries
+*.so
+*.dylib
+*.dll
+
+# Fortran module files
+*.mod
+*.smod
+
+# Compiled Static libraries
+*.lai
+*.la
+*.a
+*.lib
+
+# Executables
+*.exe
+*.out
+*.app
+### CMake template
+CMakeCache.txt
+CMakeFiles
+CMakeScripts
+Testing
+Makefile
+cmake_install.cmake
+install_manifest.txt
+compile_commands.json
+CTestTestfile.cmake
+### C template
+# Prerequisites
+
+# Object files
+*.ko
+*.elf
+
+# Linker output
+*.ilk
+*.map
+*.exp
+
+# Precompiled Headers
+
+# Libraries
+
+# Shared objects (inc. Windows DLLs)
+*.so.*
+
+# Executables
+*.i*86
+*.x86_64
+*.hex
+
+# Debug files
+*.dSYM/
+*.su
+*.idb
+*.pdb
+
+# Kernel Module Compile Results
+*.mod*
+*.cmd
+.tmp_versions/
+modules.order
+Module.symvers
+Mkfile.old
+dkms.conf
+### Autotools template
+# http://www.gnu.org/software/automake
+
+Makefile.in
+/ar-lib
+/mdate-sh
+/py-compile
+/test-driver
+/ylwrap
+
+# http://www.gnu.org/software/autoconf
+
+autom4te.cache
+/autoscan.log
+/autoscan-*.log
+/aclocal.m4
+/compile
+/config.guess
+/config.h.in
+/config.log
+/config.status
+/config.sub
+/configure
+/configure.scan
+/depcomp
+/install-sh
+/missing
+/stamp-h1
+
+# https://www.gnu.org/software/libtool/
+
+/ltmain.sh
+
+# http://www.gnu.org/software/texinfo
+
+/texinfo.tex
+
+# http://www.gnu.org/software/m4/
+
+m4/libtool.m4
+m4/ltoptions.m4
+m4/ltsugar.m4
+m4/ltversion.m4
+m4/lt~obsolete.m4

main.cpp

@@ -16,50 +16,47 @@
 #include "TH2.h"
 
 int main(int argc, char **argv) {
-    auto start = std::chrono::system_clock::now();
-    ROOT::EnableImplicitMT(2);
-    
-    if (argc < 3)
-    {
-        std::cout << "Not enough arguments! Please use:" << std::endl;
-        std::cout << "   ./main filename histoname is2d" << std::endl;
-        return -1;
-    }
-    
-    std::unique_ptr <TFile> fIn {TFile::Open(argv[1], "read")};
-    std::string outfilename = "test.root";
-    if ( strcmp( argv[3], "false") == 0 )
-    {
-        std::unique_ptr <TH1F> histo {(TH1F*) (fIn->Get(argv[2]))};
-        Centrality::BordersFinder bf;
-        bf.SetHisto(*histo);
-        bf.SetRanges( 10,0,100 );   // number of bins, min, max value
+  auto start = std::chrono::system_clock::now();
+  ROOT::EnableImplicitMT(2);
+
+  if (argc < 3) {
+    std::cout << "Not enough arguments! Please use:" << std::endl;
+    std::cout << "   ./main filename histoname is2d" << std::endl;
+    return -1;
+  }
+
+  std::unique_ptr<TFile> fIn{TFile::Open(argv[1], "read")};
+  std::string outfilename = "test.root";
+  if (strcmp(argv[3], "false") == 0) {
+    std::unique_ptr<TH1F> histo{(TH1F *) (fIn->Get(argv[2]))};
+    Centrality::BordersFinder bf;
+    bf.SetHisto(*histo);
+//        bf.SetRanges( 10,0,100 );   // number of bins, min, max value
+    bf.SetRanges({0, 5, 10, 15, 20, 25, 30, 35, 40, 45, 50, 55, 60, 70, 80, 90, 100});  // centrality bins borders with array
+    bf.SetLimits(600, 8000);
+    bf.IsSpectator(true);  // true if impact parameter b correlated with estimator (spectators eneggy),
+    // false - anticorrelated (multiplicity of produced particles)
+
+    bf.FindBorders();
+    bf.SaveBorders(outfilename, "centr_getter_1d");
+  } else if (strcmp(argv[3], "true") == 0) {
+    std::unique_ptr<TH2F> histo2d{(TH2F *) (fIn->Get(argv[2]))};
+    Centrality::BordersFinder2D bf;
+    bf.SetHisto2D(std::move(*histo2d));
+    std::unique_ptr<TH1F> h1d = bf.Convert();
+
+    bf.SetHisto(*h1d);
+    bf.SetRanges(10, 0, 100);   // number of bins, min, max value
     //   bf.SetRanges( {0,10,30,60,100} );  // centrality bins borders with array
-        bf.IsSpectator(false);  // true if impact parameter b correlated with estimator (spectators eneggy), 
-                                    // false - anticorrelated (multiplicity of produced particles) 
+    bf.IsSpectator(false);  // true if impact parameter b correlated with estimator (spectators eneggy),
+    // false - anticorrelated (multiplicity of produced particles)
+
+    bf.FindBorders();
+    bf.SaveBorders2D(outfilename, "centr_getter_2d");
+  }
 
-        bf.FindBorders();
-        bf.SaveBorders(outfilename, "centr_getter_1d");
-    }
-    else if ( strcmp( argv[3], "true") == 0 )
-    {
-        std::unique_ptr <TH2F> histo2d {(TH2F*) (fIn->Get(argv[2]))};
-        Centrality::BordersFinder2D bf;
-        bf.SetHisto2D( std::move(*histo2d) );
-        std::unique_ptr<TH1F> h1d = bf.Convert();
-        
-        bf.SetHisto(*h1d);
-        bf.SetRanges( 10,0,100 );   // number of bins, min, max value
-        //   bf.SetRanges( {0,10,30,60,100} );  // centrality bins borders with array
-        bf.IsSpectator(false);  // true if impact parameter b correlated with estimator (spectators eneggy), 
-                                    // false - anticorrelated (multiplicity of produced particles) 
-        
-        bf.FindBorders();
-        bf.SaveBorders2D(outfilename, "centr_getter_2d");
-    }
-        
-    auto end = std::chrono::system_clock::now();
-    std::chrono::duration<double> elapsed_seconds = end - start;
-    std::cout << "elapsed time: " << elapsed_seconds.count() << " s\n";
-    return 0;
+  auto end = std::chrono::system_clock::now();
+  std::chrono::duration<double> elapsed_seconds = end - start;
+  std::cout << "elapsed time: " << elapsed_seconds.count() << " s\n";
+  return 0;
 }

src/BordersFinder.cpp

@@ -3,6 +3,7 @@
 #include "BordersFinderHelper.h"
 
 #include <iostream>
+#include <TGraph.h>
 
 #include "TFile.h"
 
@@ -10,16 +11,52 @@ namespace Centrality {
 
 void BordersFinder::FindBorders()
 {
+    using namespace std;
+
     if (ranges_.size() < 2) return;
     if (norm_ == -1) norm_ = histo_.Integral(0, histo_.GetNbinsX());
-    if (!isspectator_) std::reverse(std::begin(ranges_), std::end(ranges_));
-    
+    if (!isSpectator_) std::reverse(std::begin(ranges_), std::end(ranges_));
+
+    auto axis = histo_.GetXaxis();
+
+    double xLo = (applyLimits_)? xLo_ : axis->GetXmin();
+    double xHi = (applyLimits_)? xHi_ : axis->GetXmax();
+
+
+    int n = axis->GetNbins();
+
+    double *histIntegral = histo_.GetIntegral();
+    double x[n];
+    for (int i = 0; i < n; ++i) {
+        x[i] = axis->GetBinCenter(i+1);
+    }
+
+    TGraph intVsXGraph(n, x, histIntegral);
+    intVsXGraph.SetBit(TGraph::kIsSortedX);
+
+    TGraph xVsIntGraph(n, histIntegral, x);
+    xVsIntGraph.SetBit(TGraph::kIsSortedX);
+
+    double intLo = intVsXGraph.Eval(xLo);
+    double intHi = intVsXGraph.Eval(xHi);
+    double norm = intHi - intLo;
+
+    auto cX = [=] (double x) { return 100./norm*(intVsXGraph.Eval(x) - intVsXGraph.Eval(xLo)); };
+    auto xC = [=] (double c) { return xVsIntGraph.Eval((c/100.)*norm + intLo); };
+
+    for (auto cc : ranges_) {
+        double xx = isSpectator_? xC(cc) : xC(100 - cc);
+        cout << cc << "%" << ", border:" << xx << endl;
+        borders_.push_back(xx);
+    }
+
+/*
     uint iSlice{0};
     long int entriesCurrent{0};
     
     for (int iBin=1; iBin<=histo_.GetNbinsX() && iSlice<ranges_.size() ; ++iBin)
     {
-        const float step = isspectator_ ? ranges_.at(iSlice) : 100. - ranges_.at(iSlice);
+        const float step = isSpectator_ ? ranges_.at(iSlice) : 100. - ranges_.at(iSlice);
         const long int entriesNeeeded = step/100. * norm_;
         entriesCurrent += histo_.GetBinContent(iBin);
         
@@ -34,6 +71,7 @@ void BordersFinder::FindBorders()
             iSlice++;
         }
     }
+*/
 }
 
 void BordersFinder::SaveBorders(const std::string &filename, const std::string &getter_name)

src/BordersFinder.h

@@ -11,50 +11,56 @@
 #include "TH1.h"
 
 namespace Centrality {
-    
+
 class BordersFinder {
 
-public:
-    
-    BordersFinder(){}
-    
-    void FindBorders();
-    void SaveBorders(const std::string &filename, const std::string &getter_name);
-
-
-    void SetHisto(const TH1F &h)  { histo_ = h; }
-    TH1F& GetHisto() { return histo_; }   // not const to use Draw etc
-    
-    void SetNormalization (long int norm) { norm_ = norm; }
-    long int GetNormalization () const { return norm_; }
-    
-    void SetRanges(const std::vector<float> &ranges) { ranges_ = ranges; }
-    void SetRanges(int n, float min, float max) 
-    {
-        ranges_.clear();
+ public:
+
+  BordersFinder() {}
+
+  void FindBorders();
+  void SaveBorders(const std::string &filename, const std::string &getter_name);
+
+  void SetHisto(const TH1F &h) { histo_ = h; }
+  TH1F &GetHisto() { return histo_; }   // not const to use Draw etc
+
+  void SetNormalization(long int norm) { norm_ = norm; }
+  long int GetNormalization() const { return norm_; }
+
+  void SetRanges(const std::vector<float> &ranges) { ranges_ = ranges; }
+  void SetRanges(int n, float min, float max) {
+    ranges_.clear();
 //         ranges_.reserve(n+1);
-        
-        for (int i=0; i<=n; ++i)
-            ranges_.push_back( min + i*(max-min)/n );
-    }
-    
-    void IsSpectator(bool is=true) { isspectator_ = is; }
-    
-    const std::vector<float>& GetRanges() const { return ranges_; }
-    const std::vector<double>& GetBorders() const { return borders_; }
-    bool GetIsSpectator() const { return isspectator_; }
-    
-    
-private:
-    
-    TH1F histo_;
-    long int norm_{-1};
-        
-    std::vector<float> ranges_{};
-    std::vector<double> borders_{};
-    
-    bool isspectator_{false};
-    
+
+    for (int i = 0; i <= n; ++i)
+      ranges_.push_back(min + i * (max - min) / n);
+  }
+
+  void SetLimits(double xLo, double xHi) {
+    xLo_ = xLo;
+    xHi_ = xHi;
+    applyLimits_ = true;
+  }
+
+  void IsSpectator(bool is = true) { isSpectator_ = is; }
+
+  const std::vector<float> &GetRanges() const { return ranges_; }
+  const std::vector<double> &GetBorders() const { return borders_; }
+  bool GetIsSpectator() const { return isSpectator_; }
+
+ private:
+
+  TH1F histo_;
+  Double_t norm_{-1};
+
+  std::vector<float> ranges_{};
+  std::vector<double> borders_{};
+
+  bool isSpectator_{false};
+
+  bool applyLimits_{false};
+  double xLo_{-1};
+  double xHi_{-1};
 };
 }