C++ TRandom3::Integer方法代码示例

void TestGraphRandBetter() {// named with the title so it will run automatically

   TCanvas *c1 = new TCanvas("c1","My Awesome Test Graph!!",200,10,700,500);

   c1->SetFillColor(70);
  c1->SetGrid();
  c1->GetFrame()->SetFillColor(100);
  c1->GetFrame()->SetBorderSize(1000);

   const int n = 10;

TRandom3 jrand;

vector<double> xVector;
vector<double> yVector;
vector<double> xErrorVector;
vector<double> yErrorVector;

for(int i =0; i<n; i=i+1){
xVector.push_back(jrand.Integer(20));
yVector.push_back(jrand.Integer(20));
xErrorVector.push_back(jrand.Integer(5));
yErrorVector.push_back(jrand.Integer(5));
}// I changed my parameters so that the errors would be out of a smaller number than the x and y values. 

sort(xVector.begin(), xVector.end());
sort(yVector.begin(), yVector.end());
// This puts my vectors in order so that when I graph it, the points will be plotted from least to greatest. 

   TGraphErrors *gr = new TGraphErrors(n,&xVector[0],&yVector[0],&xErrorVector[0],&yErrorVector[0]);
 
   gr->SetTitle("Test Graph With Random Data and Error Bars");
   gr->SetMarkerColor(2);
   gr->SetMarkerStyle(21);
   gr->Draw("ALP");


   c1->Update();

}

/**
 *  Factory: DD4hep_ConditionExample_stress
 *
 *  \author  M.Frank
 *  \version 1.0
 *  \date    01/12/2016
 */
static int condition_example (Detector& description, int argc, char** argv)  {
  string input;
  int    num_iov = 10, num_runs = 10;
  bool   arg_error = false;
  for(int i=0; i<argc && argv[i]; ++i)  {
    if ( 0 == ::strncmp("-input",argv[i],4) )
      input = argv[++i];
    else if ( 0 == ::strncmp("-iovs",argv[i],4) )
      num_iov = ::atol(argv[++i]);
    else if ( 0 == ::strncmp("-runs",argv[i],4) )
      num_runs = ::atol(argv[++i]);
    else
      arg_error = true;
  }
  if ( arg_error || input.empty() )   {
    /// Help printout describing the basic command line interface
    cout <<
      "Usage: -plugin <name> -arg [-arg]                                             \n"
      "     name:   factory name     DD4hep_AlignmentExample1                        \n"
      "     -input   <string>        Geometry file                                   \n"
      "     -iovs    <number>        Number of parallel IOV slots for processing.    \n"
      "     -runs    <number>        Number of collision loads to be performed.      \n"
      "\tArguments given: " << arguments(argc,argv) << endl << flush;
    ::exit(EINVAL);
  }

  // First we load the geometry
  description.fromXML(input);

  /******************** Initialize the conditions manager *****************/
  ConditionsManager manager = installManager(description);
  const IOVType*    iov_typ = manager.registerIOVType(0,"run").second;
  if ( 0 == iov_typ )  {
    except("ConditionsPrepare","++ Unknown IOV type supplied.");
  }

  /******************** Now as usual: create the slice ********************/
  shared_ptr<ConditionsContent> content(new ConditionsContent());
  shared_ptr<ConditionsSlice>   slice(new ConditionsSlice(manager,content));
  Scanner(ConditionsKeys(*content,INFO),description.world());
  Scanner(ConditionsDependencyCreator(*content,DEBUG),description.world());

  TStatistic cr_stat("Creation"), acc_stat("Access");
  /******************** Populate the conditions store *********************/
  // Have 10 run-slices [11,20] .... [91,100]
  size_t total_created = 0;
  for(int i=0; i<num_iov; ++i)  {
    TTimeStamp start;
    IOV iov(iov_typ, IOV::Key(1+i*10,(i+1)*10));
    ConditionsPool*   iov_pool = manager.registerIOV(*iov.iovType, iov.key());
    // Create conditions with all deltas.  Use a generic creator
    int count = Scanner().scan(ConditionsCreator(*slice, *iov_pool, DEBUG), description.world());
    TTimeStamp stop;
    cr_stat.Fill(stop.AsDouble()-start.AsDouble());
    printout(INFO,"Example", "Setup %ld conditions for IOV:%s [%8.3f sec]",
             count, iov.str().c_str(), stop.AsDouble()-start.AsDouble());
    total_created += count;
  }

  // ++++++++++++++++++++++++ Now compute the conditions for each of these IOVs
  TRandom3 random;
  ConditionsManager::Result total;
  for(int i=0; i<num_runs; ++i)  {
    TTimeStamp start;
    unsigned int rndm = 1+random.Integer(num_iov*10);
    IOV req_iov(iov_typ,rndm);
    // Attach the proper set of conditions to the user pool
    ConditionsManager::Result res = manager.prepare(req_iov,*slice);
    TTimeStamp stop;
    total += res;
    acc_stat.Fill(stop.AsDouble()-start.AsDouble());
    // Now compute the tranformation matrices
    printout(INFO,"Prepare","Total %ld conditions (S:%6ld,L:%6ld,C:%6ld,M:%ld) of type %s [%8.3f sec]",
             res.total(), res.selected, res.loaded, res.computed, res.missing,
             req_iov.str().c_str(), stop.AsDouble()-start.AsDouble());
  }
  printout(INFO,"Statistics","+======= Summary: # of IOV: %3d  # of Runs: %3d ===========================", num_iov, num_runs);
  printout(INFO,"Statistics","+  %-12s:  %11.5g +- %11.4g  RMS = %11.5g  N = %lld",
           cr_stat.GetName(), cr_stat.GetMean(), cr_stat.GetMeanErr(), cr_stat.GetRMS(), cr_stat.GetN());
  printout(INFO,"Statistics","+  %-12s:  %11.5g +- %11.4g  RMS = %11.5g  N = %lld",
           acc_stat.GetName(), acc_stat.GetMean(), acc_stat.GetMeanErr(), acc_stat.GetRMS(), acc_stat.GetN());
  printout(INFO,"Statistics","+  Accessed a total of %ld conditions (S:%6ld,L:%6ld,C:%6ld,M:%ld). Created:%ld",
           total.total(), total.selected, total.loaded, total.computed, total.missing, total_created);
  printout(INFO,"Statistics","+=========================================================================");
  // All done.
  return 1;
}

