C++ TParticle::Print方法代码示例

void testTDime(Int_t nev = 100) {

  gSystem->Load("libEVGEN");
  gSystem->Load("libTDime");
  gSystem->Load("libdime");

  TDime* dime = new TDime();
  dime->SetEnergyCMS(7000.0);
  dime->SetYRange(-2.0, 2.0);   // Set rapidity range of mesons
  dime->SetMinPt(0.1);          // Minimum pT of mesons
  dime->Initialize();

  // (pi+pi-) histograms
  TH1* hM = new TH1D("hM", "DIME #pi^{+}#pi^{-};M_{#pi^{+}#pi^{-}} #[]{GeV/#it{c}^{2}}", 100,  0.0, 5.0);

  TClonesArray* particles = new TClonesArray("TParticle", 6);
  TParticle* part = NULL;
  TLorentzVector v[2];
  TLorentzVector vSum;

  // Event loop
  for (Int_t i = 0; i < nev; ++i) {

    dime->GenerateEvent();
    Int_t np = dime->ImportParticles(particles, "All");
    printf("\n DIME Event %d: Imported %3d particles \n", i, np);

    Int_t nPrimary = 0;

    // Loop over pion (j = 4,5) tracks
    for (Int_t j = 4; j < 6; ++j) {
      part = (TParticle*) particles->At(j); // Choose the particle
      part->Print();
      part->Momentum(v[nPrimary]);          // Copy content to v
      nPrimary++;
    }
    //particles.Clear();

    // 4-vector sum
    vSum = v[0] + v[1];

    // Fill pi+pi- histograms
    hM->Fill(vSum.M());
  }

  // Save plots as pdf
  hM->Draw();    c1->SaveAs("massTDime.pdf");

}

void test(const char * sdir ="signal",
	  const char * bdir ="backgr") {

  TStopwatch timer;
  timer.Start();

  TString name;

  // Signal file, tree, and branch
  name = sdir;
  name += "/IlcESDs.root";
  TFile * fSig = TFile::Open(name.Data());
  TTree * tSig = (TTree*)fSig->Get("esdTree");

  IlcESDEvent * esdSig = new IlcESDEvent();// The signal ESD object is put here
  esdSig->ReadFromTree(tSig);

  // Run loader (signal events)
  name = sdir;
  name += "/gilc.root";
  IlcRunLoader* rlSig = IlcRunLoader::Open(name.Data());

  // Run loader (underlying events)
  name = bdir;
  name += "/gilc.root";
  IlcRunLoader* rlUnd = IlcRunLoader::Open(name.Data(),"Underlying");

  // gIlc
  rlSig->LoadgIlc();
  rlUnd->LoadgIlc();
  gIlc = rlSig->GetIlcRun();

  // Now load kinematics and event header
  rlSig->LoadKinematics();
  rlSig->LoadHeader();
  rlUnd->LoadKinematics();
  rlUnd->LoadHeader();

  // Loop on events: check that MC and data contain the same number of events
  Long64_t nevSig = rlSig->GetNumberOfEvents();
  Long64_t nevUnd = rlUnd->GetNumberOfEvents();
  Long64_t nSigPerUnd = nevSig/nevUnd;

  cout << nevSig << " signal events" << endl;
  cout << nevUnd << " underlying events" << endl;
  cout << nSigPerUnd << " signal events per one underlying" << endl;

  for (Int_t iev=0; iev<nevSig; iev++) {
    cout << "Signal event " << iev << endl;
    Int_t ievUnd = iev/nSigPerUnd;
    cout << "Underlying event " << ievUnd << endl;

    // Get signal ESD
    tSig->GetEntry(iev);
    // Get signal kinematics
    rlSig->GetEvent(iev);
    // Get underlying kinematics
    rlUnd->GetEvent(ievUnd);

    // Particle stack
    IlcStack * stackSig = rlSig->Stack();
    Int_t nPartSig = stackSig->GetNtrack();
    IlcStack * stackUnd = rlUnd->Stack();
    Int_t nPartUnd = stackUnd->GetNtrack();

    Int_t nrec = esdSig->GetNumberOfTracks();
    cout << nrec << " reconstructed tracks" << endl;
    for(Int_t irec=0; irec<nrec; irec++) {
      IlcESDtrack * track = esdSig->GetTrack(irec);
      UInt_t label = TMath::Abs(track->GetTPCLabel());
      if (label>=10000000) {
	// Underlying event. 10000000 is the
	// value of fkMASKSTEP in IlcRunDigitizer
// 	cout << " Track from the underlying event" << endl;
	label %=10000000;
	if (label>=nPartUnd) continue;
	TParticle * part = stackUnd->Particle(label);
 	if(part) part->Print();
      }
      else {
	cout << " Track " << label << " from the signal event" << endl;
	if (label>=nPartSig) {
	  cout <<"Strange, label outside the range "<< endl;
	  continue;
	}
	TParticle * part = stackSig->Particle(label);
	if(part) part->Print();
      }

    }

  }

  fSig->Close();

  timer.Stop();
  timer.Print();
}

