Changeset 2d494a6 in git for modules/PileUpMerger.cc

Timestamp:

Jan 26, 2015, 3:13:02 PM (10 years ago)

Author:

Pavel Demin <pavel.demin@…>

Branches:

ImprovedOutputFile, Timing, dual_readout, llp, master

Children:

3db5282

Parents:

b62c2da

Message:

add code from CMS to PileUpMerger and PileUpMergerPythia8

File:

• modules/PileUpMerger.cc (modified) (7 diffs)

modules/PileUpMerger.cc

@@ -80 +80 @@
   fTVertexSpread = GetDouble("TVertexSpread", 1.5E-09);
 
+  fInputBeamSpotX = GetDouble("InputBeamSpotX", 0.0);
+  fInputBeamSpotY = GetDouble("InputBeamSpotY", 0.0);
+  fOutputBeamSpotX = GetDouble("OutputBeamSpotX", 0.0);
+  fOutputBeamSpotY = GetDouble("OutputBeamSpotY", 0.0);
+
   // read vertex smearing formula
 
@@ -111 +116 @@
   TParticlePDG *pdgParticle;
   Int_t pid;
-  Float_t x, y, z, t;
+  Float_t x, y, z, t, vx, vy;
   Float_t px, py, pz, e;
   Double_t dz, dphi, dt;
-  Int_t numberOfEvents, event;
+  Int_t numberOfEvents, event, numberOfParticles;
   Long64_t allEntries, entry;
-  Candidate *candidate, *vertexcandidate;
+  Candidate *candidate, *vertex;
   DelphesFactory *factory;
 
@@ -123 +128 @@
   fItInputArray->Reset();
 
-  // --- Deal with Primary vertex first  ------
+  // --- Deal with primary vertex first  ------
 
   fFunction->GetRandom2(dz, dt);
@@ -129 +134 @@
   dt *= c_light*1.0E3; // necessary in order to make t in mm/c
   dz *= 1.0E3; // necessary in order to make z in mm
-
+  vx = 0.0;
+  vy = 0.0;
+  numberOfParticles = fInputArray->GetEntriesFast();
   while((candidate = static_cast<Candidate*>(fItInputArray->Next())))
   {
+    vx += candidate->Position.X();
+    vy += candidate->Position.Y();
     z = candidate->Position.Z();
     t = candidate->Position.T();
@@ -139 +148 @@
   }
 
+  if(numberOfParticles > 0)
+  {
+    vx /= numberOfParticles;
+    vy /= numberOfParticles;
+  }
+
   factory = GetFactory();
 
-  vertexcandidate = factory->NewCandidate();
-  vertexcandidate->Position.SetXYZT(0.0, 0.0, dz, dt);
-  fVertexOutputArray->Add(vertexcandidate);
+  vertex = factory->NewCandidate();
+  vertex->Position.SetXYZT(vx, vy, dz, dt);
+  fVertexOutputArray->Add(vertex);
 
   // --- Then with pile-up vertices  ------
@@ -181 +196 @@
     dphi = gRandom->Uniform(-TMath::Pi(), TMath::Pi());
 
-    vertexcandidate = factory->NewCandidate();
-    vertexcandidate->Position.SetXYZT(0.0, 0.0, dz, dt);
-    vertexcandidate->IsPU = 1;
-
-    fVertexOutputArray->Add(vertexcandidate);
-
+    vx = 0.0;
+    vy = 0.0;
+    numberOfParticles = 0;
     while(fReader->ReadParticle(pid, x, y, z, t, px, py, pz, e))
     {
@@ -204 +216 @@
       candidate->Momentum.RotateZ(dphi);
 
+      x -= fInputBeamSpotX;
+      y -= fInputBeamSpotY;
       candidate->Position.SetXYZT(x, y, z + dz, t + dt);
       candidate->Position.RotateZ(dphi);
+      candidate->Position += TLorentzVector(fOutputBeamSpotX, fOutputBeamSpotY, 0.0, 0.0);
+
+      vx += candidate->Position.X();
+      vy += candidate->Position.Y();
+      ++numberOfParticles;
 
       fParticleOutputArray->Add(candidate);
     }
-  }
-}
-
-//------------------------------------------------------------------------------
-
+
+    if(numberOfParticles > 0)
+    {
+      vx /= numberOfParticles;
+      vy /= numberOfParticles;
+    }
+
+    vertex = factory->NewCandidate();
+    vertex->Position.SetXYZT(vx, vy, dz, dt);
+    vertex->IsPU = 1;
+
+    fVertexOutputArray->Add(vertex);
+  }
+}
+
+//------------------------------------------------------------------------------