//.........这里部分代码省略.........
		num_fadc_samples = 0;
		num_fadc_samples = fadc->GetNumFadcSamples(chan, ievent);
	      }
	      if (num_fadc_events > 0) {
	      	for (Int_t jevent = 0; jevent < num_fadc_events; jevent++) {
		  // Debug output
		  if ((fadc_mode == 1 || fadc_mode == 8) && num_fadc_samples > 0) 
		    if (debugfile) *debugfile << "FADC EMULATED PI DATA = " << fadc->GetEmulatedPulseIntegralData(chan) << endl;
		  if (fadc_mode == 7 || fadc_mode == 8 || fadc_mode == 9 || fadc_mode == 10) {
		    if (fadc_mode != 8) {if (debugfile) *debugfile << "FADC PI DATA = " << fadc->GetPulseIntegralData(chan, jevent) << endl;}
		    if (debugfile) *debugfile << "FADC PT DATA = " << fadc->GetPulseTimeData(chan, jevent) << endl;
		    if (debugfile) *debugfile << "FADC PPED DATA = " << fadc->GetPulsePedestalData(chan, jevent) << endl;
		    if (debugfile) *debugfile << "FADC PPEAK DATA = " << fadc->GetPulsePeakData(chan, jevent) << endl;
		  }
		  // Fill histos
	      	  if ((fadc_mode == 1 || fadc_mode == 8) && num_fadc_samples > 0) h_pinteg[islot][chan]->Fill(fadc->GetEmulatedPulseIntegralData(chan));
		  else if (fadc_mode == 7 || fadc_mode == 8 || fadc_mode == 9 || fadc_mode == 10) {
		    if (fadc_mode != 8) {h_pinteg[islot][chan]->Fill(fadc->GetPulseIntegralData(chan, jevent));}
		    h_ptime[islot][chan]->Fill(fadc->GetPulseTimeData(chan, jevent));
		    h_pped[islot][chan]->Fill(fadc->GetPulsePedestalData(chan, jevent));
		    h_ppeak[islot][chan]->Fill(fadc->GetPulsePeakData(chan, jevent));
		  }
		  // Raw sample events
	      	  if (raw_mode && num_fadc_samples > 0) {
		    // Debug output
	      	    if (debugfile) *debugfile << "NUM FADC SAMPLES = " << num_fadc_samples << endl;
	      	    if (debugfile) *debugfile << "=============================" << endl;
		    // Populate raw sample plots
	      	    if (ievent < NUMRAWEVENTS) {
	      	      // Acquire the raw samples vector and populate graphs
	      	      raw_samples_vector[islot][chan] = fadc->GetPulseSamplesVector(chan);
	      	      for (Int_t sample_num = 0; sample_num < Int_t (raw_samples_vector[islot][chan].size()); sample_num++) 
	      	      	g_psamp_event[islot][chan][ievent]->SetPoint(sample_num, sample_num + 1, Int_t (raw_samples_vector[islot][chan][sample_num]));
	      	      mg_psamp[islot][chan]->Add(g_psamp_event[islot][chan][ievent], "ACP");
	      	    }  // NUMRAWEVENTS condition
	      	  }  // Raw mode condition
	      	}  //  FADC event loop
	      }  // Number of FADC events condition
	    }  //FADC channel loop
	  }  // FADC module found condition
	  else 
	    if (debugfile) *debugfile << "FADC MODULE NOT FOUND!!!" << endl;
	}  // Number raw words condition
      }  // Slot loop
    }  // CODA file read status condition
    if (iievent % 1000 == 0)
      //if (iievent % 1 == 0)
      cout << iievent << " events have been processed" << endl;
    iievent++;
    if (status == EOF) break;
  }  // Event loop 

  // Populate waveform graphs and multigraphs and write to root file
  TRandom3 *rand = new TRandom3();
  for(uint32_t islot = 3; islot < NUMSLOTS; islot++) {
    for (uint32_t chan = 0; chan < NADCCHAN; chan++) {
      for( uint32_t ievent = 0; ievent < NUMRAWEVENTS; ievent++) {
	// Raw sample plots
	if (g_psamp_event[islot][chan][ievent] != NULL) {
	  UInt_t rand_int = 1 + rand->Integer(9);
	  hfile->cd(Form("/mode_%d_data/slot_%d/chan_%d/raw_samples", fadc_mode_const, islot, chan));
	  c_psamp[islot][chan]->cd(ievent + 1);
	  g_psamp_event[islot][chan][ievent]->Draw("ACP");
	  g_psamp_event[islot][chan][ievent]->SetTitle(Form("FADC Mode %d Pulse Peak Data Slot %d Channel %d Event %d", fadc_mode_const, islot, chan, ievent));
	  g_psamp_event[islot][chan][ievent]->GetXaxis()->SetTitle("Sample Number");
	  g_psamp_event[islot][chan][ievent]->GetYaxis()->SetTitle("Sample Value");
	  g_psamp_event[islot][chan][ievent]->SetLineColor(rand_int);
	  g_psamp_event[islot][chan][ievent]->SetMarkerColor(rand_int);
	  g_psamp_event[islot][chan][ievent]->SetMarkerStyle(20);
	  g_psamp_event[islot][chan][ievent]->SetDrawOption("ACP");
	}  // Graph condition
      }  // Raw event loop
      // Write the canvas to file
      if (c_psamp[islot][chan] != NULL) c_psamp[islot][chan]->Write();
      // Write the multigraphs to file
      if (mg_psamp[islot][chan] != NULL) {
        hfile->cd(Form("/mode_%d_data/slot_%d/chan_%d", fadc_mode_const, islot, chan));
        mg_psamp[islot][chan]->Draw("ACP");
        //mg_psamp[islot][chan]->SetTitle(Form("%d Raw Events of FADC Mode %d Pulse Peak Data Slot %d Channel %d", NUMRAWEVENTS, fadc_mode_const, islot, chan));
        //mg_psamp[islot][chan]->GetXaxis()->SetTitle("Sample Number");
        //mg_psamp[islot][chan]->GetYaxis()->SetTitle("Sample Value");
        mg_psamp[islot][chan]->Write(); 
      }  // Mulitgraph condition
    }  // Channel loop
  }  // Slot loop

  cout << "***************************************" << endl;
  cout << iievent - 1 << " events were processed" << endl;
  cout << "***************************************" << endl;
  
  // Write and clode the data file
  hfile->Write();
  hfile->Close();

  // Calculate the analysis rate
  t = clock() - t;
  printf ("The analysis took %.1f seconds \n", ((float) t) / CLOCKS_PER_SEC);
  printf ("The analysis event rate is %.1f Hz \n", (iievent - 1) / (((float) t) / CLOCKS_PER_SEC));

}  // main()

