13 #include <JANA/JApplication.h>
14 #include <JANA/JFactory.h>
46 gDirectory->mkdir(
"FCAL_Pi0HFA");
47 gDirectory->cd(
"FCAL_Pi0HFA");
48 hCurrentGainConstants =
new TProfile(
"CurrentGainConstants",
"Current Gain Constants", 2800, -0.5, 2799.5);
62 vector< double > raw_gains;
64 eventLoop->GetCalib(
"/FCAL/gains", raw_gains);
65 for (
unsigned int i=0; i<raw_gains.size(); i++){
92 vector<const DNeutralParticle *> neutralParticleVector;
93 loop->Get(neutralParticleVector);
94 vector< const DNeutralShower* > neutralShowers;
95 loop->Get( neutralShowers );
96 vector< const DFCALShower* > locFCALShowers;
97 loop->Get(locFCALShowers);
98 vector< const DFCALCluster* > locFCALClusters;
99 loop->Get(locFCALClusters);
102 if (neutralParticleVector.size() > 6 || neutralParticleVector.size() < 2)
return NOERROR;
104 for (
unsigned int i=0; i< neutralParticleVector.size() - 1; i++){
143 for (
unsigned int j=i+1; j< neutralParticleVector.size(); j++){
164 vector<const DFCALCluster*> associated_clusters1;
165 fcalShower2->Get(associated_clusters1);
166 vector<const DFCALCluster*> associated_clusters2;
167 fcalShower2->Get(associated_clusters2);
171 double radiusShower2=
sqrt(pow(xShower2,2)+pow(yShower2,2));
175 double radiusShower1=
sqrt(pow(xShower1,2)+pow(yShower1,2));
179 for(
unsigned int loc_j = 0; loc_j < associated_clusters1.size(); loc_j++)
181 for(
unsigned int loc_jj = 0; loc_jj < associated_clusters2.size(); loc_jj++)
184 vector< DFCALCluster::DFCALClusterHit_t > hits1 = associated_clusters1[loc_j]->GetHits();
185 vector< DFCALCluster::DFCALClusterHit_t > hits2 = associated_clusters2[loc_jj]->GetHits();
187 Int_t numhits_per_cluster1 = associated_clusters1[loc_j]->GetNHits();
188 Int_t numhits_per_cluster2 = associated_clusters2[loc_jj]->GetNHits();
189 double hitEnergyMax1=associated_clusters1[loc_j]->getEmax();
190 double hitEnergyMax2=associated_clusters1[loc_j]->getEmax();
192 double frac1 = associated_clusters1[loc_j]->getEmax()/associated_clusters1[loc_j]->getEnergy();
193 if(associated_clusters1[loc_j]->getEnergy() < 0.8)
continue;
195 double frac2 = associated_clusters2[loc_jj]->getEmax()/associated_clusters2[loc_jj]->getEnergy();
196 if(associated_clusters2[loc_jj]->getEnergy() < 0.8)
continue;
198 if (numhits_per_cluster1<1)
continue;
200 double energry1[numhits_per_cluster1];
201 double fEmax1=0, eMax1=0;
202 int chMax1=0, fChannelEmax1=0;
203 for(
int i = 0; i < numhits_per_cluster1; ++i ){
204 energry1[i]=hits1[i].E;
207 if (fabs(eMax1-fEmax1) > 0.001) {
209 fChannelEmax1 = chMax1;
213 "#pi^{0} Mass; #pi^{0} Mass;",
217 double energry2[numhits_per_cluster2];
218 double fEmax2=0, eMax2=0;
219 int chMax2=0, fChannelEmax2=0;
220 for(
int i = 0; i < numhits_per_cluster2; ++i ){
221 energry2[i]=hits2[i].E;
224 if (fabs(eMax2-fEmax2) > 0.001) {
226 fChannelEmax2 = chMax2;
230 "#pi^{0} Mass; #pi^{0} Mass;",
236 double avgE = 0.5*associated_clusters1[loc_j]->getEnergy() + 0.5*associated_clusters2[loc_jj]->getEnergy();
238 if(radiusShower1>20.785 || radiusShower2>20.785){
241 "#pi^{0} Mass; #pi^{0} Mass;",
248 "#pi^{0} Mass Vs. Channel Number; CCDB Index; #pi^{0} Mass",
249 2800, -0.5, 2799.5, 200, 0.00, 0.5);
251 if(radiusShower2<108.4239 && radiusShower2>20.785){
254 "#pi^{0} Mass Vs. Channel Number; CCDB Index; #pi^{0} Mass",
255 2800, -0.5, 2799.5, 200, 0.00, 0.5);
256 if(radiusShower1<108.4239 && radiusShower1>20.785){
259 "#pi^{0} Mass Vs. Channel Number; CCDB Index; #pi^{0} Mass",
260 2800, -0.5, 2799.5, 200, 0.00, 0.5);
263 if(radiusShower1<20.785 || radiusShower2<20.785){
266 "#pi^{0} Mass Vs. Channel Number; CCDB Index; #pi^{0} Mass",
267 2800, -0.5, 2799.5, 200, 0.00, 0.5);
274 "#pi^{0} Mass Vs. Channel Number; CCDB Index; #pi^{0} Mass",
275 2800, -0.5, 2799.5, 200, 0.00, 0.5);
276 if(radiusShower1<108.4239 && radiusShower1>20.785){
279 "#pi^{0} Mass Vs. Channel Number; CCDB Index; #pi^{0} Mass",
280 2800, -0.5, 2799.5, 200, 0.00, 0.5);
281 if(radiusShower2<108.4239 && radiusShower2>20.785){
284 "#pi^{0} Mass Vs. Channel Number; CCDB Index; #pi^{0} Mass",
285 2800, -0.5, 2799.5, 200, 0.00, 0.5);
288 if(radiusShower1<20.785 || radiusShower2<20.785){
291 "#pi^{0} Mass Vs. Channel Number; CCDB Index; #pi^{0} Mass",
292 2800, -0.5, 2799.5, 200, 0.00, 0.5);
297 for(
auto hit : associated_clusters1[loc_j]->GetHits()){
299 hit.ch, pi0Mass, hit.E / associated_clusters1[loc_j]->getEnergy(),
300 "#pi^{0} Mass Vs. Channel Number Weighted; CCDB Index; #pi^{0} Mass",
301 2800, -0.5, 2799.5, 200, 0.00, 0.5);
303 hit.ch, pi0Mass, (hit.E / associated_clusters1[loc_j]->getEnergy())*(hit.E / associated_clusters1[loc_j]->getEnergy()),
304 "#pi^{0} Mass Vs. Channel Number Weighted; CCDB Index; #pi^{0} Mass",
305 2800, -0.5, 2799.5, 200, 0.00, 0.5);
306 if(radiusShower2<108.4239 && radiusShower2>20.785){
308 hit.ch, pi0Mass, hit.E / associated_clusters1[loc_j]->getEnergy(),
309 "#pi^{0} Mass Vs. Channel Number Weighted; CCDB Index; #pi^{0} Mass",
310 2800, -0.5, 2799.5, 200, 0.00, 0.5);
312 hit.ch, pi0Mass, (hit.E / associated_clusters1[loc_j]->getEnergy())*(hit.E / associated_clusters1[loc_j]->getEnergy()),
313 "#pi^{0} Mass Vs. Channel Number Weighted; CCDB Index; #pi^{0} Mass",
314 2800, -0.5, 2799.5, 200, 0.00, 0.5);
315 if(radiusShower1<108.4239 && radiusShower1>20.785){
317 hit.ch, pi0Mass, hit.E / associated_clusters1[loc_j]->getEnergy(),
318 "#pi^{0} Mass Vs. Channel Number Weighted; CCDB Index; #pi^{0} Mass",
319 2800, -0.5, 2799.5, 200, 0.00, 0.5);
321 hit.ch, pi0Mass, (hit.E / associated_clusters1[loc_j]->getEnergy())*(hit.E / associated_clusters1[loc_j]->getEnergy()),
322 "#pi^{0} Mass Vs. Channel Number Weighted; CCDB Index; #pi^{0} Mass",
323 2800, -0.5, 2799.5, 200, 0.00, 0.5);
329 for(
auto hit : associated_clusters2[loc_jj]->GetHits()){
331 hit.ch, pi0Mass, hit.E / associated_clusters2[loc_jj]->getEnergy(),
332 "#pi^{0} Mass Vs. Channel Number Weighted; CCDB Index; #pi^{0} Mass",
333 2800, -0.5, 2799.5, 200, 0.00, 0.5);
335 hit.ch, pi0Mass, (hit.E / associated_clusters2[loc_jj]->getEnergy())*(hit.E / associated_clusters2[loc_jj]->getEnergy()),
336 "#pi^{0} Mass Vs. Channel Number Weighted; CCDB Index; #pi^{0} Mass",
337 2800, -0.5, 2799.5, 200, 0.00, 0.5);
338 if(radiusShower1<108.4239 && radiusShower1>20.785){
340 hit.ch, pi0Mass, hit.E / associated_clusters2[loc_jj]->getEnergy(),
341 "#pi^{0} Mass Vs. Channel Number Weighted; CCDB Index; #pi^{0} Mass",
342 2800, -0.5, 2799.5, 200, 0.00, 0.5);
344 hit.ch, pi0Mass, (hit.E / associated_clusters2[loc_jj]->getEnergy())*(hit.E / associated_clusters2[loc_jj]->getEnergy()),
345 "#pi^{0} Mass Vs. Channel Number Weighted; CCDB Index; #pi^{0} Mass",
346 2800, -0.5, 2799.5, 200, 0.00, 0.5);
347 if(radiusShower2<108.4239 && radiusShower2>20.785){
349 hit.ch, pi0Mass, hit.E / associated_clusters2[loc_jj]->getEnergy(),
350 "#pi^{0} Mass Vs. Channel Number Weighted; CCDB Index; #pi^{0} Mass",
351 2800, -0.5, 2799.5, 200, 0.00, 0.5);
353 hit.ch, pi0Mass, (hit.E / associated_clusters2[loc_jj]->getEnergy())*(hit.E / associated_clusters2[loc_jj]->getEnergy()),
354 "#pi^{0} Mass Vs. Channel Number Weighted; CCDB Index; #pi^{0} Mass",
355 2800, -0.5, 2799.5, 200, 0.00, 0.5);
361 if (fabs(associated_clusters1[loc_j]->getEnergy() - associated_clusters2[loc_jj]->getEnergy()) < 0.25){
364 "#pi^{0} Mass Vs. Average Shower Energy; Avg. Cluster Energy; #pi^{0} Mass",
365 100, 0.0, 10.0, 100, 0.05, 0.25);
366 if(radiusShower2<108.4239 && radiusShower2>20.785){
369 "#pi^{0} Mass Vs. Average Shower Energy; Avg. Cluster Energy; #pi^{0} Mass",
370 100, 0.0, 10.0, 100, 0.05, 0.25);
371 if(radiusShower1<108.4239 && radiusShower1>20.785){
374 "#pi^{0} Mass Vs. Average Shower Energy; Avg. Cluster Energy; #pi^{0} Mass",
375 100, 0.0, 10.0, 100, 0.05, 0.25);
380 if (fabs(associated_clusters1[loc_j]->getEnergy() - associated_clusters2[loc_jj]->getEnergy()) < 0.5){
383 "#pi^{0} Mass Vs. Average Shower Energy; Cluster Energy; #pi^{0} Mass",
384 100, 0.0, 10.0, 100, 0.05, 0.25);
386 if(radiusShower2<108.4239 && radiusShower2>20.785){
389 "#pi^{0} Mass Vs. Average Shower Energy; Avg. Cluster Energy; #pi^{0} Mass",
390 100, 0.0, 10.0, 100, 0.05, 0.25);
391 if(radiusShower1<108.4239 && radiusShower1>20.785){
394 "#pi^{0} Mass Vs. Average Shower Energy; Avg. Cluster Energy; #pi^{0} Mass",
395 100, 0.0, 10.0, 100, 0.05, 0.25);
400 if (fabs(associated_clusters1[loc_j]->getEnergy() - associated_clusters2[loc_jj]->getEnergy()) < 0.1){
403 "#pi^{0} Mass Vs. Average Shower Energy; Cluster Energy; #pi^{0} Mass",
404 100, 0.0, 10.0, 100, 0.05, 0.25);
406 if(radiusShower2<108.4239 && radiusShower2>20.785){
409 "#pi^{0} Mass Vs. Average Shower Energy; Avg. Cluster Energy; #pi^{0} Mass",
410 100, 0.0, 10.0, 100, 0.05, 0.25);
411 if(radiusShower1<108.4239 && radiusShower1>20.785){
414 "#pi^{0} Mass Vs. Average Shower Energy; Avg. Cluster Energy; #pi^{0} Mass",
415 100, 0.0, 10.0, 100, 0.05, 0.25);
419 if (fabs(associated_clusters1[loc_j]->getEnergy() - associated_clusters2[loc_jj]->getEnergy()) < 0.05){
422 "#pi^{0} Mass Vs. Average Shower Energy; Cluster Energy; #pi^{0} Mass",
423 100, 0.0, 10.0, 100, 0.05, 0.25);
424 if(radiusShower2<108.4239 && radiusShower2>20.785){
427 "#pi^{0} Mass Vs. Average Shower Energy; Avg. Cluster Energy; #pi^{0} Mass",
428 100, 0.0, 10.0, 100, 0.05, 0.25);
429 if(radiusShower1<108.4239 && radiusShower1>20.785){
432 "#pi^{0} Mass Vs. Average Shower Energy; Avg. Cluster Energy; #pi^{0} Mass",
433 100, 0.0, 10.0, 100, 0.05, 0.25);
jerror_t erun(void)
Called everytime run number changes, provided brun has been called.
DetectorSystem_t dDetectorSystem
jerror_t init(void)
Called once at program start.
jerror_t evnt(jana::JEventLoop *eventLoop, uint64_t eventnumber)
Called every event.
jerror_t brun(jana::JEventLoop *eventLoop, int32_t runnumber)
Called everytime a new run number is detected.
~JEventProcessor_FCAL_Pi0HFA()
DLorentzVector lorentzMomentum(void) const
JEventProcessor_FCAL_Pi0HFA()
jerror_t fini(void)
Called after last event of last event source has been processed.
DVector3 getPosition() const
TProfile * hCurrentGainConstants