Add together an arbitrary number of signals plus background

Newbie
1

Hello,

The title explains what I want to do. If it matters, I am fitting Brillouin spectra and don’t always know how many peaks there will be or which I want to fit.

This what I have so far:

void region_1() {

  TH1F* hist = ReadBS("Region_1.csv");
  if ( hist->GetEntries() == 0 ) // Empty histo.
    return;

  // hist->Draw("AP");
  
  // Import into RooDataHist
  double xlow, xhigh;
  xlow = hist->GetXaxis()->GetXmin();
  xhigh = hist->GetXaxis()->GetXmax();
  
  RooRealVar x("x","x",xlow,xhigh);
  RooDataHist data("data","dataset with x",x,hist);
  RooPlot* frame = x.frame() ;
  data.plotOn(frame) ;
  frame->Draw();

  std::vector<lorentzian> lorentzians;
  
  // Lorentzian
  // A/pi ((G/2)/((x-P)^2 + (G/2)^2))
  // RooRealVar A1("A1","Amplitude of Lorentzian", 0.5, 0,100);
  // RooRealVar P1("P1", "Peak of Lorentzian", -74,-77, -72);
  // RooRealVar G1("G1", "Width of lorentzian",0.1,0,2);
  // RooGenericPdf L1{"L1", "A1 * ((G1/2)/((x-P1)^2 + (G1/2)^2))", RooArgList(A1,P1,G1,x)};

  // The info used here should come from TSpectrum
  for (int i = 0; i < 2; ++i) {
    lorentzian l;
    if ( i == 0 ) {
    l.A = 0.5; l.Alow = 0.0; l.Ahigh = 100.0;
    l.P = -74; l.Plow = -77; l.Phigh = -72;
    l.G = 0.1; l.Glow = 0.0; l.Ghigh = 2.0;
    } else {// -77,-80, -75
    l.A = 0.5; l.Alow = 0.0; l.Ahigh = 100.0;
    l.P = -77; l.Plow = -80; l.Phigh = -75;
    l.G = 0.1; l.Glow = 0.0; l.Ghigh = 2.0;
    }
    lorentzians.push_back(l);
  }
  std::cout << "lorentsians: " << lorentzians.size() << "\n";

  // DEBUG
  // for (int i = 0; i < 10; ++i)
  //   std::cout << "counter: " << UNIQUE_ID(poo) << "\n";
  
  // Now construct a generic PDF for our lorentzian
  std::vector<RooGenericPdf> gpdfs;
  std::vector<RooArgList> args;
  for (int i = 0; i < lorentzians.size(); ++i) {
    RooRealVar A("A", "Amplitude of Lorentzian", lorentzians.at(i).A,lorentzians.at(i).Alow,lorentzians.at(i).Ahigh);
    RooRealVar P("P", "Peak of Lorentzian", lorentzians.at(i).P,lorentzians.at(i).Plow,lorentzians.at(i).Phigh);
    RooRealVar G("G", "FWHM of Lorentzian", lorentzians.at(i).G,lorentzians.at(i).Glow,lorentzians.at(i).Ghigh);
    RooArgList arg;
    arg.add(A);
    arg.add(P);
    arg.add(G);
    arg.add(x);
    RooGenericPdf L{"L", "A * ((G/2)/((x-P)^2 + (G/2)^2))", arg};
    // RooGenericPdf LL(L);
    gpdfs.push_back(L);
  }
  std::cout << "gpdfs: " << gpdfs.size() << "\n";
  std::cout << "args: " << args.size() << "\n";
  
  // Need to sum L1 and L2 and fit that.
  RooRealVar f("f","f",0,1);
  RooAddPdf S1("S1", "Sum of Lorentians", args,f);

  S1.fitTo(data);
  // S1.plotOn(frame);
  // frame->Draw();

  // How would one assemble an arbitrary number of peaks plus a background?
  delete hist;
}

It won’t compile because RooAddPdf does not take a std::vector of RooArgList. I have looked through some of the examples in the tutorials but I can’t find one that clearly does what I want.

I appreciate any help, and by all means point me to the tutorial that I may have missed.

sprock

Please read tips for efficient and successful posting and posting code

ROOT Version: Not Provided
Platform: Not Provided
Compiler: Not Provided

2

(post deleted by author)