void flagMultCands( TString fname, TString tname ) {

  TFile *inFile = new TFile( fname, "UPDATE" );
  TTree *tree = (TTree*)inFile->Get( tname );

  if ( tree->FindBranch( "pass_multcand" ) ) {
    cout << "pass_multcand branch already exists so I won't add it" << endl;
    delete inFile;
    return;
  }

  ULong64_t  eventNumber;
  int        itype;
  bool       pass_bdt;
  bool       pass_pid;
  bool       pass_rhokst;
  bool       pass_massveto;
  UInt_t     nCandidate;
  ULong64_t  totCandidates;
  tree->SetBranchAddress(  "eventNumber"                 , &eventNumber                 );
  tree->SetBranchAddress(  "itype"                       , &itype                       );
  tree->SetBranchAddress(  "pass_bdt"                    , &pass_bdt                    );
  tree->SetBranchAddress(  "pass_pid"                    , &pass_pid                    );
  tree->SetBranchAddress(  "pass_rhokst"                 , &pass_rhokst                 );
  tree->SetBranchAddress(  "pass_massveto"               , &pass_massveto               );
  tree->SetBranchAddress(  "totCandidates"               , &totCandidates               );
  tree->SetBranchAddress(  "nCandidate"                  , &nCandidate                  );

  map< ULong64_t, ULong64_t > multCandEventNumbers;

  cout << "Finding multiple candidates" << endl;
  // first loop tree and save number of all multiple candidates
  for ( int ev=0; ev<tree->GetEntries(); ev++ ) {
    tree->GetEntry(ev);
    if ( ev%10000==0 ) cout << ev << "/" << tree->GetEntries() << endl;
    //if ( itype>0 && pass_bdt && pass_pid && (!pass_rhokst) && (!pass_massveto) && totCandidates > 1 ) {
    if ( itype>0 && pass_bdt && pass_pid && (!pass_massveto) && totCandidates > 1 ) {
      multCandEventNumbers[eventNumber] = totCandidates;
    }
  }

  // now randomly select which one to keep
  map< ULong64_t, UInt_t > eventToKeep;
  TRandom3 rand;
  rand.SetSeed(2016);
  for ( map<ULong64_t,ULong64_t>::iterator it=multCandEventNumbers.begin(); it!=multCandEventNumbers.end(); it++) {
    UInt_t keep = rand.Integer( it->second );
    eventToKeep[it->first] = keep;
  }

  cout << "Adding flag pass_multcand to tree" << endl;
  // then loop again and write the random choice back in
  bool pass_multcand;
  TBranch *bpt = tree->Branch( "pass_multcand", &pass_multcand, "pass_multcand/O" );

  for ( int ev=0; ev<tree->GetEntries(); ev++ ) {
    tree->GetEntry(ev);
    if ( ev%10000==0 ) cout << ev << "/" << tree->GetEntries() << endl;
    pass_multcand = true;
    if ( totCandidates > 1 && eventToKeep[eventNumber] != nCandidate ) pass_multcand = false;
    bpt->Fill();
  }
  tree->Write();
  delete inFile;
}

