start to work on a more efficient version of MM PDepthProfile user function.

src/external/DepthProfile/inc/PDepthProfile.h

@@ -53,6 +53,7 @@ class PDepthProfileGlobal
     mutable std::vector<Double_t> fPreviousParam;
 
     PRgeHandler *fRgeHandler{nullptr};
+    PRgeDataList fCfd;
 
   ClassDef(PDepthProfileGlobal, 1)
 };

src/external/DepthProfile/src/PDepthProfile.cpp

@@ -2,8 +2,8 @@
 
   PDepthProfile.cpp
 
-  Author: Andreas Suter
+  Authors: Maria Martins, Andreas Suter
-  e-mail: andreas.suter@psi.ch
+  e-mail: maria.martins@psi.ch, andreas.suter@psi.ch
 
 ***************************************************************************/
 
@@ -54,6 +54,20 @@ PDepthProfileGlobal::PDepthProfileGlobal() {
     std::cout << std::endl;
     fValid = false;
   }
+
+  // calculate cumulative frequency distribution of all the rge-files
+  PRgeDataList rgeData = fRgeHandler->GetRgeData();
+  fCfd.resize(fRgeHandler->GetNoOfRgeDataSets());
+  for (unsigned int i=0; i<fCfd.size(); i++) {
+    fCfd[i].energy = rgeData[i].energy;
+    fCfd[i].depth = rgeData[i].depth;
+    fCfd[i].nn.resize(rgeData[i].nn.size());
+    double dval=0.0;
+    for (unsigned int j=0; j<fCfd[i].nn.size(); j++) {
+      dval += rgeData[i].nn[j];
+      fCfd[i].nn[j] = dval;
+    }
+  }
 }
 
 //--------------------------------------------------------------------------
@@ -80,7 +94,7 @@ ClassImp(PDepthProfile)
  * <p>Clean up the global part.
  */
 PDepthProfile::~PDepthProfile() {
-    if ((fDepthProfileGlobal != 0) && fInvokedGlobal) {
+  if ((fDepthProfileGlobal != nullptr) && fInvokedGlobal) {
     delete fDepthProfileGlobal;
     fDepthProfileGlobal = nullptr;
   }
@@ -140,7 +154,13 @@ Bool_t PDepthProfile::GlobalPartIsValid() const {
 //--------------------------------------------------------------------------
 // GetStoppingProbability()
 //--------------------------------------------------------------------------
-
+/**
+ * @brief PDepthProfileGlobal::GetStoppingProbability
+ * @param a
+ * @param b
+ * @param energy
+ * @return
+ */
 double PDepthProfileGlobal::GetStoppingProbability(double a, double b, Double_t energy) const {
 
     // calculation of stopping probability for a given  z interval and experimental energy
@@ -175,12 +195,18 @@ double PDepthProfileGlobal::GetStoppingProbability(double a, double b, Double_t
 //--------------------------------------------------------------------------
 // operator()
 //--------------------------------------------------------------------------
-
+/**
+ * @brief PDepthProfile::operator ()
+ * @param t
+ * @param param
+ * @return
+ */
 Double_t PDepthProfile::operator()(Double_t t, const std::vector <Double_t> &param) const {
-    //verify number of parameters: 2n+1
+  // verify number of parameters: 2n+1, i.e. it has to be an odd number of parameters
-    // parameters: {E,f1, f2, ..., f_n, x1, ..., x_(n-1)}
+  // parameters: {f1, f2, ..., f_n, x1, ..., x_(n-1)}
   assert(param.size() > 2);
   assert(((param.size() - 1) % 2) == 0);
+
   //number of steps: n+1
   int n = (param.size() - 1) / 2;
   std::vector<double> parameters;
@@ -190,7 +216,7 @@ Double_t PDepthProfile::operator()(Double_t t, const std::vector <Double_t> &par
       fPreviousParam.push_back(param[i]);
       parameters.push_back(param[i]);
     }
-    }else{
+  } else {
     for (int i = 0;  i <(n+1); i++) {
       if (param[i]>=0 & param[i] <=1) {
         parameters.push_back(param[i]);
@@ -198,24 +224,24 @@ Double_t PDepthProfile::operator()(Double_t t, const std::vector <Double_t> &par
         parameters.push_back(fPreviousParam[i]);
       }
     }
-        if (n==1){
+    if (n==1) {
-            if(param[n + 1]>0) {
+      if (param[n + 1]>0) {
         parameters.push_back(param[n+1]);
-            } else{
+      } else {
         parameters.push_back(fPreviousParam[n+1]);
       }
     } else {
-            for (int i = n + 1; i < param.size() ; i++) {
+      for (int i=n+1; i<param.size(); i++) {
         //std::cout << "param[i]: " << param[i] << "param[i+1]: " << param[i + 1] << std::endl;
-                if (i!=(param.size() -1)){
+        if (i != (param.size() -1)){
-                    if (param[i] > param[i + 1]) {
+          if (param[i] > param[i+1]) {
             parameters.push_back(fPreviousParam[i+1]);
           } else {
             //std::cout<<"bem" << std::endl;
             parameters.push_back(param[i]);
           }
-                } else{
+        } else {
-                    if (param[i] < param[i - 1])  {
+          if (param[i] < param[i-1]) {
             parameters.push_back(fPreviousParam[i-1]);
           } else {
             //std::cout<<"bem" << std::endl;
@@ -235,7 +261,7 @@ Double_t PDepthProfile::operator()(Double_t t, const std::vector <Double_t> &par
 
   // get stopping probability
   //int l=0;
-    for (int i = 0; i < n + 1; i++) {
+  for (int i=0; i<n+1; i++) {
     double probE;
     double a;
     double b;
@@ -245,9 +271,9 @@ Double_t PDepthProfile::operator()(Double_t t, const std::vector <Double_t> &par
     } else if (i == n) { //prob between x_(n-1) and inf
       a = parameters[2 * n];
       b = 179;
-        } else {//prob between x_1-x_2, ..., x_n-x_(n-1)
+    } else { //prob between x_1-x_2, ..., x_n-x_(n-1)
-            a = parameters[n + i];
+      a = parameters[n+i];
-            b = parameters[n + 1 + i];
+      b = parameters[n+1+i];
     }
     probE = fDepthProfileGlobal->GetStoppingProbability(a, b, energy);
     probability.push_back(probE);
@@ -255,14 +281,11 @@ Double_t PDepthProfile::operator()(Double_t t, const std::vector <Double_t> &par
 
   double fit=0;
 
-    for (int j = 0; j < n + 1; j++) {
+  for (int j=0; j<n+1; j++) {
     fit = fit + parameters[j] * probability[j];
   }
 
-   /* std::cout << "Energy  " << energy << std::endl;
-    std::cout << "FRACTION  " << fit << std::endl;*/
   return fit;
-
 }