void plots() {

  gSystem->Load("libEVGEN"); // Needs to be!

  AliRunLoader* rl = AliRunLoader::Open("galice.root");
  rl->LoadKinematics();
  rl->LoadHeader();

  // 4pi histograms
  TH1* hM    = new TH1D("hM",  "DIME #rho#rho;M_{4#pi} #(){GeV/#it{c}^{2}}", 100,  1.0, 3.0);
  TH1* hPt   = new TH1D("hPt", "DIME #rho#rho;p_{T}#(){4#pi} #(){GeV/#it{c}}", 100,  0.0, 3.0);

  // pi+- histograms
  TH1* hPt1   = new TH1D("hPt1", "DIME #rho#rho;p_{T}#(){#pi^{#pm}} #(){Gev/#it{c}}",  100, 0.0, 3.0);

  AliStack* stack = NULL;
  TParticle* part = NULL;
  TLorentzVector v[4];
  TLorentzVector vSum;

  // Loop over events
  for (Int_t i = 0; i < rl->GetNumberOfEvents(); ++i) {

    rl->GetEvent(i);
    stack = rl->Stack();
    Int_t nPrimary = 0;

    // Loop over all particles
    for (Int_t j = 0; j < stack->GetNtrack(); ++j) {
      part = stack->Particle(j);           // Get particle
      part->Print();                       // Print contents

      if (abs(part->GetPdgCode()) == 211   // Is pi+ or pi-
          & part->GetStatusCode() == 1     // Is stable final state
          & stack->IsPhysicalPrimary(j)) { // Is from generator level

        part->Momentum(v[nPrimary]);       // Set content of v
        ++nPrimary;
      }
    }
    if (nPrimary != 4) {
      printf("Error: nPrimary=%d != 4 \n", nPrimary);
      continue;
    }

    // 4-vector sum
    vSum = v[0] + v[1] + v[2] + v[3];

    // Fill 4pi histograms
    hM->Fill(vSum.M());
    hPt->Fill(vSum.Perp());

    // Fill pi+- histograms
    for (Int_t k = 0; k < 4; ++k) {
      hPt1->Fill(v[k].Perp());
    }
    printf("\n");
  }

  // Save plots as pdf
  hM->Draw();    c1->SaveAs("plotM.pdf");
  hPt->Draw();   c1->SaveAs("plotPt.pdf");
  hPt1->Draw();  c1->SaveAs("plotPt1.pdf");

}

//.........这里部分代码省略.........
    py->SetMDME(d,1,0);
  }

  //
  //                        Event Loop
  //
  
  TStopwatch timer;
  timer.Start();
  for (Int_t iev = 0; iev < nev; iev++) {
    
    cout <<"Event number "<< iev << endl;
    
    //  Initialize event
    header->Reset(0,iev);
    rl->SetEventNumber(iev);
    stack->Reset();
    rl->MakeTree("K");
    
    //  Generate event
    Int_t nprim = 0;
    Int_t ntrial = 0;
    Int_t ndstar = 0;
   
   

    //-------------------------------------------------------------------------------------

    while(!ndstar) {
      // Selection of events with D*
      stack->Reset();
      stack->ConnectTree(rl->TreeK());
      gener->Generate();
      ntrial++;
      nprim = stack->GetNprimary();
      
      for(Int_t ipart =0; ipart < nprim; ipart++){
        TParticle * part = stack->Particle(ipart);
        if(part)    {
          
          if (TMath::Abs(part->GetPdgCode())== 413) {

	    TArrayI daughtersId;

	    GetFinalDecayProducts(ipart,*stack,daughtersId);

	    Bool_t kineOK = kTRUE;

	    Double_t thetaMin = TMath::Pi()/4;
	    Double_t thetaMax = 3*TMath::Pi()/4;

	    for (Int_t id=1; id<=daughtersId[0]; id++) {
	      TParticle * daughter = stack->Particle(daughtersId[id]);
	      if (!daughter) {
		kineOK = kFALSE;
		break;
	      }

	      Double_t theta = daughter->Theta();
	      if (theta<thetaMin || theta>thetaMax) {
		kineOK = kFALSE;
		break;
	      }
	    }

	    if (!kineOK) continue;

            part->Print();
            ndstar++;     
	    
          }
        }
      }   
    }   
      
    cout << "Number of particles " << nprim << endl;
    cout << "Number of trials " << ntrial << endl;
    
    //  Finish event
    header->SetNprimary(stack->GetNprimary());
    header->SetNtrack(stack->GetNtrack());  
    
    //      I/O
    stack->FinishEvent();
    header->SetStack(stack);
    rl->TreeE()->Fill();
    rl->WriteKinematics("OVERWRITE");
    
  } // event loop
  timer.Stop();
  timer.Print();
  
  //                         Termination
  //  Generator
  gener->FinishRun();
  //  Write file
  rl->WriteHeader("OVERWRITE");
  gener->Write();
  rl->Write();
}