//this code create the so-called 'PrintOuts' as done with the MT2tree::PrintOut() function
void PrintOutEvents(){

	bool savetofile = true;//save the printout to the file PrintOut.log, default = true
	bool onlyData = true;  //run only over data, default = true

//	TString  samples = "/shome/haweber/MT2Analysis_8TeV/Code/MT2AnalysisCode/RootMacros/samples/samples_MET_filter.dat";//only dummy
//	TString  samples = "/shome/haweber/MT2Analysis_8TeV/Code/MT2AnalysisCode/RootMacros/samples/samples_HT_filter.dat";//only dummy
	TString samples = "samples/samples_HTandMET_filter.dat";

	//event selection for which you want to have the printout - complicated stuff can also be done later
	std::ostringstream cutStream;
	std::ostringstream cutStreamBase;
	cutStream << " " 
	//	<< "NBJets40CSVM >= 3"		<< "&&"
		<< "NJetsIDLoose40 >= 2"		<< "&&"
		<< "NEles+NMuons+NTausIDLoose3Hits  ==0"<< "&&"
	//	<< "misc.MT2>=250"                      << "&&"
		<< "misc.Jet0Pass ==1"                      << "&&"
		<< "misc.Jet1Pass ==1"                      << "&&"
		<< "misc.Vectorsumpt < 70";
		cutStream  << "&& misc.MinMetJetDPhi4Pt40 >0.3";
//	cutStream << "&&misc.HT<750&&misc.HT>=450&&misc.MET>200";
//	cutStream << "&&misc.HT>=750&&misc.MET>30";
	cutStream << "&&((misc.HT<750&&misc.HT>=450&&misc.MET>200)||(misc.HT>=750&&misc.MET>30))";

	//Higgs
//	cutStream <<"&&NJetsIDLoose40 >= 4"<<"&&"
  //  		  <<"NBJetsCSVM>=2"<<"&&"
	//	  <<"GetSelBBMinv()>250";
	
	cutStreamBase << " " 
      << "misc.PassJet40ID ==1"                      << "&&"
      << "(misc.HBHENoiseFlag == 0 || misc.ProcessID==10)"  << "&&" // for rare SM samples
      << "misc.CSCTightHaloIDFlag == 0"             << "&&"
      << "misc.trackingFailureFlag==0"              << "&&"
      << "misc.eeBadScFlag==0"                      << "&&"
      << "misc.EcalDeadCellTriggerPrimitiveFlag==0" << "&&"
      << "misc.TrackingManyStripClusFlag==0"             << "&&"
      << "misc.TrackingTooManyStripClusFlag==0"          << "&&"
      << "misc.TrackingLogErrorTooManyClustersFlag==0"   << "&&"
      << "misc.CrazyHCAL==0";
	cutStreamBase << "&&misc.MET/misc.CaloMETRaw<=2.";
	cutStreamBase << "&&TOBTECTagger<=8&&ExtraBeamHaloFilter==0";// ONLY ON DATA
	TString cuts = cutStream.str().c_str();
	TString basecuts = cutStreamBase.str().c_str();

	load(samples.Data());

	//if have extremely many events can only select a random subset of those
	TRandom3 *random = new TRandom3(0);//somehow truly random with this seed
	vector<int> randvec; randvec.clear();
	for(int i =0; i<40;++i)
		randvec.push_back(random->Integer(7750));
	delete random;
	sort(randvec.begin(), randvec.end());
	//for(int i=0; i<randvec.size();++i) cout << "randvec["<<i<<"]="<<randvec[i] << endl;
	int counter1 = 0; int counter2 = 0;

	//don't run events twice - this vector makes that sure (in case event is in both MET and HT trigger stream)
	   vector<pair<int,pair<int,int> > > rls; rls.clear();


   	for(size_t i = 0; i < fSamples.size(); ++i){
        
   	    if(onlyData && fSamples[i].type!="data") continue;

	    MT2tree* fMT2tree = new MT2tree();
	    fSamples[i].tree->SetBranchAddress("MT2tree", &fMT2tree);
	    Long64_t nentries =  fSamples[i].tree->GetEntries();
	    Long64_t nbytes = 0, nb = 0;
	    int nev =0;
        
	    TString myCuts = cuts + "&&" + basecuts;
        
	    //if( fSamples[i].type=="data") myCuts += " && " + trigger; //cuts to be aplied only on data
        
	    cout << "sample " << fSamples[i].name << endl;
   	    cout << "Cuts for Flow: " << myCuts << endl;
   	    fSamples[i].tree->Draw(">>selList", myCuts);
        
	    TEventList *myEvtList = (TEventList*)gDirectory->Get("selList");
	    fSamples[i].tree->SetEventList(myEvtList);
	    int counter=0;
	    cout << "Filtering done, size=" <<myEvtList->GetN()  << endl;
	    if(myEvtList->GetSize()==0) continue;
		//run over selected events
        while(myEvtList->GetEntry(counter++) !=-1){	
		int jentry = myEvtList->GetEntry(counter-1);
            
		nb =  fSamples[i].tree->GetEntry(jentry);   nbytes += nb;
		fSamples[i].tree->SetBranchAddress("MT2tree", &fMT2tree);
		//implementation of running over random subset
//		if(counter1==randvec[counter2]){
//			++counter2;
//			fMT2tree->PrintOut(savetofile);
//		}

		pair <int,int> dummy1(fMT2tree->misc.LumiSection,fMT2tree->misc.Event);
		pair <int,pair<int,int> > dummy(fMT2tree->misc.Run, dummy1);
//.........这里部分代码省略.........