C++ TGraph::SetLineColor方法代码示例

void Sens_Matrix_LAr1ND_200m_T600_on_axis_Shape_and_Rate_5sigma(TCanvas* c)
{
   c->cd();
   TGraph *graph = new TGraph(501);
   graph->SetName("Graph0");
   graph->SetTitle("Graph");
   graph->SetFillColor(1);
   graph->SetLineColor(9);
   graph->SetLineStyle(2);
   graph->SetPoint(0,1,0.01);
   graph->SetPoint(1,1,0.01018591388);
   graph->SetPoint(2,1,0.01037528416);
   graph->SetPoint(3,1,0.01056817509);
   graph->SetPoint(4,1,0.01076465214);
   graph->SetPoint(5,1,0.01096478196);
   graph->SetPoint(6,1,0.01116863248);
   graph->SetPoint(7,1,0.01137627286);
   graph->SetPoint(8,1,0.01158777356);
   graph->SetPoint(9,1,0.01180320636);
   graph->SetPoint(10,1,0.01202264435);
   graph->SetPoint(11,1,0.01224616199);
   graph->SetPoint(12,1,0.01247383514);
   graph->SetPoint(13,1,0.01270574105);
   graph->SetPoint(14,1,0.01294195841);
   graph->SetPoint(15,1,0.01318256739);
   graph->SetPoint(16,1,0.01342764961);
   graph->SetPoint(17,1,0.01367728826);
   graph->SetPoint(18,1,0.01393156803);
   graph->SetPoint(19,1,0.01419057522);
   graph->SetPoint(20,1,0.01445439771);
   graph->SetPoint(21,1,0.01472312502);
   graph->SetPoint(22,1,0.01499684836);
   graph->SetPoint(23,1,0.01527566058);
   graph->SetPoint(24,1,0.01555965632);
   graph->SetPoint(25,1,0.01584893192);
   graph->SetPoint(26,1,0.01614358557);
   graph->SetPoint(27,1,0.01644371723);
   graph->SetPoint(28,1,0.01674942876);
   graph->SetPoint(29,1,0.01706082389);
   graph->SetPoint(30,1,0.01737800829);
   graph->SetPoint(31,1,0.01770108958);
   graph->SetPoint(32,1,0.01803017741);
   graph->SetPoint(33,1,0.01836538343);
   graph->SetPoint(34,1,0.0187068214);
   graph->SetPoint(35,1,0.01905460718);
   graph->SetPoint(36,1,0.01940885878);
   graph->SetPoint(37,1,0.0197696964);
   graph->SetPoint(38,1,0.0201372425);
   graph->SetPoint(39,1,0.02051162179);
   graph->SetPoint(40,1,0.02089296131);
   graph->SetPoint(41,1,0.02128139046);
   graph->SetPoint(42,1,0.02167704105);
   graph->SetPoint(43,1,0.02208004733);
   graph->SetPoint(44,1,0.02249054606);
   graph->SetPoint(45,1,0.02290867653);
   graph->SetPoint(46,1,0.02333458062);
   graph->SetPoint(47,1,0.02376840287);
   graph->SetPoint(48,1,0.02421029047);
   graph->SetPoint(49,1,0.02466039337);
   graph->SetPoint(50,1,0.02511886432);
   graph->SetPoint(51,1,0.02558585887);
   graph->SetPoint(52,1,0.0260615355);
   graph->SetPoint(53,1,0.02654605562);
   graph->SetPoint(54,1,0.02703958364);
   graph->SetPoint(55,1,0.02754228703);
   graph->SetPoint(56,1,0.02805433638);
   graph->SetPoint(57,1,0.02857590543);
   graph->SetPoint(58,1,0.02910717118);
   graph->SetPoint(59,1,0.0296483139);
   graph->SetPoint(60,1,0.0301995172);
   graph->SetPoint(61,1,0.03076096815);
   graph->SetPoint(62,1,0.03133285724);
   graph->SetPoint(63,1,0.03191537855);
   graph->SetPoint(64,1,0.03250872974);
   graph->SetPoint(65,1,0.03311311215);
   graph->SetPoint(66,1,0.03372873087);
   graph->SetPoint(67,1,0.03435579479);
   graph->SetPoint(68,1,0.0349945167);
   graph->SetPoint(69,1,0.03564511334);
   graph->SetPoint(70,1,0.03630780548);
   graph->SetPoint(71,1,0.03698281798);
   graph->SetPoint(72,1,0.0376703799);
   graph->SetPoint(73,1,0.03837072455);
   graph->SetPoint(74,1,0.03908408958);
   graph->SetPoint(75,1,0.03981071706);
   graph->SetPoint(76,1,0.04055085354);
   graph->SetPoint(77,1,0.0413047502);
   graph->SetPoint(78,1,0.04207266284);
   graph->SetPoint(79,1,0.04285485204);
   graph->SetPoint(80,1,0.04365158322);
   graph->SetPoint(81,1,0.04446312675);
   graph->SetPoint(82,1,0.04528975799);
   graph->SetPoint(83,1,0.04613175746);
   graph->SetPoint(84,1,0.04698941086);
   graph->SetPoint(85,1,0.04786300923);
   graph->SetPoint(86,1,0.04875284901);
   graph->SetPoint(87,1,0.04965923215);
   graph->SetPoint(88,1,0.0505824662);
   graph->SetPoint(89,1,0.05152286446);
   graph->SetPoint(90,1,0.05248074602);
//.........这里部分代码省略.........

TGraph* GetStat(TString className="CVHMV0M-B-NOPF-CENTNOTRD",Bool_t lumi=1, Bool_t goodOnly=0, TString part="PHYSICS_1"){
  Double_t stat[100000];
  Double_t time_stamp[100000];
  TH1D* hStat      = new TH1D("hStat","",1,0,1);
  TH1D* hTimeStart = new TH1D("hTimeStart","",1,0,1);
  TH1D* hTimeEnd   = new TH1D("hTimeEnd","",1,0,1);
  TString classNameRun;

  for (Int_t r=0;r<t->GetEntries();r++){
    t->GetEntry(r);
    if (!partition->String().Contains(part.Data())) continue;
    if (!lhcState->String().Contains("STABLE")) continue;
    if (!lhcPeriod->String().Contains("LHC15m")) continue;
    if (run<225000) continue;
    hTimeStart->Fill(Form("%i",run),timeStart);
    hTimeEnd->Fill(Form("%i",run),timeEnd);
    classNameRun = TString(className);
    if (className.EqualTo("CEMC7-B-NOPF-CENTNOTRD")                  && run<=236224) classNameRun="CEMC7-ABCE-NOPF-ALLNOTRD";
    if (className.EqualTo("CDMC7-B-NOPF-CENTNOTRD")                  && run<=236224) classNameRun="CDMC7-ABCE-NOPF-ALLNOTRD";
    if (className.EqualTo("CMUL7-B-NOPF-MUFAST")                     && run<=226606) classNameRun="CMUL7-B-NOPF-ALLNOTRD";
    if (className.EqualTo("CVHMV0M-B-NOPF-CENTNOTRD")                && run<=236221) classNameRun="CVHMV0M-B-NOPF-CENT";
    if (className.EqualTo("CVHMSH2-B-NOPF-CENTNOTRD")                && run<=236226) classNameRun="CVHMSH2-B-NOPF-CENT";
    if (className.EqualTo("CVHMV0M-B-NOPF-CENTNOTRD") && run>=238432 && run<=239144) classNameRun="CVHMV0M-B-SPD1-CENTNOTRD"; 
    if (className.EqualTo("CVHMSH2-B-NOPF-CENTNOTRD") && run>=238432               ) classNameRun="CVHMSH2-B-SPD1-CENTNOTRD";
    if (className.EqualTo("CINT7-B-NOPF-CENT")        && run>=225000 && run<=228935) classNameRun="CINT7-B-NOPF-ALLNOTRD";
    if (className.EqualTo("CINT7-B-NOPF-CENT")        && run>=228936 && run<=229893) classNameRun="CINT7-B-NOPF-CENTNOTRD";
    if (className.EqualTo("CINT7-B-NOPF-CENT")        && run>=229894 && run<=229899) classNameRun="CINT7-B-NOPF-ALLNOTRD";
    if (className.EqualTo("CINT7-B-NOPF-CENT")        && run>=229900 && run<=233911) classNameRun="CINT7-B-NOPF-CENT";
    if (className.EqualTo("CINT7-B-NOPF-CENT")        && run>=233912 && run<=234050) classNameRun="CINT7-B-NOPF-ALLNOTRD";
    if (className.EqualTo("CINT7-B-NOPF-CENT")        && run>=238890 && run<=239144) classNameRun="CINT7-I-NOPF-CENTNOTRD";
    if (className.EqualTo("CEMC7-B-NOPF-CENTNOTRD")                  && run<=236224) classNameRun="CEMC7-ABCE-NOPF-ALLNOTRD";
    if (className.EqualTo("CDMC7-B-NOPF-CENTNOTRD")                  && run<=236224) classNameRun="CDMC7-ABCE-NOPF-ALLNOTRD";

    AliTriggerClass* cl = (AliTriggerClass*) classes->FindObject(classNameRun.Data());
    if (!cl) { 
      hStat->Fill(Form("%i",run),0.); 
      continue; 
    }
    
    if (goodOnly){
      if (classNameRun.Contains("CINT7") || classNameRun.Contains("CVHM")){
        Int_t good=0;
        for (Int_t i=0;i<nGoodRuns;i++) good|= (goodRuns[i]==run);
        if (!good) {
          hStat->Fill(Form("%i",run),0.);
          continue;
        }
      }

      if (className.Contains("MC7")){
        if (run<235709 || run==236855 || run==236858 || run==236861 || fill==4440) {
          hStat->Fill(Form("%i",run),0.);
          continue;
        }
      }

      if (classNameRun.Contains("CMUL")){
        if (!activeDetectors->String().Contains("MUONTRK") ||
            !activeDetectors->String().Contains("MUONTRG") ||
            !activeDetectors->String().Contains("T0") ||
            !activeDetectors->String().Contains("VZERO") ||
            !activeDetectors->String().Contains("ITSSPD")
        ){
          hStat->Fill(Form("%i",run),0.);
          continue;
        }
      }
    }
//    printf("%i\n",classes->IndexOf(cl));
    Double_t x = lumi ? class_lumi[classes->IndexOf(cl)] : class_l2a[classes->IndexOf(cl)];
    hStat->Fill(Form("%i",run),x);
  }
  hStat->GetXaxis()->LabelsOption("a");
  hTimeEnd->GetXaxis()->LabelsOption("a");
  hTimeStart->GetXaxis()->LabelsOption("a");
  hStat->LabelsDeflate("x");
  hTimeStart->LabelsDeflate("x");
  hTimeEnd->LabelsDeflate("x");
  
  Int_t n = 2*hStat->GetNbinsX();

  for (Int_t i=0;i<hStat->GetNbinsX();i++){
    time_stamp[2*i  ] = hTimeStart->GetBinContent(i+1);
    time_stamp[2*i+1] = hTimeEnd->GetBinContent(i+1);
    stat[2*i]         = i>0?stat[2*i-1]:0;
    stat[2*i+1]       = stat[2*i]+hStat->GetBinContent(i+1)/1000000.;
  }

  delete hStat;
  delete hTimeStart;
  delete hTimeEnd;
  TGraph* gStat = new TGraph(n,time_stamp,stat);
  gStat->SetLineColor(kBlue);
  gStat->SetLineWidth(2);
  return gStat;
}

//.........这里部分代码省略.........
  Double_t *x1bins = new Double_t[NP1X1Bins];
  Double_t *eRms   = new Double_t[NP1X1Bins];
  for(Int_t i=1;i<=hP1X1prof->GetNbinsX();i++) {
    x1bins[i] = hP1X1prof->GetBinCenter(i);
    eRms[i] = 100 * hP1X1prof->GetBinError(i) / hP1X1prof->GetBinContent(i);
  }
  gErmsB = new TGraph(NP1X1Bins,x1bins,eRms);
  
  // Vertical Energy histogram:
  // --------------------------------------------------------------------------------   
  TGraph *gP1left = NULL;
  if(hP1) {
    Double_t *yarray   = new Double_t[yNbin];
    Double_t *xarray   = new Double_t[yNbin];
    
    // This is for the right side:
    // Double_t xMax = x1Min + (x1Max-x1Min) * 0.9;
    // Double_t xMin = x1Max;
    // And this for left:
    Double_t xMin = hX1->GetXaxis()->GetXmin();
    Double_t xMax = hX1->GetXaxis()->GetXmin() + (hX1->GetXaxis()->GetXmax()
						  -hX1->GetXaxis()->GetXmin()) * 0.2;
    Double_t EneMax = hP1->GetMaximum();
    // cout << Form("  EneMax = %f ", EneMax) << endl;
 
    for(Int_t j=0; j<yNbin; j++) {
      yarray[j] = hP1->GetBinCenter(j+1);
      xarray[j] = ((xMax-xMin)/EneMax)*hP1->GetBinContent(j+1) + xMin;

      // cout << Form("  x = %f  y = %f ", xarray[j],yarray[j]) << endl;
    }

    gP1left = new TGraph(yNbin,xarray,yarray);
    gP1left->SetLineColor(PlasmaGlob::elecLine);
    gP1left->SetLineWidth(2);
    gP1left->SetFillStyle(1001);
    gP1left->SetFillColor(PlasmaGlob::elecFill);
       
  }

  

  // Plotting
  // -----------------------------------------------
    
  // Canvas setup
  // Create the canvas and the pads before the Frame loop
  // Resolution:
  Int_t sizex = 800;
  Int_t sizey = 600;
  if(opt.Contains("hres")) {
    Int_t sizex = 1600;
    Int_t sizey = 1200;    
  }
  
  TCanvas *C = new TCanvas("C1","Evolution of Injection",sizex,sizey);
  C->cd();

  // Set palette:
  PPalette * pPalette = (PPalette*) gROOT->FindObject("electron");
  pPalette->cd();

  // Float_t Max  = hP1X1->GetMaximum();
  // Float_t Min  = hP1X1->GetMinimum();
  
  // hP1X1->GetZaxis()->SetRangeUser(Min,Max); 

//____________________________________________________________________________________
// Plot blue band (SM Higgs curve)
void blueBand( )
{

  double zerror = 4;
  double xmin = 76;

  int npoints = 31;
  double x[] = { 10.000000, 11.392849, 12.979700, 14.787576, 16.847262, 19.193831, 21.867241, 
                 24.913018, 28.383024, 32.336350, 36.840315, 41.971614, 47.817625, 54.477897,
                 62.065844, 70.710678, 76, 80.559606, 91.780341, 104.563955, 119.128133, 135.720881,
                 154.624747, 176.161637, 200.698289, 228.652526, 260.500365, 296.784127, 
                 338.121669, 385.216905, 438.871795, 500.000000 };

  double y[] = { 20.929260, 18.757168, 16.647426, 14.611732, 12.661571, 10.808107, 9.062193,
                 7.434501, 5.935701, 4.576664, 3.368675, 2.323612, 1.454084, 0.776111,
                 0.300018, 0.042673, 0.0001, 0.020469, 0.251087, 0.750519, 1.538231, 2.633789,
                 4.056692, 5.826546, 7.962673, 10.483595, 13.406320, 16.745357, 20.511276,
                 24.717583, 29.345917, 34.380810 };


  TGraph* smGraph = new TGraph( npoints, x, y );

  if ( npoints > MAXPOINTS ) {
    std::cerr << "FATAL: npoints = " << npoints << " is greater than MAXPOINTS" << std::endl;
    exit(-1);
  }
  double x1[2*MAXPOINTS], y1[2*MAXPOINTS];
  double x2[MAXPOINTS], y2[MAXPOINTS];

  // Find mH at minimum
  int i1=0; // iterator for upper curve
  int i2=0; // iterator for lower curve
  double x0,y0;
  double step = 5;
  smGraph->GetPoint( 0, x0, y0 );
  for ( int i=0; i<100; ++i )
    {
      double ix = x0+i*step;
      double iy = smGraph->Eval( ix, 0, "S" );
      double error = TMath::Log( ix/12. )*3.;
      if ( ix<=xmin-error ) { 
        x1[i1] = ix+error; y1[i1]=iy; ++i1; 
        x2[i2] = ix-error; y2[i2]=iy; ++i2;
      }
      else if ( ix>=xmin+error) { 
        x1[i1] = ix-error; y1[i1]=iy; ++i1; 
        x2[i2] = ix+error; y2[i2]=iy; ++i2;
      } else {//if ( iy>0. ) {
        x2[i2] = (ix<=xmin)?ix-error:ix+error; 
        y2[i2]=iy; 
        ++i2;
      }
      
    }
  // Merge arrays
  for ( int i=0; i<i2; ++i ) {
    x1[i+i1] = x2[i2-i-1];
    y1[i+i1] = y2[i2-i-1];
  }
  if ( doBands ) {
    TGraph* blueBand = new TGraph( i1+i2, x1, y1 );
    blueBand->SetFillColor(7);
    blueBand->SetLineColor(7);
    blueBand->Draw("LF");
    blueBand->Draw("C");
  } else {
    smGraph->SetLineWidth(5);
  }

  smGraph->SetLineColor(4);
  smGraph->Draw("C");
  double chi2limit = smGraph->Eval( LEPLIMIT, 0, "S" );
  std::cout << "Value at limit (" << LEPLIMIT << "): " << chi2limit << std::endl;

  TGraph* lowband = new TGraph( i2, x2, y2 );
  std::cout << "Value of (low) blueband at limit: " << lowband->Eval( LEPLIMIT, 0, "S" )
            << std::endl;

}

TCanvas *timeonaxis()
{

//Author:

   time_t script_time;
   script_time = time(0);
   script_time = 3600*(int)(script_time/3600);

// The time offset is the one that will be used by all graphs.
// If one changes it, it will be changed even on the graphs already defined
   gStyle->SetTimeOffset(script_time);

   TCanvas *ct = new TCanvas("ct","Time on axis",10,10,700,900);
   ct->Divide(1,3);
   ct->SetFillColor(28);

   int i;

//======= Build a signal : noisy damped sine ======
//        Time interval : 30 minutes

   gStyle->SetTitleH(0.08);
   float noise;
   TH1F *ht = new TH1F("ht","Love at first sight",3000,0.,2000.);
   for (i=1;i<3000;i++) {
      noise = gRandom->Gaus(0,120);
      if (i>700) {
         noise += 1000*sin((i-700)*6.28/30)*exp((double)(700-i)/300);
      }
      ht->SetBinContent(i,noise);
   }
   ct->cd(1);
   gPad->SetFillColor(41);
   gPad->SetFrameFillColor(33);
   ht->SetLineColor(2);
   ht->GetXaxis()->SetLabelSize(0.05);
   ht->Draw();
// Sets time on the X axis
// The time used is the one set as time offset added to the value
// of the axis. This is converted into day/month/year hour:min:sec and
// a reasonable tick interval value is chosen.
   ht->GetXaxis()->SetTimeDisplay(1);

//======= Build a simple graph beginning at a different time ======
//        Time interval : 5 seconds

   float x[100], t[100];
   for (i=0;i<100;i++) {
      x[i] = sin(i*4*3.1415926/50)*exp(-(double)i/20);
      t[i] = 6000+(double)i/20;
   }
   TGraph *gt = new TGraph(100,t,x);
   gt->SetTitle("Politics");
   ct->cd(2);
   gPad->SetFillColor(41);
   gPad->SetFrameFillColor(33);
   gt->SetFillColor(19);
   gt->SetLineColor(5);
   gt->SetLineWidth(2);
   gt->Draw("AL");
   gt->GetXaxis()->SetLabelSize(0.05);
// Sets time on the X axis
   gt->GetXaxis()->SetTimeDisplay(1);
   gPad->Modified();

//======= Build a second simple graph for a very long time interval ======
//        Time interval : a few years

   float x2[10], t2[10];
   for (i=0;i<10;i++) {
      x2[i] = gRandom->Gaus(500,100)*i;
      t2[i] = i*365*86400;
   }
   TGraph *gt2 = new TGraph(10,t2,x2);
   gt2->SetTitle("Number of monkeys on the moon");
   ct->cd(3);
   gPad->SetFillColor(41);
   gPad->SetFrameFillColor(33);
   gt2->SetFillColor(19);
   gt2->SetMarkerColor(4);
   gt2->SetMarkerStyle(29);
   gt2->SetMarkerSize(1.3);
   gt2->Draw("AP");
   gt2->GetXaxis()->SetLabelSize(0.05);
// Sets time on the X axis
   gt2->GetXaxis()->SetTimeDisplay(1);
//
// One can choose a different time format than the one chosen by default
// The time format is the same as the one of the C strftime() function
// It's a string containing the following formats :
//    for date :
//      %a abbreviated weekday name
//      %b abbreviated month name
//      %d day of the month (01-31)
//      %m month (01-12)
//      %y year without century
//      %Y year with century
//
//    for time :
//.........这里部分代码省略.........

/**
The purpose of this function is to conveniently plot an event. Each event is plotted on a TCanvas. Divided into as many pads as there are panels. Each pad is a TMultiGraph with a TLegend, and shows the channels of the digitizer connected to the lines of that panel. The grouping into pads must not necessarily be according to panels, but in any other prefered way.

@param a_channels - a vector of a vector of samples, containing all 32 channels
@ param a_channelsToPadsAssociation - a map from std::string, which is the name of the panel (or group of channels) to be assigned to each pad to a vector if integers, which is the list of channels indices corresponding to indices of channels in the paramater a_channels to associate to a pad
@param sEventTitle - a string containing the title of the event (for example, the time stamp)
*/
void RangePlotter::PlotRanges(Channels_t& a_channels, Range_t& a_channelsToPadsAssociation, std::string sEventTitle)
{	
//	printf("Plotting\n");
	//m_pCanvas->Clear();
	m_pCanvas->SetTitle(sEventTitle.c_str());
	
	
	int iPadCounter = 0;
	if(0 == m_vpMultiGraph.size())
	{
		MakePads(a_channelsToPadsAssociation.size());	
		for (auto& rangeIt: a_channelsToPadsAssociation)
		{	
			TMultiGraph* pMg = new TMultiGraph();	
			m_vpMultiGraph.push_back(std::unique_ptr<TMultiGraph>(pMg));

	//		m_pCanvas->cd(iPadCounter + 1);
			ChangePad(iPadCounter);
			int i = 0;

			auto legend = new TLegend(0.8,0.8,1,1, "Channels");
			m_vpLegends.push_back(std::unique_ptr<TLegend>(legend));
			printf("printing panel %s\n", rangeIt.first.c_str());
			for (auto& chanIt: rangeIt.second)
			{
				int iNumOfSamples = a_channels[chanIt].size();
				TGraph* pGr = new TGraph(iNumOfSamples);
				std::vector<float> vTimeSeq = CommonUtils::GenerateTimeSequence(iNumOfSamples, m_fSamplingFreqGHz);
				for (int counter = 0; counter < iNumOfSamples; counter++)
				{
					pGr->SetPoint(counter, vTimeSeq[counter], TransformToVoltage(a_channels[chanIt][counter]));
				}

				m_vpGraph[chanIt] = pGr;
				pGr->SetLineColor(m_colors[i%(sizeof(m_colors)/sizeof(int))]);
				pGr->SetName((m_sInstanceName + std::string("Pan_") + rangeIt.first + std::string("chan_") + std::to_string(chanIt)).c_str());
				std::string sGraphTitle = std::string("Channel ") + std::to_string(chanIt);
				pGr->SetTitle(sGraphTitle.c_str());
				legend->AddEntry(pGr,std::to_string(chanIt).c_str(), "l");
		
				pMg->Add(pGr);
				i++;
			}

			if(Configuration::Instance().ShowTriggerInWaveformsStep())
			{
				int iNumOfSamples = a_channels[a_channels.size() - 1].size();
				m_vpGraphPrecisionTrigger = new TGraph(iNumOfSamples);
				std::vector<float> vTimeSeq = CommonUtils::GenerateTimeSequence(iNumOfSamples, m_fSamplingFreqGHz);
				for (int counter = 0; counter < iNumOfSamples; counter++)
				{
					m_vpGraphPrecisionTrigger->SetPoint(counter, vTimeSeq[counter], TransformToVoltage(a_channels[a_channels.size() - 1][counter]));
				}
				
				m_vpGraphPrecisionTrigger->SetName((m_sInstanceName + std::string("Pan_") + rangeIt.first + "_Trig").c_str());	
				m_vpGraphPrecisionTrigger->SetTitle("Trigger");
				legend->AddEntry(m_vpGraphPrecisionTrigger,"Trigger", "l");
				pMg->Add(m_vpGraphPrecisionTrigger);
			}
			
			std::string sMultiGraphTitle = std::string("Panel ") + rangeIt.first;
			pMg->SetTitle(sMultiGraphTitle.c_str());

			pMg->Draw("AL");
			pMg->GetXaxis()->SetTitle("Time [nanoseconds]");
			pMg->GetXaxis()->CenterTitle();
			pMg->GetYaxis()->SetTitle("Voltage [volts]");
			pMg->GetYaxis()->CenterTitle();

			gPad->Modified();

			pMg->SetMinimum(m_fMinVoltage);
			pMg->SetMaximum(m_fMaxVoltage);
			legend->Draw();
			iPadCounter++;
		}
		m_pCanvas->Update();
	}
	else
	{
//		printf("Plottin again\n");
		for (auto& rangeIt: a_channelsToPadsAssociation)
		{
			printf("Panel %s\n", rangeIt.first.c_str());
			m_pCanvas->cd(iPadCounter + 1);
			for (auto& chanIt: rangeIt.second)
			{
//				printf("Chanenl %d\n", chanIt);
				//TODO: num of samples is constant per run at least!
				m_vpGraph[chanIt]->SetLineWidth(1);
				int iNumOfSamples = a_channels[chanIt].size();	
				std::vector<float> vTimeSeq = CommonUtils::GenerateTimeSequence(iNumOfSamples, m_fSamplingFreqGHz);	
				for (int counter = 0; counter < iNumOfSamples; counter++)
//.........这里部分代码省略.........

// Function for the computation of channeling efficiency at various incoming angle
Int_t AnalyseChannelingEfficiency(TTree *fTree,Float_t fChannelingMinimum = 35., Float_t fChannelingMaximum = 70.){
    //**//Channeling Gaussian Fit Function
    TF1 *vChanneling = new TF1("vChanneling","gaus",fChannelingMinimum,fChannelingMaximum);
    vChanneling->SetParNames("Const","Mean","Sigma");
    vChanneling->SetLineColor(4);
    vChanneling->SetLineStyle(2);

    TH2D *hChannelingPlot = new TH2D("hChannelingPlot","Deflection Angle vs. Incoming Angle;Horizontal Incoming Angle [#murad];Horizontal Deflection Angle [#murad]",21,-10.5,10.5,256,-127.5,128.5);
    
    TH1F *hChannelingEfficiency = new TH1F("hChannelingEfficiency","G4Channeling;Horizontal Incoming Angle [#murad];Efficiency [%]",21,-10.5,10.5);

    fTree->Draw("-(angXout-angXin):-angXin>>hChannelingPlot");
    
    Double_t vNormalizationToAmorphous = 0.965; // Normalization for channeling efficiency, see PRSTAB 11, 063501 (2008)
    
    for(int i=2;i<=21;i++){
        TH1D* h1 = hChannelingPlot->ProjectionY("h1",i,i);
        h1->Fit(vChanneling,"QR");
        Double_t *vChannelingParameters;
        vChannelingParameters = vChanneling->GetParameters();
        hChannelingEfficiency->SetBinContent(i,ComputeEfficiency(h1,vChannelingParameters)/vNormalizationToAmorphous);
        h1->Delete();
    }
    hChannelingEfficiency->SetLineColor(3);
    hChannelingEfficiency->SetLineStyle(4);
    hChannelingEfficiency->SetMarkerColor(3);
    hChannelingEfficiency->SetFillStyle(0);
    hChannelingEfficiency->SetMarkerStyle(20);
    hChannelingEfficiency->Draw("PL");

    TGraph* gRoughExperimentalData = new TGraph(11);
    gRoughExperimentalData->SetPoint(	0	,	-10	,	20	);
    gRoughExperimentalData->SetPoint(	1	,	-8	,	38	);
    gRoughExperimentalData->SetPoint(	2	,	-6	,	56	);
    gRoughExperimentalData->SetPoint(	3	,	-4	,	72	);
    gRoughExperimentalData->SetPoint(	4	,	-2	,	80	);
    gRoughExperimentalData->SetPoint(	5	,	0	,	84	);
    gRoughExperimentalData->SetPoint(	6	,	2	,	82	);
    gRoughExperimentalData->SetPoint(	7	,	4	,	78	);
    gRoughExperimentalData->SetPoint(	8	,	6	,	66	);
    gRoughExperimentalData->SetPoint(	9	,	8	,	52	);
    gRoughExperimentalData->SetPoint(	10	,	10	,	37	);

    gRoughExperimentalData->SetLineColor(4);
    gRoughExperimentalData->SetLineStyle(3);
    gRoughExperimentalData->SetFillStyle(0);
    gRoughExperimentalData->SetFillColor(0);
    gRoughExperimentalData->SetMarkerColor(4);
    gRoughExperimentalData->SetMarkerStyle(21);
    gRoughExperimentalData->SetTitle("Phys. Lett. B 680, 129");

    gRoughExperimentalData->Draw("sameCP");

    TLegend *aLegend = new TLegend(0.30,0.15,0.55,0.3);
    aLegend->AddEntry(hChannelingEfficiency);
    aLegend->AddEntry(gRoughExperimentalData);
    aLegend->SetFillStyle(0);
    aLegend->SetLineColor(0);
    aLegend->Draw();

    return 0;
}

//.........这里部分代码省略.........
   exp1[7] = 0.565467804363;
   exp1[22] = 0.425115372019;
   exp2[7] = 0.672064509509;
   exp2[22] = 0.4120593318;
   exp[8] = 0.563947415859;
   obs[8] = 0.536922403468;
   exp1[8] = 0.672650081221;
   exp1[21] = 0.534487182571;
   exp2[8] = 0.808457381294;
   exp2[21] = 0.514649041607;
   exp[9] = 0.700579859433;
   obs[9] = 0.676280384249;
   exp1[9] = 0.798482092783;
   exp1[20] = 0.674378686191;
   exp2[9] = 1.0120004248;
   exp2[20] = 0.663954563503;
   exp[10] = 0.918970653344;
   obs[10] = 0.889845853577;
   exp1[10] = 1.04923440155;
   exp1[19] = 0.858375432352;
   exp2[10] = 1.24717913491;
   exp2[19] = 0.746731022724;
   exp[11] = 1.19793401961;
   obs[11] = 1.17653186674;
   exp1[11] = 1.39512272299;
   exp1[18] = 1.16185610477;
   exp2[11] = 1.62812001722;
   exp2[18] = 1.01020656441;
   exp[12] = 1.58254134376;
   obs[12] = 1.56834374153;
   exp1[12] = 1.84607923047;
   exp1[17] = 1.51628586669;
   exp2[12] = 2.18103685499;
   exp2[17] = 1.31304946585;
   exp[13] = 2.30981854338;
   obs[13] = 2.30117260535;
   exp1[13] = 2.42189837009;
   exp1[16] = 2.13395645197;
   exp2[13] = 2.92307163497;
   exp2[16] = 2.04337780276;
   exp[14] = 3.19375101546;
   obs[14] = 3.16458908761;
   exp1[14] = 3.28113723449;
   exp1[15] = 2.99139611866;
   exp2[14] = 4.07139223438;
   exp2[15] = 2.84407249777;

  TGraph* grExp = new TGraph(15, masses, exp);
  TGraph* grObs = new TGraph(15, masses, obs);
  TGraph* grExp1 = new TGraph(30, masses2D, exp1);
  TGraph* grExp2 = new TGraph(30, masses2D, exp2);

  grExp->SetLineColor(kBlack);
  grExp->SetMarkerColor(kBlack);
  grExp->SetLineStyle(kDashed);
  grObs->SetLineColor(kBlack);
  grObs->SetMarkerColor(kBlack);
  grObs->SetMarkerStyle(20);
  grExp1->SetFillStyle(1001);
  grExp1->SetFillColor(kGreen);
  grExp2->SetFillStyle(1001);
  grExp2->SetFillColor(kYellow);

  mg->Add(grExp2, "F");
  mg->Add(grExp1, "F");
  mg->Add(grExp,  "L");
  mg->Add(grObs,  "P");
  mg->SetMinimum(0.1);
  mg->SetMaximum(10.);
  mg->Draw("a");

  TLine* line = new TLine(200, 1, 1500, 1);
  line->Draw();

  TLegend* leg = new TLegend(0.2, 0.62, 0.58, 0.92);
  leg->SetTextSize(0.05);
  leg->SetTextFont(42);
  leg->AddEntry(grObs, "Obs. Limit", "p");
  leg->AddEntry(grExp, "Exp. Limit", "l");
  leg->AddEntry(grExp1, "Exp. #pm 1#sigma", "f");
  leg->AddEntry(grExp2, "Exp. #pm 2#sigma", "f");
  leg->AddEntry(line, "W'_{SSM}#rightarrowWZ Coupling", "l");
  leg->SetBorderSize(0);
  leg->SetFillColor(0);
  leg->Draw();

  TLatex latexLabel;
  latexLabel.SetNDC();
  latexLabel.SetTextSize(0.05);
  latexLabel.SetTextFont(42);
  latexLabel.DrawLatex(0.33, 0.96, "CMS Preliminary 2011");
  latexLabel.DrawLatex(0.67, 0.38, "#sqrt{s} = 7 TeV");
  latexLabel.DrawLatex(0.62, 0.25, Form("#intL dt = %.2f fb^{-1}",4.98));


  c1->RedrawAxis();
  c1->Print("wprimewz_limits.C");
  c1->Print("wprimewz_limits.pdf");
  c1->Print("wprimewz_limits.png");
}

void BackgroundPrediction(std::string pname,int rebin_factor,int model_number = 0,int imass=750, bool plotBands = false)
{
    rebin = rebin_factor;
    std::string fname = std::string("../fitFilesMETPT34/") + pname + std::string("/histos_bkg.root");
    
    stringstream iimass ;
    iimass << imass;
    std::string dirName = "info_"+iimass.str()+"_"+pname;
    
    
    gStyle->SetOptStat(000000000);
    gStyle->SetPadGridX(0);
    gStyle->SetPadGridY(0);
    
    setTDRStyle();
    gStyle->SetPadGridX(0);
    gStyle->SetPadGridY(0);
    gStyle->SetOptStat(0000);
    
    writeExtraText = true;       // if extra text
    extraText  = "Preliminary";  // default extra text is "Preliminary"
    lumi_13TeV  = "2.7 fb^{-1}"; // default is "19.7 fb^{-1}"
    lumi_7TeV  = "4.9 fb^{-1}";  // default is "5.1 fb^{-1}"
    
    
    double ratio_tau=-1;
    
    TFile *f=new TFile(fname.c_str());
    TH1F *h_mX_CR_tau=(TH1F*)f->Get("distribs_18_10_1")->Clone("CR_tau");
    TH1F *h_mX_SR=(TH1F*)f->Get("distribs_18_10_0")->Clone("The_SR");
    double maxdata = h_mX_SR->GetMaximum();
    double nEventsSR = h_mX_SR->Integral(600,4000);
    ratio_tau=(h_mX_SR->GetSumOfWeights()/(h_mX_CR_tau->GetSumOfWeights()));
    //double nEventsSR = h_mX_SR->Integral(600,4000);
    
    std::cout<<"ratio tau "<<ratio_tau<<std::endl;
    
    TH1F *h_SR_Prediction;
    TH1F *h_SR_Prediction2;
    
    if(blind) {
        h_SR_Prediction2 = (TH1F*)h_mX_CR_tau->Clone("h_SR_Prediction2");
        h_mX_CR_tau->Rebin(rebin);
        h_mX_CR_tau->SetLineColor(kBlack);
        h_SR_Prediction=(TH1F*)h_mX_CR_tau->Clone("h_SR_Prediction");
    } else {
        h_SR_Prediction2=(TH1F*)h_mX_SR->Clone("h_SR_Prediction2");
        h_mX_SR->Rebin(rebin);
        h_mX_SR->SetLineColor(kBlack);
        h_SR_Prediction=(TH1F*)h_mX_SR->Clone("h_SR_Prediction");
        
    }
    h_SR_Prediction->SetMarkerSize(0.7);
    h_SR_Prediction->GetYaxis()->SetTitleOffset(1.2);
    h_SR_Prediction->Sumw2();
    
    /*TFile *f_sig = new TFile((dirName+"/w_signal_"+iimass.str()+".root").c_str());
    RooWorkspace* xf_sig = (RooWorkspace*)f_sig->Get("Vg");
    RooAbsPdf *xf_sig_pdf = (RooAbsPdf *)xf_sig->pdf((std::string("signal_fixed_")+pname).c_str());
    
    RooWorkspace w_sig("w");
    w_sig.import(*xf_sig_pdf,RooFit::RenameVariable((std::string("signal_fixed_")+pname).c_str(),(std::string("signal_fixed_")+pname+std::string("low")).c_str()),RooFit::RenameAllVariablesExcept("low","x"));
    xf_sig_pdf = w_sig.pdf((std::string("signal_fixed_")+pname+std::string("low")).c_str());
   
    RooArgSet* biasVars = xf_sig_pdf->getVariables();
    TIterator *it = biasVars->createIterator();
    RooRealVar* var = (RooRealVar*)it->Next();
    while (var) {
        var->setConstant(kTRUE);
        var = (RooRealVar*)it->Next();
    }
    */
    RooRealVar x("x", "m_{X} (GeV)", SR_lo, SR_hi);
    
    RooRealVar nBackground((std::string("bg_")+pname+std::string("_norm")).c_str(),"nbkg",h_mX_SR->GetSumOfWeights());
    RooRealVar nBackground2((std::string("alt_bg_")+pname+std::string("_norm")).c_str(),"nbkg",h_mX_SR->GetSumOfWeights());
    std::string blah = pname;
    //pname=""; //Antibtag=tag to constrain b-tag to the anti-btag shape
    
    
    /* RooRealVar bg_p0((std::string("bg_p0_")+pname).c_str(), "bg_p0", 4.2, 0, 200.);
     RooRealVar bg_p1((std::string("bg_p1_")+pname).c_str(), "bg_p1", 4.5, 0, 300.);
     RooRealVar bg_p2((std::string("bg_p2_")+pname).c_str(), "bg_p2", 0.000047, 0, 10.1);
     RooGenericPdf bg_pure = RooGenericPdf((std::string("bg_pure_")+blah).c_str(),"(pow([email protected]/13000,@1)/pow(@0/13000,@[email protected]*log(@0/13000)))",RooArgList(x,bg_p0,bg_p1,bg_p2));
   */
    RooRealVar bg_p0((std::string("bg_p0_")+pname).c_str(), "bg_p0", 0., -1000, 200.);
    RooRealVar bg_p1((std::string("bg_p1_")+pname).c_str(), "bg_p1", -13, -1000, 1000.);
    RooRealVar bg_p2((std::string("bg_p2_")+pname).c_str(), "bg_p2", -1.4, -1000, 1000.);
    bg_p0.setConstant(kTRUE);
    //RooGenericPdf bg_pure = RooGenericPdf((std::string("bg_pure_")+blah).c_str(),"(pow(@0/13000,@[email protected]*log(@0/13000)))",RooArgList(x,bg_p1,bg_p2));
    RooGenericPdf bg = RooGenericPdf((std::string("bg_")+blah).c_str(),"(pow(@0/13000,@[email protected]*log(@0/13000)))",RooArgList(x,bg_p1,bg_p2));
  

    /*TF1* biasFunc = new TF1("biasFunc","(0.63*x/1000-1.45)",1350,3600);
    TF1* biasFunc2 = new TF1("biasFunc2","TMath::Min(2.,2.3*x/1000-3.8)",1350,3600);
    double bias_term_s = 0;
    if ((imass > 2450 && blah == "antibtag") || (imass > 1640 && blah == "btag")) {
        if (blah == "antibtag") {
            bias_term_s = 2.7*biasFunc->Eval(imass);
        } else {
//.........这里部分代码省略.........

//.........这里部分代码省略.........
   381.3793,
   381.9813,
   382.5189,
   382.9921,
   383.401,
   383.7455,
   384.0256,
   384.2414,
   384.3928,
   384.4798,
   384.5024,
   384.4607,
   384.3546,
   384.1841,
   383.9493,
   383.65,
   383.2864,
   382.8585,
   382.3661,
   381.8094,
   381.1884,
   380.5029,
   379.7531,
   378.9389,
   378.9389};
   TGraph *graph = new TGraph(102,curvetot_fx1,curvetot_fy1);
   graph->SetName("curvetot");
   graph->SetTitle("Projection of totalPdf");
   graph->SetFillColor(1);

   Int_t ci;      // for color index setting
   TColor *color; // for color definition with alpha
   ci = TColor::GetColor("#0000ff");
   graph->SetLineColor(ci);
   graph->SetLineWidth(3);
   graph->SetMarkerStyle(20);
   
   TH1F *Graph_curvetot1 = new TH1F("Graph_curvetot1","Projection of totalPdf",102,306,474);
   Graph_curvetot1->SetMinimum(117.439);
   Graph_curvetot1->SetMaximum(408.7809);
   Graph_curvetot1->SetDirectory(0);
   Graph_curvetot1->SetStats(0);
   Graph_curvetot1->SetLineWidth(2);
   Graph_curvetot1->SetMarkerStyle(20);
   Graph_curvetot1->GetXaxis()->SetNdivisions(505);
   Graph_curvetot1->GetXaxis()->SetLabelFont(132);
   Graph_curvetot1->GetXaxis()->SetLabelOffset(0.01);
   Graph_curvetot1->GetXaxis()->SetLabelSize(0.06);
   Graph_curvetot1->GetXaxis()->SetTitleSize(0.072);
   Graph_curvetot1->GetXaxis()->SetTitleOffset(0.95);
   Graph_curvetot1->GetXaxis()->SetTitleFont(132);
   Graph_curvetot1->GetYaxis()->SetLabelFont(132);
   Graph_curvetot1->GetYaxis()->SetLabelOffset(0.01);
   Graph_curvetot1->GetYaxis()->SetLabelSize(0.06);
   Graph_curvetot1->GetYaxis()->SetTitleSize(0.072);
   Graph_curvetot1->GetYaxis()->SetTitleOffset(0.95);
   Graph_curvetot1->GetYaxis()->SetTitleFont(132);
   Graph_curvetot1->GetZaxis()->SetLabelFont(132);
   Graph_curvetot1->GetZaxis()->SetLabelSize(0.06);
   Graph_curvetot1->GetZaxis()->SetTitleSize(0.072);
   Graph_curvetot1->GetZaxis()->SetTitleOffset(1.2);
   Graph_curvetot1->GetZaxis()->SetTitleFont(132);
   graph->SetHistogram(Graph_curvetot1);
   
   graph->Draw("l");
   

void draw_width() {


    TCanvas *c1 = new TCanvas("c1","Heaviest H width, M_{A} = 15 GeV",200,10,700,500);
    c1->SetLogy();
    c1->SetFillColor(0);
    c1->SetGrid();

    const Int_t n = 8;
    Double_t M[n]= {126, 200, 350, 500, 600, 700, 800, 900};
    Double_t W_15[n]= {8.12477269e-02, 1.36517184e+00, 1.19121726e+01, 4.40044052e+01, 7.74130955e+01, 1.22668305e+02, 1.81925299e+02, 2.57306382e+02};

    const Int_t n2 = 6;
    Double_t M2[n2] = {350, 500, 600, 700, 800, 900};
    Double_t W_200[n2] = {2.16288379e+00, 2.69462533e+01, 5.58588957e+01, 9.67469376e+01, 1.51721744e+02, 2.22878605e+02};

    const Int_t n3 = 4;
    Double_t M3[n3] = {600, 700, 800, 900};
    Double_t W_500[n3] = {9.10391108e+00, 2.11743549e+01, 4.72556689e+01, 8.96368631e+01};

    const Int_t n4 = 2;
    Double_t M4[n4] = {800, 900};
    Double_t W_700[n4] = { 1.47759627e+01, 2.74841698e+01};


    const Int_t nA = 14;
    Double_t M_A[nA]= {10, 11, 15, 23, 53, 142, 248, 329, 378, 435, 500, 662, 875, 1006 };
    Double_t W_A[nA]= { 3.59374088e-04, 5.89460809e-04, 1.07254526e-03, 1.68222604e-03, 3.39350655e-03, 7.98661276e-03, 1.50107862e-02, 4.75666163e-02, 5.56925012e+00, 8.12049501e+00, 1.02416303e+01, 1.44832862e+01, 1.92577931e+01, 2.19907445e+01};


    TCanvas *c1 = new TCanvas("c1","c1",1000,900);

    gPad->SetLogy();
    gPad->SetGridx();
    gPad->SetGridy();
    gStyle->SetOptTitle(0);

    gPad->SetBottomMargin(0.15);
    gPad->SetLeftMargin(0.15);

    TGraph *gr = new TGraph(n,M,W_15);
    gr->SetLineColor(kRed-8);
    gr->SetLineWidth(4);
    gr->SetMarkerColor(kRed-8);
    gr->SetMarkerStyle(21);
    gr->SetTitle("Heaviest H width, M_{A} = 15 GeV");
    gr->GetXaxis()->SetTitle("M_{H}");
    gr->GetYaxis()->SetTitle("#Gamma_{Tot}^{H}");
    gr->GetXaxis()->SetTitleSize(0.05);
    gr->GetYaxis()->SetTitleSize(0.05);
    gr->GetXaxis()->SetLabelSize(0.05);
    gr->GetYaxis()->SetLabelSize(0.05);
    gr->GetXaxis()->SetTitleOffset(1.2);
    gr->GetYaxis()->SetTitleOffset(1.2);
    gr->GetXaxis()->SetNdivisions(5);


    gr->Draw("ACP");


    TGraph *gr2 = new TGraph(n2,M2,W_200);
    gr2->SetLineColor(kRed-6);
    gr2->SetLineWidth(4);
    gr2->SetMarkerColor(kRed-6);
    gr2->SetMarkerStyle(21);
    //gr->SetTitle("Heaviest H width, M_{A} = 15 GeV");
    //gr->GetXaxis()->SetTitle("M_{H}");
    //gr->GetYaxis()->SetTitle("#Gamma_{Tot}");
    gr2->Draw("CP");


    TGraph *gr3 = new TGraph(n3,M3,W_500);
    gr3->SetLineColor(kRed-3);
    gr3->SetLineWidth(4);
    gr3->SetMarkerColor(kRed-3);
    gr3->SetMarkerStyle(21);
    gr3->Draw("CP");


    TGraph *gr4 = new TGraph(n4,M4,W_700);
    gr4->SetLineColor(kRed+2);
    gr4->SetLineWidth(4);
    gr4->SetMarkerColor(kRed+2);
    gr4->SetMarkerStyle(21);
    gr4->Draw("CP");

    TPaveText *pave1 = new TPaveText(0.2,0.7,0.4,0.9,"brNDC");
    pave1->SetBorderSize(0);
    pave1->SetFillStyle(0);
    pave1->SetTextAlign(12);
    pave1->SetTextFont(62);
    pave1->SetTextSize(0.04);
    pave1->AddText("#splitline{cos(#beta-#alpha) = 0.01}{tan #beta = 1.5}");
    pave1->Draw("same");

    TLegend *leg = new TLegend(0.6,0.2,0.89,0.4);
    leg->SetLineColor(0);
    leg->SetFillStyle(0);
    leg->AddEntry(gr,"M_{A} = 15 GeV","l");
    leg->AddEntry(gr2,"M_{A} = 200 GeV","l");
//.........这里部分代码省略.........

void buildPlotAxial_DD(){

	TH2F *dummyHist = new TH2F("dum","",1,1,900,1,6E-43,1E-37);
	dummyHist->GetYaxis()->SetTitle("#sigma_{SD} (cm^{2})");
	dummyHist->GetXaxis()->SetTitle("m_{DM} (GeV)");
	dummyHist->GetYaxis()->SetTitleOffset(1.5);
	dummyHist->GetYaxis()->SetTitleSize(0.046);
	dummyHist->GetXaxis()->SetTitleSize(0.046);
	dummyHist->GetYaxis()->SetTitleOffset(1.5);
	dummyHist->GetXaxis()->SetTitleOffset(1.25);
	dummyHist->GetYaxis()->SetTitleSize(0.045);
	dummyHist->GetXaxis()->SetTitleSize(0.045);
	dummyHist->GetYaxis()->SetLabelSize(0.04);
	dummyHist->GetXaxis()->SetLabelSize(0.04);


	gROOT->ProcessLine(".x paperStyle.C"); 
	gStyle->SetOptStat(0);
	gStyle->SetPalette(51);
	
	TFile *fiDD = TFile::Open("MassLimit_1_801_0_Both_updatePICO.root");
	TGraph *DDF = (TGraph*)fiDD->Get("DD");

	TFile *fi = TFile::Open("MassLimit_1_801_0_Both.root");

	TCanvas *can = new TCanvas();
	dummyHist->Draw("AXIS");
	
	TGraph *combinedD 	= makeOBV(combined_obs); combinedD->SetLineColor(1); combinedD->SetName("EXP");
	TGraph *combined_obsD 	= makeOBV(combined); combined_obsD->SetLineColor(1);combined_obsD->SetName("OBS");
	TGraph *monojetD 	= makeOBV(monojet); monojetD->SetLineColor(kBlue); monojetD->SetName("mjet");
	// Phil swapped boosted and resolved
	TGraph *boostedD 	= makeOBV(resolved); boostedD->SetLineColor(kMagenta); boostedD->SetName("bjet");
	TGraph *resolvedD 	= makeOBV(boosted); resolvedD->SetLineColor(kGreen); resolvedD->SetName("rjet");
	TGraph *monovD 		= makeOBV(monov); monovD->SetLineColor(kRed); monovD->SetName("MV");

	/* other lines */
	combinedD->Draw("lsame"); 
	combined_obsD->Draw("lsame"); 
	monojetD->Draw("lsame"); 
	boostedD->Draw("lsame"); 
	resolvedD->Draw("lsame"); 
	monovD->Draw("lsame"); 
	DDF->Draw("lsame");
	

	TLegend *leg = new TLegend(0.19,0.40,0.51,0.63,NULL,"brNDC");
	leg->SetFillStyle(0);
	leg->AddEntry(combinedD,"Median Expected","L");
	leg->AddEntry(combined_obsD,"Observed","L");
	leg->AddEntry(monojetD,"Monojet","L");
	leg->AddEntry(boostedD,"Boosted","L");
	leg->AddEntry(resolvedD,"Resolved","L");
	leg->AddEntry(monovD,"V-tagged","L");
	leg->AddEntry(DDF,"PICO-2L","L");
	leg->Draw();
	
   	TLatex *   texCMS = new TLatex(0.20,0.84,"#bf{CMS}");
	//TLegend *WHITEBOX = new TLegend(0.18,0.83,0.3,0.89);
   	//TLatex *   texCMS = new TLatex(0.20,0.84,"#bf{CMS}#it{Preliminary}");
	TLegend *WHITEBOX = new TLegend(0.18,0.83,0.42,0.89);
	WHITEBOX->SetFillColor(kWhite);
	//WHITEBOX->Draw();
	texCMS->SetNDC();
   	texCMS->SetTextFont(42);
   	texCMS->SetLineWidth(2);
   	texCMS->SetTextSize(0.042); texCMS->Draw();
	//tex->SetFillColor(kWhite);
	TLatex * tex = new TLatex();
	tex->SetNDC();
   	tex->SetTextFont(42);
   	tex->SetLineWidth(2);
   	tex->SetTextSize(0.035);
	tex->Draw();
   	tex->DrawLatex(0.69,0.94,"19.7 fb^{-1} (8 TeV)");
   	tex->DrawLatex(0.20,0.8,"g_{DM}=g_{SM}=1");
   	tex->DrawLatex(0.20,0.7,"Axial");
	
	can->SetRightMargin(0.11);
	can->SetLogx();
	can->SetLogy();
	can->RedrawAxis();	

	can->SaveAs("MassLimit_1_801_0_Both_DD.pdf");
	
        TFile *rout = new TFile("axialDD_out.root","RECREATE");
        rout->WriteTObject(combined_obsD,"combined");
        rout->WriteTObject(combinedD,"combined_expected"); 
        rout->WriteTObject(DDF,"DD_mass");
	rout->Close();
	
/*
 KEY: TCanvas	A;1	A
  KEY: TGraph	DD_mass;1	DD_mass
  KEY: TGraph	DD;1	DD
  KEY: TGraph	monojet;1	monojet
  KEY: TGraph	resolved;1	resolved
  KEY: TGraph	boosted;1	boosted
  KEY: TGraph	monov;1	monov
  KEY: TGraph	combined;1	combined
//.........这里部分代码省略.........

//.........这里部分代码省略.........
		yData_GLOB = yEvents.at(i);
		yErrorData_GLOB = yErrorEvents.at(i);

		tbstart = clock();
		vector< vector<double> > bSplineValues = bSpline(nPoints, npar, xEvents.at(i), yEvents.at(i), yErrorEvents.at(i), stepSpline, xminBSplineWorkspace, xmaxBSplineWorkspace);
		tbstop = clock();
		timeb.push_back(((float)tbstop-(float)tbstart)/ (CLOCKS_PER_SEC/1000.) );		

		std::cout<<timeb.back()<<std::endl;

		xBSplineValues.push_back(bSplineValues.at(0));
		yBSplineValues.push_back(bSplineValues.at(1));

		tcstart = clock();
		vector< vector<double> > cSplineValues = cSpline(nPoints, npar, xEvents.at(i), stepSpline, myMinuit);
		tcstop = clock();
		timec.push_back(((float)tcstop-(float)tcstart)/ (CLOCKS_PER_SEC/1000.) );		

		xCSplineValues.push_back(cSplineValues.at(0));
		yCSplineValues.push_back(cSplineValues.at(1));
	}

	//Histograms______________________________________________________________________________________

	//Time
	int nbins = 100;
	double xlow = 0;
	double xup = 1.;

	TH1D *hTimeB = new TH1D("Time","Timing; time [ms]; Number of Events", nbins, xlow, xup); 
	hTimeB->SetStats(0);
	hTimeB->SetMarkerStyle(10);
	TH1D *hTimeC = new TH1D("TimeC","Timing; time [ms]; Number of Events", nbins, xlow, xup); 
	hTimeC->SetLineColor(kRed);
	hTimeC->SetMarkerStyle(10);
	hTimeC->SetStats(0);

	for(int i=0; i<(int)timec.size(); i++) 
	{
		hTimeB->Fill(timeb.at(i));
		hTimeC->Fill(timec.at(i));	
	}

	//Interpolation
	vector <double> interpB, interpC;
	for(int i = 0; i < (int)yEvents.size(); i++)
	{
		for(int j = 0; j < (int)yEvents[i].size(); j++)
		{
			int indexForB = binarySearch(xBSplineValues[i], xEvents[i][j]);
			int indexForC = binarySearch(xCSplineValues[i], xEvents[i][j]);

			interpB.push_back( (yEvents[i][j]-yBSplineValues[i][indexForB])/yErrorEvents[i][j] );
			interpC.push_back( (yEvents[i][j]-yCSplineValues[i][indexForC])/yErrorEvents[i][j] );
		}
	}	

	//Test graphs for splines
	TGraph *GCspline = new TGraph(xCSplineValues[iEventLook].size(), &xCSplineValues[iEventLook][0], &yCSplineValues[iEventLook][0]);
	GCspline->SetLineColor(kRed);
	TGraph *GBspline = new TGraph(xBSplineValues[iEventLook].size(), &xBSplineValues[iEventLook][0], &yBSplineValues[iEventLook][0]);
	TGraph *Gdata = new TGraph(xEvents[0].size(), &xEvents[iEventLook][0], &yEvents[iEventLook][0]);
		Gdata->SetMarkerStyle(20);

	int nbinsI = 101;
	double xlowI = -0.1;

int evaluate( std::string filelist, std::string outfile )
{
  gStyle->SetOptStat(0);

  TCanvas *ctemp = new TCanvas();

  TCanvas *cres = new TCanvas("TimeDependence");
  TH1F* hres = new TH1F("hres","",100,0,650);
  hres->GetYaxis()->SetRangeUser(0,50);
  hres->SetTitle("");
  hres->GetXaxis()->SetTitle("time (s)");
  hres->GetYaxis()->SetTitle("B_{int} (mT)");
  hres->Draw();
  leg = new TLegend(0.2,0.6,0.9,0.9);
//  leg->SetHeader("The Legend Title"); // option "C" allows to center the header
  leg->SetNColumns(5);

  vector< double > v_Bint;
  vector< double > v_BintErr;
  vector< double > v_Bext;
  vector< double > v_BextErr;

  /* Loop over all lines in input file */
  std::ifstream infilelist(filelist);
  std::string line;

  unsigned colorcounter=38;

  while (std::getline(infilelist, line))
  {
    // skip lines with '#' and empty lines
    if ( line.find("#") != string::npos )
      {
        cout << "Skip line " << line << endl;
        continue;
      }

    if ( line == "" )
      continue;

    //cout << "Processing file " << line << endl;


    TString infilename("data_calib/");
    infilename.Append(line);

    TFile *fin = new TFile( infilename );
    TTree *tin = (TTree*)fin->Get("t");

    ctemp->cd();
    tin->Draw("Bi:time");
    TGraph *gtime = new TGraph(tin->GetEntries(), &(tin->GetV2()[0]), &(tin->GetV1()[0]));
    gtime->SetLineColor(colorcounter);
    colorcounter++;

    TH1F* hBext = new TH1F("hBext","",100,0,1000);
    tin->Draw("Bo >> hBext");

    cres->cd();
    gtime->Draw("lsame");

    double Bext_i = hBext->GetMean();
    double BextErr_i = hBext->GetRMS();

    double Bint_i = gtime->Eval(590);
    double BintErr_i = 0;

    /* add legend entry */
    TString legname("B_ext ~ ");
    legname += (int)Bext_i;
    leg->AddEntry(gtime,legname,"l");

    cout << "B_ext: " << Bext_i << " \t B_int: " << Bint_i << endl;

    v_Bint.push_back(Bint_i);
    v_BintErr.push_back(BintErr_i);
    v_Bext.push_back(Bext_i);
    v_BextErr.push_back(BextErr_i);

  }

  cres->cd();
  leg->Draw();

  TGraphErrors *gfinal = new TGraphErrors(v_Bext.size(), &(v_Bext[0]), &(v_Bint[0]), &(v_BextErr[0]), &(v_BintErr[0]));
  gfinal->Sort();
  gfinal->SetName("Bint_Vs_Bext");
  gfinal->SetTitle("");
  gfinal->GetXaxis()->SetTitle("B_{ext} (mT)");
  gfinal->GetYaxis()->SetTitle("B_{int} (mT)");

  TCanvas *cfinal = new TCanvas();
  gfinal->Draw("APL");

  /* Save output graph */
  TString outfilename("output/");
  outfilename.Append(outfile);
  TFile *fout = new TFile(outfilename,"RECREATE");

  cres->Write();
//.........这里部分代码省略.........

TGraph* getGraph(TGraph* gin , int ndiv = 10 , bool verbose = false ){

  const unsigned int n = gin->GetN();

  if( verbose ) cout << "npoints " << n << endl;

  float x[n];
  float y[n];

  Double_t thisx;
  Double_t thisy;

  for(int i = 0 ; i < n ; ++i){
    gin->GetPoint(i,thisx,thisy);
    x[i] = thisx;
    y[i] = thisy;
  }

  if( verbose ){
    cout << endl << "Original points" << endl;
    for( int i = 0 ; i < n ; i++ ){
      cout << i << " " << x[i] << " " << y[i] << endl;
    }
  }

  const unsigned int nnew = ndiv * (n-1);

  float xnew[nnew];
  float ynew[nnew];

  if( verbose ) cout << endl << "New points" << endl;
  for( int i = 0 ; i < n-1 ; ++i){

    float xgraph = x[i];
    float ygraph = y[i];

    float slope = (y[i+1]-y[i])/(x[i+1]-x[i]);

    if( verbose ) cout << i << " " << xgraph << " " << ygraph << " " << slope << endl;

    for( int j = 0 ; j < ndiv ; ++j ){

      //cout << "counter " << i*10+j << endl;
      xnew[i*ndiv+j] = xgraph;
      ynew[i*ndiv+j] = ygraph;

      xgraph += 20.0 / (float)ndiv;
      ygraph += (20.0/(float)ndiv)*slope;
    }
  }

  if( verbose ){
    cout << endl << "New graph points" << endl;
    for( int i = 0 ; i < nnew ; i++ ){
      cout << i << " " << xnew[i] << " " << ynew[i] << endl;
    }
  }

  TGraph *gout = new TGraph(nnew,xnew,ynew);

  gout->SetLineColor(gin->GetLineColor());
  gout->SetLineWidth(gin->GetLineWidth());
  gout->SetLineStyle(gin->GetLineStyle());
  //gout->SetMarkerColor(gin->GetLineColor());

  return gout;


}