/********************** STEP 4 **********************/
  cout<<"Step 4: Use FillDetection(DET, N , resolution, original spectrum, oscillation function) to fill histogram DET"<<endl;
  cout<<"------ where N is the number of events generated using an input spectrum function with a given resolution"<<endl;
  cout<<"------ use the reactor (original spectrum) and oscilation functions created above, and try changing N and resolution"<<endl<<endl;

  //Let's create two DET histograms: the near, ND, and the far, FD
  //Here, 32 is the chosen number of bins (should be the same for both detectors)
  TH1D *ND=new TH1D("ND","Near Detector; E [MeV]; # events",32,1.8,8.2); ND->SetLineColor(kBlue);
  TH1D *FD=new TH1D("FD","Far Detector; E [MeV]; # events",32,1.8,8.2); FD->SetLineColor(kRed);   
  ND->Sumw2(); FD->Sumw2(); //this line is necessary for root to store statistical uncertainties
  
  cout<<"Write in your macro the value for EventsInNearDetector, the resolution and calculate the value for EventRatio"<<endl<<endl<<endl;

  //Comment the line below once filled the EventsInNearDetector, resolution and EventRatio below
  return;

  double EventsInNearDetector;//= ---write-value-here---; //replace by the value you want to test. These are the events that you want to see in ND
  double resolution;//= ---write-value-here--- ; //replace by the value you want to test
  double EventRatio;//= ---write-value-here---; //This should be approximately equal to the ratio (distance_far/distance_near)^2/(mass_far/mass_near)
  // suggestion: use 40 kt far detector
  
  double EventsInFarDetector=EventsInNearDetector/EventRatio; 
  
    
  //Please see the FillDetection function at the end of the macro
  new TCanvas;
  oscil->SetParameter(0,dNear); 
  FillDetection(ND,EventsInNearDetector,resolution,reactor,oscil);
  ND->DrawClone("hist,pe");
  
  new TCanvas;
  oscil->SetParameter(0,dFar); 
  FillDetection(FD,EventsInFarDetector,resolution,reactor,oscil);
  FD->DrawClone("hist,pe");
  
  cout<<endl<<"Number of ND events "<<ND->Integral()<<" and number of FD events "<<FD->Integral()<<endl;
  cout<<"Fraction of ND events "<<ND->Integral()/EventsInNearDetector<<" and fraction of FD events "<<FD->Integral()/EventsInFarDetector<<endl<<endl<<endl;
  
  return;
  //comment '//' the 2 "new TCanvas lines;" the 2 "DET->DrawClone();" lines, and the line above and
  //paste STEP 5 to proceed with the exercise