diff --git a/macros/CCEScan/voltageFit_spline.C b/macros/CCEScan/voltageFit_spline.C index 50d3344..6396699 100755 --- a/macros/CCEScan/voltageFit_spline.C +++ b/macros/CCEScan/voltageFit_spline.C @@ -105,6 +105,8 @@ int val = 3; bool notSave = false; bool temSave = false; + bool use_height = false; + bool average_ratio = false; double radius = -10.0; @@ -130,7 +132,7 @@ } - + std::vector v_volt_val; std::vector v_volt_err; @@ -149,10 +151,7 @@ diskVar,det.Data(),lay.Data(),sector,fill)); TString fn_graph_p(Form("%s/volt_val%d", dn_graph.Data(),val)); - - - - + if (radius<0.0) { }else{ vn_result +=Form("_r%d",(int)radius); @@ -176,17 +175,17 @@ TString dn_data(Form("%s/data/ST/Aging" - ,diskVar)); - TString fn_data(Form("%s/CCEScan.root",//Michele for TESTING + ,diskVar)); + TString fn_data(Form("%s/CCEScan.root", dn_data.Data())); TString fn_dummy(Form("%s/check",dn_data.Data())); // TString dn_ratio(Form("%s/cmtuser/Vetra_v13r2/ST/STAging/data/db" // ,homeVar)); TString dn_ratio(Form("%s/data/db", homeVar)); - TString fn_ratio(Form("%s/ratio.root", - dn_ratio.Data())); - + TString fn_ratio(Form("%s/ratio_int.root", + dn_ratio.Data())); + TString fn_syst(Form("%s/systVolt.root", dn_ratio.Data())); @@ -222,14 +221,9 @@ fill)); - FILE* myfile; - TString myfile_nm(Form("%s/noSystematic_val%d.txt",dn_graph.Data(),val)); - cout << myfile_nm << endl; - myfile =fopen(myfile_nm,"w"); - fprintf(myfile, "List of the calibration steps with no estimated systematic on pulse shape\n\n"); - + + for (int i=0; i<11; i++) { - TString vn_input(Form("%s/v_pulse%d_val%d",dn_input.Data(),i,val)); if (radius>0.0) { vn_input +=Form("_r%d",(int)radius); @@ -263,10 +257,8 @@ double height_val = v_input(0); double height_err = v_input(1); - double height_sys = v_input(11);//Michele + double pulse_quality = v_input(11); - /* Michele - Use maximum value, not integral - //This is a fake if... if (true || det.EqualTo("IT")) { height_val /= 0.1306; // transform back to the amplitude height_err /= 0.1306; // transform back to the amplitude @@ -275,28 +267,24 @@ TMath::Power(v_input(3)/v_input(2),2)+ 2*v_input(8)/(height_val*v_input(2))); } - */ + - Printf(" %8.2f V: %8.4f+/-%8.4f(stat) +/-%8.4f (sys)",volt,height_val,height_err, height_sys);//MIchele - - // Load the input values in to vector - - if (height_sys == -999){ - fprintf(myfile," det=%s lay=%s sect=%d fill=%d val=%d rad=%f\n", det.Data(),lay.Data(),sector,fill,val,radius); - + Printf(" %8.2f V: %8.4f+/-%8.4f",volt,height_val,height_err); + if (pulse_quality == -1){ + cout << "Ignoring this point for the voltageFit..."<Get(Form("v_ratio_%s",lay.Data())); - TVectorD* vp_ratio_err = (TVectorD*)f_ratio->Get(Form("v_ratio_err_%s",lay.Data())); - TVectorD* vp_secRatio = (TVectorD*)f_ratio->Get(Form("v_sector_%s",lay.Data())); + // Extract the ratio values for the sampling point + TVectorD* vp_ratio = (TVectorD*)f_ratio->Get(Form("v_ratio_%s_val%d",lay.Data(),val)); + TVectorD* vp_ratio_err = (TVectorD*)f_ratio->Get(Form("v_ratio_err_%s_val%d",lay.Data(),val)); + TVectorD* vp_secRatio = (TVectorD*)f_ratio->Get(Form("v_sector_%s_val%d",lay.Data(),val)); if (!vp_ratio || !vp_secRatio || !vp_ratio_err) { Error("voltageFit_spline","Ratio vectors not available"); @@ -342,24 +330,41 @@ double ratio_err = STTool::ITratio_err; if (det.EqualTo("TT")) { + ratio_val = STTool::TTratio_val; ratio_err = STTool::TTratio_err; + + /* + if(val == 3){ + ratio_val = STTool::TTratio_val_3; + ratio_err = STTool::TTratio_err_3; + } + else if(val == 5){ + ratio_val = STTool::TTratio_val_5; + ratio_err = STTool::TTratio_err_5; + } + else if(val == 7){ + ratio_val = STTool::TTratio_val_7; + ratio_err = STTool::TTratio_err_7; + } + */ } - // Elena - Cancelled by Michele because of what said in l.159 in the Radiation Damage paper + if (false && det.EqualTo("TT")) { - - // Using ratio per sector - if (true) { - int ratioInd = getClosestIndex(v_secRatio,(double)sector);//Gets the index in v_secRatio that has a double that is closest to sector - if (!(ratioInd<0) && (int)v_secRatio(ratioInd)==sector && true) {//if this is not senseless and the sector is correct - Info("voltageFit_spline","Using specific ratio for sector %d",sector); - ratio_val = v_ratio(ratioInd); - ratio_err = v_ratio_err(ratioInd); - } + // Using ratio per sector + //This is excluded by default + //if (true) { + int ratioInd = getClosestIndex(v_secRatio,(double)sector); + if (!(ratioInd<0) && (int)v_secRatio(ratioInd)==sector && true) { + Info("voltageFit_spline","Using specific ratio for sector %d",sector); + Printf(" Individual ratio: %f",ratio_val); + ratio_val = v_ratio(ratioInd); + ratio_err = v_ratio_err(ratioInd); } } + Printf(" Ratio: %f",ratio_val); // Define fit range @@ -399,7 +404,7 @@ // Scale them Info("voltageFit_spline","#chi^{2}/ndf: %6.2f/%d",fu_volt->GetChisquare(),fu_volt->GetNDF()); - /*Michele, confirmed in 8/31 2017 + /*Michele if (chi2ndf>0.00) { Info("voltageFit_spline","scaling"); vr_adc_err_corr*=TMath::Sqrt(chi2ndf); @@ -433,7 +438,6 @@ int iFirst = -1; int iSecon = -1; // Excluding first and last point as they have been artifically added - //Looking for the two voltages that are closest to Vdepl for(int i=1; iGetParameter("A_{0}")/s5_volt->Derivative(0.0); //Wrong, sometimes this is a bad indication! + double Vdepl_sys1 = Vdepl_val-fu_volt->GetParameter("A_{0}")/s5_volt->Derivative(0.0); // Alternative determination of the systematic uncertainty by taking the // difference between the determined depletion voltage and the extracted value from the sigmoid function @@ -566,6 +570,7 @@ // Load map + /* TFile* f_map = new TFile(fn_syst.Data(),"UPDATE"); TExMap* m_lin = (TExMap*)f_map->Get(Form("m_lin_%s",lay.Data())); // Deviations from linear approach TExMap* m_sig = (TExMap*)f_map->Get(Form("m_sig_%s",lay.Data())); // Deviations from sigmoid approach @@ -586,14 +591,14 @@ m_lin->Write(Form("m_lin_%s",lay.Data())); m_sig->Write(Form("m_sig_%s",lay.Data())); f_map->Close(); - + */ // Save the result v_res(0) = Vdepl_val; v_res(1) = Vdepl_err; v_res(2) = Vdepl_sys; v_res(3) = fu_volt->GetParameter("A_{0}"); v_res(4) = fu_volt->GetParError(fu_volt->GetParNumber("A_{0}")); - + @@ -664,7 +669,7 @@ gDirectory->Delete(Form("%s;*", vn_result.Data())); v_res.Write(Form("%s", - vn_result.Data())); + vn_result.Data())); Info("voltageFit_spline","Result written to %s",fn_data.Data()); } @@ -682,4 +687,3 @@ return EXIT_SUCCESS; } -