Commit 6182c58f authored by Jérémie Dudouet's avatar Jérémie Dudouet
Browse files

Merge branch 'Cubix' into 'preprod'

Cubix

See merge request !70
parents 6aa4388f 739e0701
......@@ -264,16 +264,14 @@ void CXBgdUtility::DoSubtract()
TH1 *hist = fMainWindow->GetHisto();
if(hist == nullptr || hist->InheritsFrom("TH2"))
{
if(hist == nullptr || hist->InheritsFrom("TH2")) {
cout<<"No 1D istogram found, ignored"<<endl;
return;
}
TH1 *bkgrd = dynamic_cast<TH1*>(gPad->GetListOfPrimitives()->FindObject(Form("%s_BG",hist->GetName())));
if(bkgrd == nullptr)
{
if(bkgrd == nullptr) {
cout<<"No background histogram built in the pad for hist: "<<hist->GetName()<<endl;
return;
}
......@@ -298,20 +296,16 @@ void CXBgdUtility::ToggleBckSupp()
if(gPad == nullptr)
return;
if(fActivate->GetState() == kButtonDown)
{
if(fActivate->GetState() == kButtonDown) {
SetButtonsStatus(true);
GetParams();
}
else
{
else {
SetButtonsStatus(false);
for(int i=gPad->GetListOfPrimitives()->GetEntries()-1 ; i>=0 ; i--)
{
for(int i=gPad->GetListOfPrimitives()->GetEntries()-1 ; i>=0 ; i--) {
TObject *o = gPad->GetListOfPrimitives()->At(i);
if(o->InheritsFrom("TH1") && ((TString)o->GetName()).EndsWith("_BG"))
{
if(o->InheritsFrom("TH1") && ((TString)o->GetName()).EndsWith("_BG")) {
delete o;
break;
}
......@@ -326,8 +320,7 @@ void CXBgdUtility::SetButtonsStatus(bool on)
{
fNItersButton1D->SetState(on);
for(int i=1 ; i<fListOfButtons1D->GetEntries() ; i++)
{
for(int i=1 ; i<fListOfButtons1D->GetEntries() ; i++) {
auto *btn = dynamic_cast<TGButton*>(fListOfButtons1D->At(i));
if(!fSmoothing && i>=8 && i<=14)
......@@ -359,8 +352,7 @@ void CXBgdUtility::GetParams()
else if(fFilterOrder1D == TSpectrum::kBackOrder8) Option += "BackOrder8";
if(!fSmoothing) Option += "nosmoothing";
else
{
else {
if(fSmoothingWindow == TSpectrum::kBackSmoothing3) Option += "BackSmoothing3";
else if(fSmoothingWindow == TSpectrum::kBackSmoothing5) Option += "BackSmoothing5";
else if(fSmoothingWindow == TSpectrum::kBackSmoothing7) Option += "BackSmoothing7";
......@@ -372,10 +364,9 @@ void CXBgdUtility::GetParams()
}
if(fCompton) Option += "Compton";
TH1 *hist = dynamic_cast<TH1F*>(fMainWindow->GetHisto());
TH1 *hist = dynamic_cast<TH1*>(fMainWindow->GetHisto());
if(hist == nullptr || hist->InheritsFrom("TH2"))
{
if(hist == nullptr || hist->InheritsFrom("TH2")) {
cout<<"No 1D istogram found, ignored"<<endl;
return;
}
......@@ -409,24 +400,19 @@ void CXBgdUtility::HandleWindowButtons(Int_t id)
if (id == -1)
id = btn->WidgetId();
if(id == 0)
{
if(id == 0) {
fDirectionButton1D[1]->SetState(kButtonUp);
fDirection1D = TSpectrum::kBackIncreasingWindow;
}
if(id == 1)
{
if(id == 1) {
fDirectionButton1D[0]->SetState(kButtonUp);
fDirection1D = TSpectrum::kBackDecreasingWindow;
}
if(id == 10)
{
if(id == 10) {
fDirectionButton2D[1]->SetState(kButtonUp);
fDirection2D = TSpectrum::kBackIncreasingWindow;
}
if(id == 11)
{
if(id == 11) {
fDirectionButton2D[0]->SetState(kButtonUp);
fDirection2D = TSpectrum::kBackDecreasingWindow;
}
......@@ -465,15 +451,11 @@ void CXBgdUtility::HandleSmoothButtons(Int_t id)
if (id == -1)
id = btn->WidgetId();
if(id == -2)
{
if(btn->GetState() == kButtonDown)
{
if(id == -2) {
if(btn->GetState() == kButtonDown) {
fSmoothing = true;
for(int i=0 ; i<7 ; i++)
{
if(fSmoothingButton1D.at(i)->IsDown())
{
for(int i=0 ; i<7 ; i++) {
if(fSmoothingButton1D.at(i)->IsDown()) {
fSmoothingButton1D.at(i)->SetState(kButtonDown);
fSmoothingWindow = 3+2*i;
}
......@@ -481,8 +463,7 @@ void CXBgdUtility::HandleSmoothButtons(Int_t id)
fSmoothingButton1D.at(i)->SetState(kButtonUp);
}
}
else
{
else {
fSmoothing = false;
for(int i=0 ; i<7 ; i++)
fSmoothingButton1D.at(i)->SetState(kButtonDisabled);
......@@ -490,8 +471,7 @@ void CXBgdUtility::HandleSmoothButtons(Int_t id)
return;
}
for(int i=0 ; i<7 ; i++)
{
for(int i=0 ; i<7 ; i++) {
if(id != i)
fSmoothingButton1D.at(i)->SetState(kButtonUp);
else
......@@ -510,8 +490,7 @@ void CXBgdUtility::Do2DEvaluation(TH1 *hist_in)
if(hist && hist->InheritsFrom("TH2"))
fCurrent2DHist = dynamic_cast<TH2*>(hist);
else
{
else {
cout<<"No 2D histogram found in the current pad, ignored"<<endl;
fCurrent2DHist = nullptr;
return;
......
......@@ -180,7 +180,7 @@ void CXFit::Fit()
ROOT::Math::MinimizerOptions::SetDefaultTolerance(fPlayer->GetTolerance());
ROOT::Math::MinimizerOptions::SetDefaultPrintLevel(fPlayer->GetPrintLevel());
// fPlayer->GetFitResultsBox()->RemoveAll();
// fPlayer->GetFitResultsBox()->RemoveAll();
Double_t DefFWHM = fPlayer->fNE_FWHM[0]->GetNumber();
Double_t DefFWHM_min = fPlayer->fNE_FWHM[1]->GetNumber();
......@@ -234,7 +234,7 @@ void CXFit::Fit()
fFitFunction->SetParameter(2, 0);
fFitFunction->SetParLimits(2, -50., 0.);
fFitFunction->SetParameter(3, 0.);
// fFitFunction->SetParLimits(3, -50., 0.);
// fFitFunction->SetParLimits(3, -50., 0.);
if(!fPlayer->fUsePol1)
fFitFunction->FixParameter(2,0);
......@@ -376,9 +376,16 @@ void CXFit::Fit()
Double_t FWHM_L = peak->GetX(Max/2,fBackgd[0],Mean,1e-6);
Double_t FWHM_L_err = peak->GetX((Max-MaxErr)/2,fBackgd[0],Mean,1e-6);
Double_t F01_L = Mean-peak->GetX(Max/10.,fBackgd[0],Mean,1e-6);
Double_t F01_R = peak->GetX(Max/10.,Mean,fBackgd[1],1e-6)-Mean;
Double_t FWHM_R = peak->GetX(Max/2,Mean,fBackgd[1],1e-6);
Double_t FWHM_R_err = peak->GetX((Max-MaxErr)/2,Mean, fBackgd[1],1e-6);
Double_t LeftTailParam = F01_L/(FWHM*0.5);
Double_t RightTailParam = F01_R/(FWHM*0.5);
cout<<F01_L<<" "<<FWHM<<" "<<endl;
Double_t FWHM_Real = FWHM_R-FWHM_L;
Double_t FWHM_Real_err = (FWHM_R_err-FWHM_L_err)-FWHM_Real;
peak->SetParameter(1,1);//with backgroud
......@@ -403,12 +410,12 @@ void CXFit::Fit()
fPlayer->PrintInListBox(text.str(),kInfo);
text.str("");
text<<left<<setw(11)<<"L Tail"<<": "<<setprecision(7)<<setw(10)<<LeftT<<" ("<<setprecision(7)<<setw(10)<<LeftTErr<<")";
text<<left<<setw(11)<<"L Tail"<<": "<<setprecision(7)<<setw(10)<<LeftTailParam;
cout<<text.str()<<endl;
fPlayer->PrintInListBox(text.str(),kInfo);
text.str("");
text<<left<<setw(11)<<"R Tail"<<": "<<setprecision(7)<<setw(10)<<Right<<" ("<<setprecision(7)<<setw(10)<<RightErr<<")";
text<<left<<setw(11)<<"R Tail"<<": "<<setprecision(7)<<setw(10)<<RightTailParam;
cout<<text.str()<<endl;
fPlayer->PrintInListBox(text.str(),kInfo);
text.str("");
......
......@@ -275,7 +275,7 @@ void CXHist1DCalib::Calibrate()
vector < Fitted > FitResults = fRecalEnergy->GetFitResults();
if(fVerboseLevel->GetNumber()==0) {
if(fVerboseLevel->GetNumber()==0 && fRecalEnergy->fCalibFunction) {
cout<< left << scientific << setprecision(6);
cout<< hist->GetName()<<": ";
cout << setw(14) << fRecalEnergy->fCalibFunction->GetParameter(0);
......@@ -283,14 +283,12 @@ void CXHist1DCalib::Calibrate()
cout<<endl;
}
double xmin=-1;
double xmax=-1;
int NGoodPeak=0;
for(size_t i=0 ; i<FitResults.size() ; i++)
{
for(size_t i=0 ; i<FitResults.size() ; i++) {
Fitted FitRes = FitResults[i];
if(!FitRes.good) continue;
......@@ -345,20 +343,24 @@ void CXHist1DCalib::Calibrate()
i += NSubPeaks-1;
}
hist->GetXaxis()->SetRangeUser(xmin-(xmax-xmin)*0.1,xmax+(xmax-xmin)*0.1);
if(xmin!=-1) hist->GetXaxis()->SetRangeUser(xmin-(xmax-xmin)*0.1,xmax+(xmax-xmin)*0.1);
fMainWindow->RefreshPads();
if(fCalibCanvas != nullptr) delete fCalibCanvas;
fCalibCanvas = new TCanvas;
fCalibCanvas->SetName("CalibrationResults");
fCalibCanvas->SetTitle("Calibration Results");
fCalibCanvas->Divide(1,2,0.0001,0.0001);
fCalibCanvas->cd(1);
fRecalEnergy->fCalibGraph->Draw("ap");
fRecalEnergy->fCalibFunction->Draw("same");
fCalibCanvas->cd(2);
fRecalEnergy->fResidueGraph->Draw("ape");
if(FitResults.size()>1 && fRecalEnergy->fCalibFunction) {
if(fCalibCanvas != nullptr) delete fCalibCanvas;
fCalibCanvas = new TCanvas;
fCalibCanvas->SetName("CalibrationResults");
fCalibCanvas->SetTitle("Calibration Results");
fCalibCanvas->Divide(1,2,0.0001,0.0001);
fCalibCanvas->cd(1);
fRecalEnergy->fCalibGraph->Draw("ap");
fRecalEnergy->fCalibFunction->Draw("same");
fCalibCanvas->cd(2);
fRecalEnergy->fResidueGraph->Draw("ape");
fCalibCanvas->Update();
fCalibCanvas->Modified();
}
fMainWindow->GetCanvas()->cd();
}
......
......@@ -123,11 +123,13 @@ LevelScheme *CXLevelSchemePlayer::DrawArrows(TString ListOfNuclei,TH1 *h, TStrin
fGuiLSPlayer->GetLifeTime(fMinLifeTime,fMaxLifeTime);
fCurrentHist = h;
bool dodraw = (fCurrentHist!=nullptr);
fCurrentHist->GetYaxis()->UnZoom();
gPad->Update();
gPad->Modified();
if(dodraw) {
fCurrentHist->GetYaxis()->UnZoom();
gPad->Update();
gPad->Modified();
}
Int_t IColor = 0;
......@@ -202,7 +204,7 @@ LevelScheme *CXLevelSchemePlayer::DrawArrows(TString ListOfNuclei,TH1 *h, TStrin
return RefLevel;
}
LevelScheme *CXLevelSchemePlayer::DrawArrowsForNuc(TString NucName, TString DataSet)
LevelScheme *CXLevelSchemePlayer::DrawArrowsForNuc(TString NucName, TString DataSet, bool do_draw)
{
LevelScheme *RefLevel = nullptr;
......@@ -223,7 +225,7 @@ LevelScheme *CXLevelSchemePlayer::DrawArrowsForNuc(TString NucName, TString Data
CalcYrast(lev);
for(Int_t i=0 ; i<lev->GetLinks().GetSize() ; i++) {
const GammaLink *gamma=(GammaLink*)lev->GetLinks().At(i);
DrawArrow(fCurrentHist,gamma,nuc,"ENSDF");
if(do_draw) DrawArrow(fCurrentHist,gamma,nuc,"ENSDF");
}
}
......@@ -249,7 +251,7 @@ LevelScheme *CXLevelSchemePlayer::DrawArrowsForNuc(TString NucName, TString Data
for(Int_t i=0 ; i<lev->GetLinks().GetSize() ; i++) {
const GammaLink *gamma=(GammaLink*)lev->GetLinks().At(i);
DrawArrow(fCurrentHist,gamma,nuc,"EXP");
if(do_draw) DrawArrow(fCurrentHist,gamma,nuc,"EXP");
}
}
RefLevel = lev;
......
......@@ -73,7 +73,7 @@ public:
void CleanArrows();
LevelScheme *DrawArrows(TString ListOfNuclei, TH1 *h, TString DataSet);
LevelScheme *DrawArrowsForNuc(TString NucName, TString DataSet);
LevelScheme *DrawArrowsForNuc(TString NucName, TString DataSet, bool do_draw=true);
Bool_t CanReplot(){return fCanReplot;}
......
......@@ -139,6 +139,8 @@ void GwRecalEnergy::StartCalib()
else
np = PeakSearch2(specData, specFromDef, specToDef, specFWHMdef, specAMPLdef, specPeaks);
cout<<"yo :" <<np<<endl;
if(np) {
int fp = FitPeaks( Verbosity );
if(fp) {
......@@ -160,6 +162,7 @@ void GwRecalEnergy::StartCalib()
Peaks.clear();
if(Verbosity > 0)
cout<<Form("\n");
return;
}
if(fCalibOrder) {
......@@ -821,9 +824,9 @@ int GwRecalEnergy::xP_Next1(float *yDat, int yLen)
}
else { // transition to negative
if(widthP >= xP_nMinWidthP && // positive lobe wide enough
//widthP < 4*xP_nMinWidthP && // but not too wide
xP_1 > xP_0 && // and has a maximum
yMax > xP_minAmpl) { // which is large enough
//widthP < 4*xP_nMinWidthP && // but not too wide
xP_1 > xP_0 && // and has a maximum
yMax > xP_minAmpl) { // which is large enough
state = stateN; // initiate a negative lobe sequence
xP_3 = xP_4 = xP_5 = nn;
yMin = yInp;
......@@ -842,8 +845,8 @@ int GwRecalEnergy::xP_Next1(float *yDat, int yLen)
xP_5 = nn;
widthN++; // incr negative width and stay in this state
if(widthN > xP_nMinWidthN && // second positive lobe wide enough for a valid sequence
yInp > yMin && // and signal is beyond minimum
yMin < -xP_minAmpl) { // and minimum large enough
yInp > yMin && // and signal is beyond minimum
yMin < -xP_minAmpl) { // and minimum large enough
float cmp = abs(yMax/yMin);
if(cmp < 0.2f || cmp > 2.f)
state = stateUnknown;
......@@ -905,7 +908,7 @@ int GwRecalEnergy::xP_Next2(float *yDat, int yLen)
}
else { // transition to negative
if(widthP1 >= xP_nMinWidthP1 && // positive lobe wide enough
xP_1 > xP_0) { // and has a maximum
xP_1 > xP_0) { // and has a maximum
state = stateN; // initiate a negative lobe sequence
xP_3 = xP_4 = nn;
yMin = yInp;
......@@ -926,9 +929,9 @@ int GwRecalEnergy::xP_Next2(float *yDat, int yLen)
}
else { // transition to positive
if(widthN >= xP_nMinWidthN && // negative lobe wide enough
yMin < -xP_minAmpl && // its amplitude large enough
yMin < -1.1*yMax1 && // and minimum-to-maximum is acceptable
yMin > -4.5*yMax1) { // (exact value is -0.5*exp(1.5)) = -2.2408)
yMin < -xP_minAmpl && // its amplitude large enough
yMin < -1.1*yMax1 && // and minimum-to-maximum is acceptable
yMin > -4.5*yMax1) { // (exact value is -0.5*exp(1.5)) = -2.2408)
state = stateP2; // initiate the second positive lobe sequence
xP_5 = xP_6 = nn;
yMax2 = yInp;
......@@ -951,7 +954,7 @@ int GwRecalEnergy::xP_Next2(float *yDat, int yLen)
xP_6 = nn;
widthP2++; // incr positive width and stay in this state
if(widthP2 > xP_nMinWidthP2 && // second positive lobe wide enough for a valid sequence
yInp < yMax2) { // and signal is beyond maximum
yInp < yMax2) { // and signal is beyond maximum
return nn; // this is a good trigger
}
}
......@@ -1050,6 +1053,10 @@ Int_t GwRecalEnergy::EROOTCalibration()
{
vector <Fitted>::iterator Iter;
double bestSlope = 1.;
Float_t Xmin = 0;
Float_t Xmax = 0.;
if(Energies.size() == 1) {
// select the peak with the largest area
int pn_max = 0;
......@@ -1058,118 +1065,112 @@ Int_t GwRecalEnergy::EROOTCalibration()
pn_max = nn;
}
}
double slope_ = Energies[0] / Peaks[pn_max].posi;
bestSlope = Energies[0] / Peaks[pn_max].posi;
Peaks[pn_max].good = true;
Peaks[pn_max].erefindex = 0;
return 0;
Xmin = Peaks[pn_max].posi-100;
Xmax = Peaks[pn_max].posi+100;
}
if( Peaks.size() < 2 || Energies.size() < 2) {
else if( Peaks.size() < 2 || Energies.size() < 2){
if(Verbosity> 1) cout << "# Number of peaks (" << Peaks.size() << ") or number of energies (" << Energies.size() << ") too small" << endl;
return 0;
}
else {
#if 0
// all peaks used as pivot
int npmax = Peaks.size();
for(int np = 0; np < npmax; np++)
Peaks[np].good = true;
#else
// only the peak with largest area used as pivot
sort( Peaks.begin( ), Peaks.end( ), largerAmplitude() );
int npmax = min( min(Peaks.size(), size_t(4)), Energies.size() ); // limit the number of peaks to min(np,ne,4)
for(int np = 0; np < npmax; np++)
Peaks[np].good = true;
for(size_t np = npmax; np < Peaks.size(); np++)
Peaks[np].good = false;
#endif
// (re)order the peaks
sort( Peaks.begin( ), Peaks.end( ), smallerPosi() );
if(Verbosity> 2) {
cout<<Form("# Sorted Peak Positions = (");
for ( Iter = Peaks.begin( ) ; Iter != Peaks.end( ) ; Iter++ ) cout<<Form( " %d", int((*Iter).posi) );
cout<<Form(" )\n");
}
// only the peak with largest area used as pivot
sort( Peaks.begin( ), Peaks.end( ), largerAmplitude() );
int npmax = min( min(Peaks.size(), size_t(4)), Energies.size() ); // limit the number of peaks to min(np,ne,4)
for(int np = 0; np < npmax; np++)
Peaks[np].good = true;
for(size_t np = npmax; np < Peaks.size(); np++)
Peaks[np].good = false;
// (re)order the peaks
sort( Peaks.begin( ), Peaks.end( ), smallerPosi() );
if(Verbosity> 2) {
cout<<Form("# Sorted Peak Positions = (");
for ( Iter = Peaks.begin( ) ; Iter != Peaks.end( ) ; Iter++ ) cout<<Form( " %d", int((*Iter).posi) );
cout<<Form(" )\n");
}
// Connect every (large) peak with every energy and check how many other (peak,energy) pairs agree with this slope
// The combination with the largest number of agreeing pairs provides the reference gain for the final average
// Equal number of matches are ordered by chi2
int np_max = 0; // number of matches
int np_ind = 0; // reference peak
int np_ref = 0; // reference energy
double np_chi = 1.e40; // to distinguish among equal matchers
double dpos2max = specFWHMdef*specFWHMdef;
for(size_t np = 0; np < Peaks.size(); np++) {
if(!Peaks[np].good)
continue;
int ec_max = 0; // number of matches for peak np
int ec_ref = 0; // reference energy
double ec_chi = 1.e30; // its chi2
for(size_t ne = 0; ne < Energies.size(); ne++) { // connecting peak np with energy ne
double lgain = (Peaks[np].posi)/Energies[ne]; // gives this gain
int match = 0;
double chi2 = 0;
for(size_t nee = 0; nee < Energies.size(); nee++) { // count the number of matches for this combination
double epos = Energies[nee]*lgain; // expected position
for(size_t ip = 0; ip < Peaks.size(); ip++) {
double dpos = Peaks[ip].posi-epos;
double dpos2 = dpos*dpos;
if(dpos2 < dpos2max) {
match++;
chi2 += dpos2;
break;
// Connect every (large) peak with every energy and check how many other (peak,energy) pairs agree with this slope
// The combination with the largest number of agreeing pairs provides the reference gain for the final average
// Equal number of matches are ordered by chi2
int np_max = 0; // number of matches
int np_ind = 0; // reference peak
int np_ref = 0; // reference energy
double np_chi = 1.e40; // to distinguish among equal matchers
double dpos2max = specFWHMdef*specFWHMdef;
for(size_t np = 0; np < Peaks.size(); np++) {
if(!Peaks[np].good)
continue;
int ec_max = 0; // number of matches for peak np
int ec_ref = 0; // reference energy
double ec_chi = 1.e30; // its chi2
for(size_t ne = 0; ne < Energies.size(); ne++) { // connecting peak np with energy ne
double lgain = (Peaks[np].posi)/Energies[ne]; // gives this gain
int match = 0;
double chi2 = 0;
for(size_t nee = 0; nee < Energies.size(); nee++) { // count the number of matches for this combination
double epos = Energies[nee]*lgain; // expected position
for(size_t ip = 0; ip < Peaks.size(); ip++) {
double dpos = Peaks[ip].posi-epos;
double dpos2 = dpos*dpos;
if(dpos2 < dpos2max) {
match++;
chi2 += dpos2;
break;
}
}
}
}
if( (match > ec_max) || // more matches
if( (match > ec_max) || // more matches
((match == ec_max) && (chi2 < ec_chi)) ) { // or better chi2
ec_max = match; // record the best combination so far
ec_ref = ne;
ec_chi = chi2;
ec_max = match; // record the best combination so far
ec_ref = ne;
ec_chi = chi2;
}
}
}
if( (ec_max > np_max) || // more matches
if( (ec_max > np_max) || // more matches
((ec_max == np_max) && (ec_chi < np_chi)) ) { // or better chi2
np_max = ec_max; // record the best combination so far
np_ind = np;
np_ref = ec_ref;
np_chi = ec_chi;
np_max = ec_max; // record the best combination so far
np_ind = np;
np_ref = ec_ref;
np_chi = ec_chi;
}
//if(np_max == Energies.size())
// break; // cannot do better than this ?
}
//if(np_max == Energies.size())
// break; // cannot do better than this ?
}
double bestSlope = Energies[np_ref]/(Peaks[np_ind].posi);
if(Verbosity > 2) cout<<Form("# Best-match slope %g [p=%d e=%d] with %d values\n", bestSlope, np_ind, np_ref, np_max);
bestSlope = Energies[np_ref]/(Peaks[np_ind].posi);
if(Verbosity > 2) cout<<Form("# Best-match slope %g [p=%d e=%d] with %d values\n", bestSlope, np_ind, np_ref, np_max);
// find the good peaks
for(size_t np = 0; np < Peaks.size(); np++) {
Peaks[np].erefindex = -1;
Peaks[np].eref = 0.;
Peaks[np].good = false;
}
Float_t Xmin = 0;
Float_t Xmax = 0.;
int match = 0;
for(size_t ne = 0; ne < Energies.size(); ne++) {
double epos = Energies[ne]/bestSlope;
// find the good peaks
for(size_t np = 0; np < Peaks.size(); np++) {
if(abs(Peaks[np].posi-epos) < specFWHMdef && !Peaks[np].good) {
match++;
Peaks[np].erefindex = ne;
Peaks[np].eref = Energies[ne];
Peaks[np].good = true;
// Peaks[np].errposi = pow(Peaks[np].fwhm/s2fwhm, 2.)/Peaks[np].area;
if(Xmin==0. || Peaks[np].posi < Xmin) Xmin = Peaks[np].posi;
if(Peaks[np].posi > Xmax) Xmax = Peaks[np].posi;
break;
Peaks[np].erefindex = -1;
Peaks[np].eref = 0.;
Peaks[np].good = false;
}
int match = 0;
for(size_t ne = 0; ne < Energies.size(); ne++) {
double epos = Energies[ne]/bestSlope;
for(size_t np = 0; np < Peaks.size(); np++) {
if(abs(Peaks[np].posi-epos) < specFWHMdef && !Peaks[np].good) {
match++;
Peaks[np].erefindex = ne;
Peaks[np].eref = Energies[ne];
Peaks[np].good = true;
// Peaks[np].errposi = pow(Peaks[np].fwhm/s2fwhm, 2.)/Peaks[np].area;
if(Xmin==0. || Peaks[np].posi < Xmin) Xmin = Peaks[np].posi;
if(Peaks[np].posi > Xmax) Xmax = Peaks