libraries/ANALYSIS/DHistogramActions

Bug Summary

File:	libraries/ANALYSIS/DHistogramActions_Independent.cc
Location:	line 202, column 3
Description:	Value stored to 'locBaseFactory' is never read

Annotated Source Code

1	#include "ANALYSIS/DHistogramActions.h"
2
3	void DHistogramAction_ObjectMemory::Initialize(JEventLoop* locEventLoop)
4	{
5	string locHistName, locHistTitle;
6	vector<string> locBinLabels; //fill this
7
8	//ANALYSIS
9	dFactoryPoolBinMap["DParticleComboBlueprintStep"] = 1;
10	locBinLabels.push_back("DParticleComboBlueprintStep");
11
12	dFactoryPairsToTrack.push_back(pair<string, string>("DParticleComboBlueprint", ""));
13	dFactoryPairBinMap[dFactoryPairsToTrack.back()] = 2;
14	locBinLabels.push_back("DParticleComboBlueprint");
15
16	dFactoryPairsToTrack.push_back(pair<string, string>("DTrackTimeBased", "Combo"));
17	dFactoryPairBinMap[dFactoryPairsToTrack.back()] = 3;
18	locBinLabels.push_back("DTrackTimeBased_Combo");
19
20	dFactoryPairsToTrack.push_back(pair<string, string>("DEventRFBunch", "Combo"));
21	dFactoryPairBinMap[dFactoryPairsToTrack.back()] = 4;
22	locBinLabels.push_back("DEventRFBunch_Combo");
23
24	dFactoryPairsToTrack.push_back(pair<string, string>("DChargedTrackHypothesis", "Combo"));
25	dFactoryPairBinMap[dFactoryPairsToTrack.back()] = 5;
26	locBinLabels.push_back("DChargedTrackHypothesis_Combo");
27
28	dFactoryPairsToTrack.push_back(pair<string, string>("DNeutralParticleHypothesis", "Combo"));
29	dFactoryPairBinMap[dFactoryPairsToTrack.back()] = 6;
30	locBinLabels.push_back("DNeutralParticleHypothesis_Combo");
31
32	dFactoryPoolBinMap["DParticleComboStep_PreKinFit"] = 7;
33	locBinLabels.push_back("DParticleComboStep_PreKinFit");
34
35	dFactoryPoolBinMap["DKinematicData_ComboPreKinFit"] = 8;
36	locBinLabels.push_back("DKinematicData_ComboPreKinFit");
37
38	dFactoryPairsToTrack.push_back(pair<string, string>("DParticleCombo", "PreKinFit"));
39	dFactoryPairBinMap[dFactoryPairsToTrack.back()] = 9;
40	locBinLabels.push_back("DParticleCombo_PreKinFit");
41
42	dFactoryPoolBinMap["DKinFitParticle"] = 10;
43	locBinLabels.push_back("DKinFitParticle");
44
45	dFactoryPoolBinMap["DKinFitConstraint_Vertex"] = 11;
46	locBinLabels.push_back("DKinFitConstraint_Vertex");
47
48	dFactoryPoolBinMap["DKinFitConstraint_Spacetime"] = 12;
49	locBinLabels.push_back("DKinFitConstraint_Spacetime");
50
51	dFactoryPoolBinMap["DKinFitConstraint_P4"] = 13;
52	locBinLabels.push_back("DKinFitConstraint_P4");
53
54	dFactoryPoolBinMap["TMatrixDSym_KinFitter"] = 14;
55	locBinLabels.push_back("TMatrixDSym_KinFitter");
56
57	dFactoryPairsToTrack.push_back(pair<string, string>("DKinFitResults", ""));
58	dFactoryPairBinMap[dFactoryPairsToTrack.back()] = 15;
59	locBinLabels.push_back("DKinFitResults");
60
61	dFactoryPairsToTrack.push_back(pair<string, string>("DBeamPhoton", "KinFit"));
62	dFactoryPairBinMap[dFactoryPairsToTrack.back()] = 16;
63	locBinLabels.push_back("DBeamPhoton_KinFit");
64
65	dFactoryPairsToTrack.push_back(pair<string, string>("DChargedTrackHypothesis", "KinFit"));
66	dFactoryPairBinMap[dFactoryPairsToTrack.back()] = 17;
67	locBinLabels.push_back("DChargedTrackHypothesis_KinFit");
68
69	dFactoryPairsToTrack.push_back(pair<string, string>("DNeutralParticleHypothesis", "KinFit"));
70	dFactoryPairBinMap[dFactoryPairsToTrack.back()] = 18;
71	locBinLabels.push_back("DNeutralParticleHypothesis_KinFit");
72
73	dFactoryPoolBinMap["DParticleComboStep"] = 19;
74	locBinLabels.push_back("DParticleComboStep");
75
76	dFactoryPoolBinMap["DKinematicData_Combo"] = 20;
77	locBinLabels.push_back("DKinematicData_Combo");
78
79	dFactoryPairsToTrack.push_back(pair<string, string>("DParticleCombo", ""));
80	dFactoryPairBinMap[dFactoryPairsToTrack.back()] = 21;
81	locBinLabels.push_back("DParticleCombo");
82
83	//CREATE THE HISTOGRAMS
84	japp->RootWriteLock(); //ACQUIRE ROOT LOCK!!
85	{
86	CreateAndChangeTo_ActionDirectory();
87
88	// Total Memory
89	locHistName = "TotalMemory";
90	locHistTitle = ";Event # ;Total Memory (MB)";
91	dHist_TotalMemory = GetOrCreate_Histogram<TH1I>(locHistName, locHistTitle, dMaxNumEvents, 0.5, float(dMaxNumEvents) + 0.5);
92
93	// # Objects
94	locHistName = "NumObjects2D";
95	locHistTitle = "# Objects;Event #";
96	dHist_NumObjects = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dMaxNumEvents, 0.5, float(dMaxNumEvents) + 0.5, locBinLabels.size(), 0.5, float(locBinLabels.size()) + 0.5);
97	for(size_t loc_i = 0; loc_i < locBinLabels.size(); ++loc_i)
98	dHist_NumObjects->GetYaxis()->SetBinLabel(1 + loc_i, locBinLabels[loc_i].c_str());
99
100	// Object Memory
101	locHistName = "Memory2D";
102	locHistTitle = "Memory (Bytes);Event #";
103	dHist_Memory = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dMaxNumEvents, 0.5, float(dMaxNumEvents) + 0.5, locBinLabels.size(), 0.5, float(locBinLabels.size()) + 0.5);
104	for(size_t loc_i = 0; loc_i < locBinLabels.size(); ++loc_i)
105	dHist_Memory->GetYaxis()->SetBinLabel(1 + loc_i, locBinLabels[loc_i].c_str());
106
107	for(size_t loc_i = 0; loc_i < locBinLabels.size(); ++loc_i)
108	{
109	// # Objects
110	locHistName = string("NumObjects_") + locBinLabels[loc_i];
111	locHistTitle = locBinLabels[loc_i] + string(";Event # ;# Objects");
112	dHistMap_NumObjects[loc_i + 1] = GetOrCreate_Histogram<TH1I>(locHistName, locHistTitle, dMaxNumEvents, 0.5, float(dMaxNumEvents) + 0.5);
113
114	// # Objects
115	locHistName = string("Memory_") + locBinLabels[loc_i];
116	locHistTitle = locBinLabels[loc_i] + string(";Event # ;Memory (Bytes)");
117	dHistMap_Memory[loc_i + 1] = GetOrCreate_Histogram<TH1I>(locHistName, locHistTitle, dMaxNumEvents, 0.5, float(dMaxNumEvents) + 0.5);
118	}
119
120	//Return to the base directory
121	ChangeTo_BaseDirectory();
122	}
123	japp->RootUnLock(); //RELEASE ROOT LOCK!!
124	}
125
126	bool DHistogramAction_ObjectMemory::Perform_Action(JEventLoop* locEventLoop, const DParticleCombo* locParticleCombo)
127	{
128	if(Get_NumPreviousParticleCombos() != 0)
129	return true; //else double-counting!
130
131	if(dEventCounter > dMaxNumEvents)
132	return true;
133
134	if(locParticleCombo != NULL__null)
135	return true; // Protect against infinite recursion (see below)
136
137	//THIS IS EXTREMELY DANGEROUS, AND SHOULD BE AVOIDED UNLESS YOU KNOW !!!!EXACTLY!!!! WHAT YOU ARE DOING
138	//And if you're doing this, you probably don't
139	//This will result in a infinite-recursion crash if it is called by DAnalysisResults_factory.
140	//This action should only be used directly in an event processsor
141	//This casuses the analysis to run, generating the objects needed for histogramming below.
142	vector<const DAnalysisResults*> locAnalysisResults;
143	locEventLoop->Get(locAnalysisResults);
144
145	//FACTORIES
146	//call get-n-rows first outside of lock, just to make sure
147	map<int, size_t> locNumObjectsMap; //int is bin
148	map<int, unsigned long long> locMemoryMap; //int is bin
149	double locTotalMemory = 0;
150	for(size_t loc_i = 0; loc_i < dFactoryPairsToTrack.size(); ++loc_i)
151	{
152	string locClassName = dFactoryPairsToTrack[loc_i].first;
153	JFactory_base* locFactory = locEventLoop->GetFactory(locClassName.c_str(), dFactoryPairsToTrack[loc_i].second.c_str());
154	unsigned long long locNumObjects = locFactory->GetNrows();
155	unsigned long long locDataClassSize = locFactory->GetDataClassSize();
156
157	unsigned long long locMemory = locDataClassSize*locNumObjects;
158	if(locClassName == "DChargedTrackHypothesis")
159	locMemory += locNumObjects(778 + 55*8); //error matrices //8 = double
160	if((locClassName == "DNeutralParticleHypothesis") \|\| (locClassName == "DBeamPhoton"))
161	locMemory += locNumObjects(77*8); //error matrices //8 = double
162
163	int locBin = dFactoryPairBinMap[dFactoryPairsToTrack[loc_i]];
164	locNumObjectsMap[locBin] = locNumObjects;
165	locMemoryMap[locBin] = locMemory;
166	locTotalMemory += double(locMemory);
167	}
168
169	//RESOURCE POOLS
170	{
171	unsigned long long locMemory;
172	JFactory_base* locBaseFactory;
173	int locBin;
174
175	//DParticleComboBlueprintStep
176	locBin = dFactoryPoolBinMap["DParticleComboBlueprintStep"];
177	locBaseFactory = locEventLoop->GetFactory("DParticleComboBlueprint", "");
178	DParticleComboBlueprint_factory* locParticleComboBlueprintFactory = static_cast<DParticleComboBlueprint_factory*>(locBaseFactory);
179	locNumObjectsMap[locBin] = locParticleComboBlueprintFactory->Get_ParticleComboBlueprintStepPoolSize();
180	locMemory = sizeof(DParticleComboBlueprintStep)*locNumObjectsMap[locBin];
181	locMemoryMap[locBin] = locMemory;
182	locTotalMemory += double(locMemory);
183
184	//DParticleComboStep_PreKinFit
185	locBin = dFactoryPoolBinMap["DParticleComboStep_PreKinFit"];
186	locBaseFactory = locEventLoop->GetFactory("DParticleCombo", "PreKinFit");
187	DParticleCombo_factory_PreKinFit* locParticleComboFactory_PreKinFit = static_cast<DParticleCombo_factory_PreKinFit*>(locBaseFactory);
188	locNumObjectsMap[locBin] = locParticleComboFactory_PreKinFit->Get_ParticleComboStepPoolSize();
189	locMemory = sizeof(DParticleComboStep)*locNumObjectsMap[locBin];
190	locMemoryMap[locBin] = locMemory;
191	locTotalMemory += double(locMemory);
192
193	//DKinematicData_ComboPreKinFit
194	locBin = dFactoryPoolBinMap["DKinematicData_ComboPreKinFit"];
195	locNumObjectsMap[locBin] = locParticleComboFactory_PreKinFit->Get_KinematicDataPoolSize();
196	locMemory = sizeof(DKinematicData)*locNumObjectsMap[locBin];
197	locMemoryMap[locBin] = locMemory;
198	locTotalMemory += double(locMemory);
199
200	//DKinFitParticle
201	locBin = dFactoryPoolBinMap["DKinFitParticle"];
202	locBaseFactory = locEventLoop->GetFactory("DKinFitResults", "");
	Value stored to 'locBaseFactory' is never read
203	// DKinFitResults_factory* locKinFitResultsFactory = static_cast<DKinFitResults_factory*>(locBaseFactory);
204	// locNumObjectsMap[locBin] = locKinFitResultsFactory->Get_KinFitParticlePoolSize();
205	locNumObjectsMap[locBin] = 1;
206	locMemory = sizeof(DKinFitParticle)*locNumObjectsMap[locBin];
207	locMemoryMap[locBin] = locMemory;
208	locTotalMemory += double(locMemory);
209
210	//DKinFitConstraint_Vertex
211	locBin = dFactoryPoolBinMap["DKinFitConstraint_Vertex"];
212	// locNumObjectsMap[locBin] = locKinFitResultsFactory->Get_KinFitConstraintVertexPoolSize();
213	locNumObjectsMap[locBin] = 1;
214	locMemory = sizeof(DKinFitConstraint_Vertex)*locNumObjectsMap[locBin];
215	locMemoryMap[locBin] = locMemory;
216	locTotalMemory += double(locMemory);
217
218	//DKinFitConstraint_Spacetime
219	locBin = dFactoryPoolBinMap["DKinFitConstraint_Spacetime"];
220	// locNumObjectsMap[locBin] = locKinFitResultsFactory->Get_KinFitConstraintSpacetimePoolSize();
221	locNumObjectsMap[locBin] = 1;
222	locMemory = sizeof(DKinFitConstraint_Spacetime)*locNumObjectsMap[locBin];
223	locMemoryMap[locBin] = locMemory;
224	locTotalMemory += double(locMemory);
225
226	//DKinFitConstraint_P4
227	locBin = dFactoryPoolBinMap["DKinFitConstraint_P4"];
228	// locNumObjectsMap[locBin] = locKinFitResultsFactory->Get_KinFitConstraintP4PoolSize();
229	locNumObjectsMap[locBin] = 1;
230	locMemory = sizeof(DKinFitConstraint_P4)*locNumObjectsMap[locBin];
231	locMemoryMap[locBin] = locMemory;
232	locTotalMemory += double(locMemory);
233
234	//TMatrixDSym_KinFitter
235	locBin = dFactoryPoolBinMap["TMatrixDSym_KinFitter"];
236	// locNumObjectsMap[locBin] = locKinFitResultsFactory->Get_MatrixDSymPoolSize() + locKinFitResultsFactory->Get_LargeMatrixDSymPoolSize();
237	locNumObjectsMap[locBin] = 1;
238	// locMemory = (sizeof(TMatrixDSym) + 778)*locKinFitResultsFactory->Get_MatrixDSymPoolSize(); //assume 7x7 matrix of doubles (8)
239	// locMemory += (sizeof(TMatrixDSym) + 30308)*locKinFitResultsFactory->Get_LargeMatrixDSymPoolSize(); //assume 30x30 matrix of doubles (8)
240	locMemoryMap[locBin] = locMemory;
241	locTotalMemory += double(locMemory);
242
243	//DParticleComboStep
244	locBin = dFactoryPoolBinMap["DParticleComboStep"];
245	locBaseFactory = locEventLoop->GetFactory("DParticleCombo", "");
246	DParticleCombo_factory* locParticleComboFactory = static_cast<DParticleCombo_factory*>(locBaseFactory);
247	locNumObjectsMap[locBin] = locParticleComboFactory->Get_ParticleComboStepPoolSize();
248	locMemory = sizeof(DParticleComboStep)*locNumObjectsMap[locBin];
249	locMemoryMap[locBin] = locMemory;
250	locTotalMemory += double(locMemory);
251
252	//DKinematicData_Combo
253	locBin = dFactoryPoolBinMap["DKinematicData_Combo"];
254	locNumObjectsMap[locBin] = locParticleComboFactory->Get_KinematicDataPoolSize();
255	locMemory = sizeof(DKinematicData)*locNumObjectsMap[locBin];
256	locMemoryMap[locBin] = locMemory;
257	locTotalMemory += double(locMemory);
258	}
259	locTotalMemory /= (1024.0*1024.0); //convert to MB
260
261	map<int, size_t>::iterator locIterator = locNumObjectsMap.begin();
262	japp->RootWriteLock(); //ACQUIRE ROOT LOCK!!
263	{
264	++dEventCounter;
265	if(dEventCounter <= dMaxNumEvents)
266	{
267	for(; locIterator != locNumObjectsMap.end(); ++locIterator)
268	{
269	int locObjectBin = locIterator->first;
270
271	dHistMap_NumObjects[locObjectBin]->SetBinContent(dEventCounter, locNumObjectsMap[locObjectBin]);
272	dHist_NumObjects->SetBinContent(dEventCounter, locObjectBin, locNumObjectsMap[locObjectBin]);
273
274	dHistMap_Memory[locObjectBin]->SetBinContent(dEventCounter, locMemoryMap[locObjectBin]);
275	dHist_Memory->SetBinContent(dEventCounter, locObjectBin, locMemoryMap[locObjectBin]);
276	}
277	dHist_TotalMemory->SetBinContent(dEventCounter, locTotalMemory);
278	}
279	}
280	japp->RootUnLock(); //RELEASE ROOT LOCK!!
281
282	return true;
283	}
284
285	void DHistogramAction_Reconstruction::Initialize(JEventLoop* locEventLoop)
286	{
287	//Create any histograms/trees/etc. within a ROOT lock.
288	//This is so that when running multithreaded, only one thread is writing to the ROOT file at a time.
289
290	//When creating a reaction-independent action, only modify member variables within a ROOT lock.
291	//Objects created within a plugin (such as reaction-independent actions) can be accessed by many threads simultaneously.
292
293	string locHistName, locHistTitle;
294
295	bool locIsRESTEvent = (string(locEventLoop->GetJEvent().GetJEventSource()->className()) == string("DEventSourceREST"));
296
297	vector<const DMCThrown*> locMCThrowns;
298	locEventLoop->Get(locMCThrowns);
299
300	japp->RootWriteLock(); //ACQUIRE ROOT LOCK!!
301	{
302	//Required: Create a folder in the ROOT output file that will contain all of the output ROOT objects (if any) for this action.
303	//If another thread has already created the folder, it just changes to it.
304	CreateAndChangeTo_ActionDirectory();
305
306	//FCAL
307	locHistName = "FCALShowerYVsX";
308	dHist_FCALShowerYVsX = GetOrCreate_Histogram<TH2I>(locHistName, ";FCAL Shower X (cm);FCAL Shower Y (cm)", dNumFCALTOFXYBins, -130.0, 130.0, dNumFCALTOFXYBins, -130.0, 130.0);
309	locHistName = "FCALShowerEnergy";
310	dHist_FCALShowerEnergy = GetOrCreate_Histogram<TH1I>(locHistName, ";FCAL Shower Energy (GeV)", dNumShowerEnergyBins, dMinShowerEnergy, dMaxShowerEnergy);
311
312	//BCAL
313	locHistName = "BCALShowerEnergy";
314	dHist_BCALShowerEnergy = GetOrCreate_Histogram<TH1I>(locHistName, ";BCAL Shower Energy (GeV)", dNumShowerEnergyBins, dMinShowerEnergy, dMaxBCALP);
315	locHistName = "BCALShowerPhi";
316	dHist_BCALShowerPhi = GetOrCreate_Histogram<TH1I>(locHistName, ";BCAL Shower #phi#circ", dNumPhiBins, dMinPhi, dMaxPhi);
317	locHistName = "BCALShowerPhiVsZ";
318	dHist_BCALShowerPhiVsZ = GetOrCreate_Histogram<TH2I>(locHistName, ";BCAL Shower Z (cm);BCAL Shower #phi#circ", dNum2DBCALZBins, 0.0, 450.0, dNum2DPhiBins, dMinPhi, dMaxPhi);
319
320	//TOF
321	locHistName = "TOFPointEnergy";
322	dHist_TOFPointEnergy = GetOrCreate_Histogram<TH1I>(locHistName, ";TOF Point Energy (MeV)", dNumHitEnergyBins, dMinHitEnergy, dMaxHitEnergy);
323	locHistName = "TOFPointYVsX";
324	dHist_TOFPointYVsX = GetOrCreate_Histogram<TH2I>(locHistName, ";TOF Point X (cm);TOF Point Y (cm)", dNumFCALTOFXYBins, -130.0, 130.0, dNumFCALTOFXYBins, -130.0, 130.0);
325
326	//SC
327	locHistName = "SCHitSector";
328	dHist_SCHitSector = GetOrCreate_Histogram<TH1I>(locHistName, ";SC Hit Sector", 30, 0.5, 30.5);
329	locHistName = "SCHitEnergy";
330	dHist_SCHitEnergy = GetOrCreate_Histogram<TH1I>(locHistName, ";SC Hit Energy (MeV)", dNumHitEnergyBins, dMinHitEnergy, dMaxHitEnergy);
331	locHistName = "SCHitEnergyVsSector";
332	dHist_SCHitEnergyVsSector = GetOrCreate_Histogram<TH2I>(locHistName, ";SC Hit Sector;SC Hit Energy (MeV)", 30, 0.5, 30.5, dNum2DHitEnergyBins, dMinHitEnergy, dMaxHitEnergy);
333
334	//TRACKING
335	CreateAndChangeTo_Directory("Tracking", "Tracking");
336	locHistName = "NumDCHitsPerTrack";
337	dHist_NumDCHitsPerTrack = GetOrCreate_Histogram<TH1I>(locHistName, ";# Track Hits", 50, 0.5, 50.5);
338	locHistName = "NumDCHitsPerTrackVsTheta";
339	dHist_NumDCHitsPerTrackVsTheta = GetOrCreate_Histogram<TH2I>(locHistName, ";#theta#circ;# Track Hits", dNum2DThetaBins, dMinTheta, dMaxTheta, 46, 4.5, 50.5);
340	locHistName = "TrackingFOM_WireBased";
341	dHist_TrackingFOM_WireBased = GetOrCreate_Histogram<TH1I>(locHistName, ";Confidence Level", dNumFOMBins, 0.0, 1.0);
342	locHistName = "TrackingFOM";
343	dHist_TrackingFOM = GetOrCreate_Histogram<TH1I>(locHistName, ";Confidence Level", dNumFOMBins, 0.0, 1.0);
344	locHistName = "TrackingFOMVsP";
345	dHist_TrackingFOMVsP = GetOrCreate_Histogram<TH2I>(locHistName, ";p (GeV/c);Confidence Level", dNum2DPBins, dMinP, dMaxP, dNum2DFOMBins, 0.0, 1.0);
346	locHistName = "TrackingFOMVsTheta";
347	dHist_TrackingFOMVsTheta = GetOrCreate_Histogram<TH2I>(locHistName, ";#theta#circ;Confidence Level", dNum2DThetaBins, dMinTheta, dMaxTheta, dNum2DFOMBins, 0.0, 1.0);
348	locHistName = "TrackingFOMVsNumHits";
349	dHist_TrackingFOMVsNumHits = GetOrCreate_Histogram<TH2I>(locHistName, ";# Track Hits;Confidence Level", 46, 4.5, 50.5, dNum2DFOMBins, 0.0, 1.0);
350
351	if(!locIsRESTEvent)
352	{
353	locHistName = "CDCRingVsTheta_Candidates";
354	dHist_CDCRingVsTheta_Candidates = GetOrCreate_Histogram<TH2I>(locHistName, "Hits on Track Candidates;#theta#circ;CDC Ring", dNum2DThetaBins, dMinTheta, dMaxTheta, 28, 0.5, 28.5);
355	locHistName = "CDCRingVsTheta_WireBased";
356	dHist_CDCRingVsTheta_WireBased = GetOrCreate_Histogram<TH2I>(locHistName, "Hits on Wire-Based Tracks;#theta#circ;CDC Ring", dNum2DThetaBins, dMinTheta, dMaxTheta, 28, 0.5, 28.5);
357	}
358
359	locHistName = "CDCRingVsTheta_TimeBased";
360	dHist_CDCRingVsTheta_TimeBased = GetOrCreate_Histogram<TH2I>(locHistName, "Hits on Time-Based Tracks;#theta#circ;CDC Ring", dNum2DThetaBins, dMinTheta, dMaxTheta, 28, 0.5, 28.5);
361	locHistName = "CDCRingVsTheta_TimeBased_GoodTrackFOM";
362	dHist_CDCRingVsTheta_TimeBased_GoodTrackFOM = GetOrCreate_Histogram<TH2I>(locHistName, "Hits on Good FOM Time-Based Tracks;#theta#circ;CDC Ring", dNum2DThetaBins, dMinTheta, dMaxTheta, 28, 0.5, 28.5);
363
364	if(!locIsRESTEvent)
365	{
366	locHistName = "FDCPlaneVsTheta_Candidates";
367	dHist_FDCPlaneVsTheta_Candidates = GetOrCreate_Histogram<TH2I>(locHistName, "Hits on Track Candidates;p (GeV/c);FDC Plane", dNum2DThetaBins, dMinTheta, dMaxTheta, 24, 0.5, 24.5);
368	locHistName = "FDCPlaneVsTheta_WireBased";
369	dHist_FDCPlaneVsTheta_WireBased = GetOrCreate_Histogram<TH2I>(locHistName, "Hits on Wire-Based Tracks;p (GeV/c);FDC Plane", dNum2DThetaBins, dMinTheta, dMaxTheta, 24, 0.5, 24.5);
370	}
371
372	locHistName = "FDCPlaneVsTheta_TimeBased";
373	dHist_FDCPlaneVsTheta_TimeBased = GetOrCreate_Histogram<TH2I>(locHistName, "Hits on Time-Based Tracks;p (GeV/c);FDC Plane", dNum2DThetaBins, dMinTheta, dMaxTheta, 24, 0.5, 24.5);
374	locHistName = "FDCPlaneVsTheta_TimeBased_GoodTrackFOM";
375	dHist_FDCPlaneVsTheta_TimeBased_GoodTrackFOM = GetOrCreate_Histogram<TH2I>(locHistName, "Hits on Good FOM Time-Based Tracks;p (GeV/c);FDC Plane", dNum2DThetaBins, dMinTheta, dMaxTheta, 24, 0.5, 24.5);
376
377	if(!locMCThrowns.empty())
378	{
379	locHistName = "MCMatchedHitsVsTheta";
380	dHist_MCMatchedHitsVsTheta = GetOrCreate_Histogram<TH2I>(locHistName, "Fraction of Track Hits Matched to MC;#theta#circ;", dNum2DThetaBins, dMinTheta, dMaxTheta, 100, 0.0, 1.0);
381	locHistName = "MCMatchedHitsVsP";
382	dHist_MCMatchedHitsVsP = GetOrCreate_Histogram<TH2I>(locHistName, "Fraction of Track Hits Matched to MC;p (GeV/c);", dNum2DPBins, dMinP, dMaxP, 100, 0.0, 1.0);
383	}
384
385	for(int locCharge = -1; locCharge <= 1; locCharge += 2)
386	{
387	string locParticleROOTName = (locCharge == -1) ? "#it{q}^{-}" : "#it{q}^{+}";
388	string locParticleName = (locCharge == -1) ? "q-" : "q+";
389	CreateAndChangeTo_Directory(locParticleName, locParticleName);
390
391	if(!locIsRESTEvent)
392	{
393	// PVsTheta Track Candidates
394	locHistName = string("PVsTheta_Candidates_") + locParticleName;
395	locHistTitle = locParticleROOTName + string(" Track Candidates;#theta#circ;p (GeV/c)");
396	dHistMap_PVsTheta_Candidates[locCharge] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DThetaBins, dMinTheta, dMaxTheta, dNum2DPBins, dMinP, dMaxP);
397
398	// PVsTheta Wire-Based Tracks
399	locHistName = string("PVsTheta_WireBased_") + locParticleName;
400	locHistTitle = locParticleROOTName + string(" Wire-Based Tracks;#theta#circ;p (GeV/c)");
401	dHistMap_PVsTheta_WireBased[locCharge] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DThetaBins, dMinTheta, dMaxTheta, dNum2DPBins, dMinP, dMaxP);
402	}
403
404	// PVsTheta Time-Based Tracks
405	locHistName = string("PVsTheta_TimeBased_") + locParticleName;
406	locHistTitle = locParticleROOTName + string(" Time-Based Tracks;#theta#circ;p (GeV/c)");
407	dHistMap_PVsTheta_TimeBased[locCharge] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DThetaBins, dMinTheta, dMaxTheta, dNum2DPBins, dMinP, dMaxP);
408
409	// PVsTheta Time-Based Tracks Good Track FOM
410	locHistName = string("PVsTheta_TimeBased_GoodTrackFOM_") + locParticleName;
411	locHistTitle = locParticleROOTName + string(" Time-Based Tracks, Good Tracking FOM;#theta#circ;p (GeV/c)");
412	dHistMap_PVsTheta_TimeBased_GoodTrackFOM[locCharge] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DThetaBins, dMinTheta, dMaxTheta, dNum2DPBins, dMinP, dMaxP);
413
414	// PVsTheta Time-Based Tracks Low Track FOM
415	locHistName = string("PVsTheta_TimeBased_LowTrackFOM_") + locParticleName;
416	locHistTitle = locParticleROOTName + string(" Time-Based Tracks, Low Tracking FOM;#theta#circ;p (GeV/c)");
417	dHistMap_PVsTheta_TimeBased_LowTrackFOM[locCharge] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DThetaBins, dMinTheta, dMaxTheta, dNum2DPBins, dMinP, dMaxP);
418
419	// PVsTheta Time-Based Tracks High Track FOM
420	locHistName = string("PVsTheta_TimeBased_HighTrackFOM_") + locParticleName;
421	locHistTitle = locParticleROOTName + string(" Time-Based Tracks, High Tracking FOM;#theta#circ;p (GeV/c)");
422	dHistMap_PVsTheta_TimeBased_HighTrackFOM[locCharge] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DThetaBins, dMinTheta, dMaxTheta, dNum2DPBins, dMinP, dMaxP);
423
424	if(!locIsRESTEvent)
425	{
426	// PVsTheta: Good Wire-Based, Good Time-Based
427	locHistName = string("PVsTheta_GoodWireBased_GoodTimeBased_") + locParticleName;
428	locHistTitle = locParticleROOTName + string(" Good Wire-Based, Good Time-Based;#theta#circ;p (GeV/c)");
429	dHistMap_PVsTheta_GoodWireBased_GoodTimeBased[locCharge] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DThetaBins, dMinTheta, dMaxTheta, dNum2DPBins, dMinP, dMaxP);
430
431	// PVsTheta: Good Wire-Based, Bad Time-Based
432	locHistName = string("PVsTheta_GoodWireBased_BadTimeBased_") + locParticleName;
433	locHistTitle = locParticleROOTName + string(" Good Wire-Based, Bad Time-Based;#theta#circ;p (GeV/c)");
434	dHistMap_PVsTheta_GoodWireBased_BadTimeBased[locCharge] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DThetaBins, dMinTheta, dMaxTheta, dNum2DPBins, dMinP, dMaxP);
435	}
436
437	gDirectory(TDirectory::CurrentDirectory())->cd(".."); //end of charge
438	}
439	gDirectory(TDirectory::CurrentDirectory())->cd(".."); //End of "Tracking"
440
441	//Return to the base directory
442	ChangeTo_BaseDirectory();
443	}
444	japp->RootUnLock(); //RELEASE ROOT LOCK!!
445	}
446
447	bool DHistogramAction_Reconstruction::Perform_Action(JEventLoop* locEventLoop, const DParticleCombo* locParticleCombo)
448	{
449	//Expect locParticleCombo to be NULL since this is a reaction-independent action.
450
451	//Optional: Quit the action if it has already been executed this event (else may result in double-counting when filling histograms)
452	if(Get_NumPreviousParticleCombos() != 0)
453	return true;
454
455	bool locIsRESTEvent = (string(locEventLoop->GetJEvent().GetJEventSource()->className()) == string("DEventSourceREST"));
456
457	vector<const DBCALShower*> locBCALShowers;
458	locEventLoop->Get(locBCALShowers);
459
460	vector<const DFCALShower*> locFCALShowers;
461	locEventLoop->Get(locFCALShowers);
462
463	vector<const DTOFPoint*> locTOFPoints;
464	locEventLoop->Get(locTOFPoints);
465
466	vector<const DSCHit*> locSCHits;
467	locEventLoop->Get(locSCHits);
468
469	const DDetectorMatches* locDetectorMatches = NULL__null;
470	locEventLoop->GetSingle(locDetectorMatches);
471
472	vector<const DTrackTimeBased*> locTrackTimeBasedVector;
473	locEventLoop->Get(locTrackTimeBasedVector);
474
475	vector<const DMCThrownMatching*> locMCThrownMatchingVector;
476	locEventLoop->Get(locMCThrownMatchingVector);
477
478	vector<const DMCThrown*> locMCThrowns;
479	locEventLoop->Get(locMCThrowns, "FinalState");
480
481	const DParticleID* locParticleID = NULL__null;
482	locEventLoop->GetSingle(locParticleID);
483
484	const DDetectorMatches* locDetectorMatches_WireBased = NULL__null;
485	vector<const DTrackCandidate*> locTrackCandidates;
486	vector<const DTrackWireBased*> locTrackWireBasedVector;
487	if(!locIsRESTEvent)
488	{
489	vector<const DDetectorMatches*> locDetectorMatchesVector_WireBased;
490	locEventLoop->Get(locDetectorMatchesVector_WireBased, "WireBased");
491	if(!locDetectorMatchesVector_WireBased.empty())
492	locDetectorMatches_WireBased = locDetectorMatchesVector_WireBased[0];
493	locEventLoop->Get(locTrackCandidates);
494	locEventLoop->Get(locTrackWireBasedVector);
495	}
496
497	//select the best DTrackWireBased for each track: use best tracking FOM
498	map<JObject::oid_t, const DTrackWireBased*> locBestTrackWireBasedMap; //lowest tracking FOM for each candidate id
499	for(size_t loc_i = 0; loc_i < locTrackWireBasedVector.size(); ++loc_i)
500	{
501	JObject::oid_t locCandidateID = locTrackWireBasedVector[loc_i]->candidateid;
502	if(locBestTrackWireBasedMap.find(locCandidateID) == locBestTrackWireBasedMap.end())
503	locBestTrackWireBasedMap[locCandidateID] = locTrackWireBasedVector[loc_i];
504	else if(locTrackWireBasedVector[loc_i]->FOM > locBestTrackWireBasedMap[locCandidateID]->FOM)
505	locBestTrackWireBasedMap[locCandidateID] = locTrackWireBasedVector[loc_i];
506	}
507
508	//select the best DTrackTimeBased for each track: use best tracking FOM
509	//also, make map from WBT -> TBT (if not rest)
510	map<JObject::oid_t, const DTrackTimeBased*> locBestTrackTimeBasedMap; //lowest tracking FOM for each candidate id
511	map<const DTrackWireBased, const DTrackTimeBased> locWireToTimeBasedTrackMap;
512	for(size_t loc_i = 0; loc_i < locTrackTimeBasedVector.size(); ++loc_i)
513	{
514	JObject::oid_t locCandidateID = locTrackTimeBasedVector[loc_i]->candidateid;
515	if(locBestTrackTimeBasedMap.find(locCandidateID) == locBestTrackTimeBasedMap.end())
516	locBestTrackTimeBasedMap[locCandidateID] = locTrackTimeBasedVector[loc_i];
517	else if(locTrackTimeBasedVector[loc_i]->FOM > locBestTrackTimeBasedMap[locCandidateID]->FOM)
518	locBestTrackTimeBasedMap[locCandidateID] = locTrackTimeBasedVector[loc_i];
519	if(locIsRESTEvent)
520	continue;
521	const DTrackWireBased* locTrackWireBased = NULL__null;
522	locTrackTimeBasedVector[loc_i]->GetSingle(locTrackWireBased);
523	locWireToTimeBasedTrackMap[locTrackWireBased] = locTrackTimeBasedVector[loc_i];
524	}
525
526	//Fill Histograms
527	japp->RootWriteLock(); //ACQUIRE ROOT LOCK!!
528	{
529	for(size_t loc_i = 0; loc_i < locFCALShowers.size(); ++loc_i)
530	{
531	dHist_FCALShowerEnergy->Fill(locFCALShowers[loc_i]->getEnergy());
532	dHist_FCALShowerYVsX->Fill(locFCALShowers[loc_i]->getPosition().X(), locFCALShowers[loc_i]->getPosition().Y());
533	}
534
535	for(size_t loc_i = 0; loc_i < locBCALShowers.size(); ++loc_i)
536	{
537	dHist_BCALShowerEnergy->Fill(locBCALShowers[loc_i]->E);
538
539	DVector3 locBCALPosition(locBCALShowers[loc_i]->x, locBCALShowers[loc_i]->y, locBCALShowers[loc_i]->z);
540	double locBCALPhi = locBCALPosition.Phi()*180.0/TMath::Pi();
541	dHist_BCALShowerPhi->Fill(locBCALPhi);
542	dHist_BCALShowerPhiVsZ->Fill(locBCALPosition.Z(), locBCALPhi);
543	}
544
545	for(size_t loc_i = 0; loc_i < locTOFPoints.size(); ++loc_i)
546	{
547	dHist_TOFPointEnergy->Fill(locTOFPoints[loc_i]->dE*1.0E3);
548	dHist_TOFPointYVsX->Fill(locTOFPoints[loc_i]->pos.X(), locTOFPoints[loc_i]->pos.Y());
549	}
550
551	for(size_t loc_i = 0; loc_i < locSCHits.size(); ++loc_i)
552	{
553	dHist_SCHitSector->Fill(locSCHits[loc_i]->sector);
554	dHist_SCHitEnergy->Fill(locSCHits[loc_i]->dE*1.0E3);
555	dHist_SCHitEnergyVsSector->Fill(locSCHits[loc_i]->sector, locSCHits[loc_i]->dE*1.0E3);
556	}
557
558	for(size_t loc_i = 0; loc_i < locTrackCandidates.size(); ++loc_i)
559	{
560	int locCharge = (locTrackCandidates[loc_i]->charge() > 0.0) ? 1 : -1;
561	double locTheta = locTrackCandidates[loc_i]->momentum().Theta()*180.0/TMath::Pi();
562	double locP = locTrackCandidates[loc_i]->momentum().Mag();
563	dHistMap_PVsTheta_Candidates[locCharge]->Fill(locTheta, locP);
564
565	set<int> locCDCRings;
566	locParticleID->Get_CDCRings(locTrackCandidates[loc_i]->dCDCRings, locCDCRings);
567	for(set<int>::iterator locIterator = locCDCRings.begin(); locIterator != locCDCRings.end(); ++locIterator)
568	dHist_CDCRingVsTheta_Candidates->Fill(locTheta, *locIterator);
569
570	set<int> locFDCPlanes;
571	locParticleID->Get_FDCPlanes(locTrackCandidates[loc_i]->dFDCPlanes, locFDCPlanes);
572	for(set<int>::iterator locIterator = locFDCPlanes.begin(); locIterator != locFDCPlanes.end(); ++locIterator)
573	dHist_FDCPlaneVsTheta_Candidates->Fill(locTheta, *locIterator);
574	}
575
576	map<JObject::oid_t, const DTrackWireBased*>::iterator locWireBasedIterator = locBestTrackWireBasedMap.begin();
577	for(; locWireBasedIterator != locBestTrackWireBasedMap.end(); ++locWireBasedIterator)
578	{
579	const DTrackWireBased* locTrackWireBased = locWireBasedIterator->second;
580	int locCharge = (locTrackWireBased->charge() > 0.0) ? 1 : -1;
581	double locTheta = locTrackWireBased->momentum().Theta()*180.0/TMath::Pi();
582	double locP = locTrackWireBased->momentum().Mag();
583	dHistMap_PVsTheta_WireBased[locCharge]->Fill(locTheta, locP);
584
585	set<int> locCDCRings;
586	locParticleID->Get_CDCRings(locTrackWireBased->dCDCRings, locCDCRings);
587	for(set<int>::iterator locIterator = locCDCRings.begin(); locIterator != locCDCRings.end(); ++locIterator)
588	dHist_CDCRingVsTheta_WireBased->Fill(locTheta, *locIterator);
589
590	set<int> locFDCPlanes;
591	locParticleID->Get_FDCPlanes(locTrackWireBased->dFDCPlanes, locFDCPlanes);
592	for(set<int>::iterator locIterator = locFDCPlanes.begin(); locIterator != locFDCPlanes.end(); ++locIterator)
593	dHist_FDCPlaneVsTheta_WireBased->Fill(locTheta, *locIterator);
594
595	dHist_TrackingFOM_WireBased->Fill(locTrackWireBased->FOM);
596	}
597
598	map<JObject::oid_t, const DTrackTimeBased*>::iterator locTimeBasedIterator = locBestTrackTimeBasedMap.begin();
599	for(; locTimeBasedIterator != locBestTrackTimeBasedMap.end(); ++locTimeBasedIterator)
600	{
601	const DTrackTimeBased* locTrackTimeBased = locTimeBasedIterator->second;
602	int locCharge = (locTrackTimeBased->charge() > 0.0) ? 1 : -1;
603	double locTheta = locTrackTimeBased->momentum().Theta()*180.0/TMath::Pi();
604	double locP = locTrackTimeBased->momentum().Mag();
605
606	dHistMap_PVsTheta_TimeBased[locCharge]->Fill(locTheta, locP);
607	dHist_NumDCHitsPerTrack->Fill(locTrackTimeBased->Ndof + 5);
608	dHist_NumDCHitsPerTrackVsTheta->Fill(locTheta, locTrackTimeBased->Ndof + 5);
609
610	dHist_TrackingFOM->Fill(locTrackTimeBased->FOM);
611	dHist_TrackingFOMVsTheta->Fill(locTheta, locTrackTimeBased->FOM);
612	dHist_TrackingFOMVsP->Fill(locP, locTrackTimeBased->FOM);
613	dHist_TrackingFOMVsNumHits->Fill(locTrackTimeBased->Ndof + 5, locTrackTimeBased->FOM);
614
615	set<int> locCDCRings;
616	locParticleID->Get_CDCRings(locTrackTimeBased->dCDCRings, locCDCRings);
617	for(set<int>::iterator locIterator = locCDCRings.begin(); locIterator != locCDCRings.end(); ++locIterator)
618	{
619	dHist_CDCRingVsTheta_TimeBased->Fill(locTheta, *locIterator);
620	if(locTrackTimeBased->FOM > dGoodTrackFOM)
621	dHist_CDCRingVsTheta_TimeBased_GoodTrackFOM->Fill(locTheta, *locIterator);
622	}
623
624	set<int> locFDCPlanes;
625	locParticleID->Get_FDCPlanes(locTrackTimeBased->dFDCPlanes, locFDCPlanes);
626	for(set<int>::iterator locIterator = locFDCPlanes.begin(); locIterator != locFDCPlanes.end(); ++locIterator)
627	{
628	dHist_FDCPlaneVsTheta_TimeBased->Fill(locTheta, *locIterator);
629	if(locTrackTimeBased->FOM > dGoodTrackFOM)
630	dHist_FDCPlaneVsTheta_TimeBased_GoodTrackFOM->Fill(locTheta, *locIterator);
631	}
632
633	if(locTrackTimeBased->FOM > dGoodTrackFOM)
634	dHistMap_PVsTheta_TimeBased_GoodTrackFOM[locCharge]->Fill(locTheta, locP);
635	else
636	dHistMap_PVsTheta_TimeBased_LowTrackFOM[locCharge]->Fill(locTheta, locP);
637	if(locTrackTimeBased->FOM > dHighTrackFOM)
638	dHistMap_PVsTheta_TimeBased_HighTrackFOM[locCharge]->Fill(locTheta, locP);
639	}
640
641	// If "Good" WBT, see if TBT is good
642	locWireBasedIterator = locBestTrackWireBasedMap.begin();
643	for(; locWireBasedIterator != locBestTrackWireBasedMap.end(); ++locWireBasedIterator)
644	{
645	if(locDetectorMatches_WireBased == NULL__null)
646	continue;
647	const DTrackWireBased* locTrackWireBased = locWireBasedIterator->second;
648	if(locTrackWireBased->FOM < dGoodTrackFOM)
649	continue; //no good
650	if(!locDetectorMatches_WireBased->Get_IsMatchedToHit(locTrackWireBased))
651	continue; //no good
652
653	int locCharge = (locTrackWireBased->charge() > 0.0) ? 1 : -1;
654	double locTheta = locTrackWireBased->momentum().Theta()*180.0/TMath::Pi();
655	double locP = locTrackWireBased->momentum().Mag();
656
657	map<const DTrackWireBased, const DTrackTimeBased>::iterator locReconIterator = locWireToTimeBasedTrackMap.find(locTrackWireBased);
658	if(locReconIterator == locWireToTimeBasedTrackMap.end())
659	{
660	dHistMap_PVsTheta_GoodWireBased_BadTimeBased[locCharge]->Fill(locTheta, locP);
661	continue; //no time-based
662	}
663
664	const DTrackTimeBased* locTrackTimeBased = locReconIterator->second;
665	if((locTrackTimeBased->FOM < dGoodTrackFOM) \|\| (!locDetectorMatches->Get_IsMatchedToHit(locTrackTimeBased)))
666	dHistMap_PVsTheta_GoodWireBased_BadTimeBased[locCharge]->Fill(locTheta, locP);
667	else
668	dHistMap_PVsTheta_GoodWireBased_GoodTimeBased[locCharge]->Fill(locTheta, locP);
669	}
670
671	for(size_t loc_i = 0; loc_i < locMCThrowns.size(); ++loc_i)
672	{
673	if(fabs(locMCThrowns[loc_i]->charge()) < 0.9)
674	continue;
675
676	double locMatchFOM;
677	const DChargedTrackHypothesis* locChargedTrackHypothesis = locMCThrownMatchingVector[0]->Get_MatchingChargedHypothesis(locMCThrowns[loc_i], locMatchFOM);
678	if(locChargedTrackHypothesis == NULL__null)
679	continue;
680
681	const DTrackTimeBased* locTrackTimeBased = NULL__null;
682	locChargedTrackHypothesis->GetSingle(locTrackTimeBased);
683
684	double locHitFraction = 1.0*locTrackTimeBased->dNumHitsMatchedToThrown/(locTrackTimeBased->Ndof + 5);
685	dHist_MCMatchedHitsVsTheta->Fill(locTrackTimeBased->momentum().Theta()*180.0/TMath::Pi(), locHitFraction);
686	dHist_MCMatchedHitsVsP->Fill(locTrackTimeBased->momentum().Mag(), locHitFraction);
687	}
688	}
689	japp->RootUnLock(); //RELEASE ROOT LOCK!!
690
691	return true; //return false if you want to use this action to apply a cut (and it fails the cut!)
692	}
693
694	void DHistogramAction_DetectorMatching::Initialize(JEventLoop* locEventLoop)
695	{
696	//Create any histograms/trees/etc. within a ROOT lock.
697	//This is so that when running multithreaded, only one thread is writing to the ROOT file at a time.
698
699	//When creating a reaction-independent action, only modify member variables within a ROOT lock.
700	//Objects created within a plugin (such as reaction-independent actions) can be accessed by many threads simultaneously.
701	string locHistName, locHistTitle;
702
703	bool locIsRESTEvent = (string(locEventLoop->GetJEvent().GetJEventSource()->className()) == string("DEventSourceREST"));
704
705	japp->RootWriteLock(); //ACQUIRE ROOT LOCK!!
706	{
707	//Required: Create a folder in the ROOT output file that will contain all of the output ROOT objects (if any) for this action.
708	//If another thread has already created the folder, it just changes to it.
709	CreateAndChangeTo_ActionDirectory();
710
711	//Loop over filling for time-based and wire-based tracks
712	for(unsigned int locDummy = 0; locDummy < 2; ++locDummy)
713	{
714	bool locIsTimeBased = (locDummy == 0);
715	if(locIsRESTEvent && (!locIsTimeBased))
716	continue;
717
718	string locDirectoryName = locIsTimeBased ? "TimeBased" : "WireBased";
719	CreateAndChangeTo_Directory(locDirectoryName.c_str(), locDirectoryName.c_str());
720	string locTrackString = locIsTimeBased ? "Time-Based Tracks" : "Wire-Based Tracks";
721
722	//Kinematics of has (no) hit
723	vector<DetectorSystem_t> locDetectorSystems;
724	locDetectorSystems.push_back(SYS_START); locDetectorSystems.push_back(SYS_BCAL);
725	locDetectorSystems.push_back(SYS_TOF); locDetectorSystems.push_back(SYS_FCAL);
726	for(size_t loc_i = 0; loc_i < locDetectorSystems.size(); ++loc_i)
727	{
728	DetectorSystem_t locSystem = locDetectorSystems[loc_i];
729
730	double locMaxTheta = ((locSystem == SYS_FCAL) \|\| (locSystem == SYS_TOF)) ? 12.0 : dMaxTheta;
731	double locMaxP = (locSystem == SYS_BCAL) ? 3.0 : dMaxP;
732
733	string locSystemName = SystemName(locSystem);
734	if(locSystemName == "ST")
735	locSystemName = "SC";
736	string locDirName = locSystemName;
737	if(locSystemName == "TOF")
738	locDirName = "TOFPoint";
739
740	CreateAndChangeTo_Directory(locDirName, locDirName);
741
742	// PVsTheta Has Hit
743	locHistName = "PVsTheta_HasHit";
744	locHistTitle = locTrackString + string(", Has Other Match, ") + locSystemName + string(" Has Hit;#theta#circ;p (GeV/c)");
745	dHistMap_PVsTheta_HasHit[locSystem][locIsTimeBased] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DThetaBins, dMinTheta, locMaxTheta, dNum2DPBins, dMinP, locMaxP);
746
747	// PVsTheta Has No Hit
748	locHistName = "PVsTheta_NoHit";
749	locHistTitle = locTrackString + string(", Has Other Match, ") + locSystemName + string(" No Hit;#theta#circ;p (GeV/c)");
750	dHistMap_PVsTheta_NoHit[locSystem][locIsTimeBased] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DThetaBins, dMinTheta, locMaxTheta, dNum2DPBins, dMinP, locMaxP);
751
752	// PhiVsTheta Has Hit
753	locHistName = "PhiVsTheta_HasHit";
754	locHistTitle = locTrackString + string(", Has Other Match, ") + locSystemName + string(" Has Hit;#theta#circ;#phi#circ");
755	dHistMap_PhiVsTheta_HasHit[locSystem][locIsTimeBased] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DThetaBins, dMinTheta, locMaxTheta, dNum2DPhiBins, dMinPhi, dMaxPhi);
756
757	// PhiVsTheta Has No Hit
758	locHistName = "PhiVsTheta_NoHit";
759	locHistTitle = locTrackString + string(", Has Other Match, ") + locSystemName + string(" No Hit;#theta#circ;#phi#circ");
760	dHistMap_PhiVsTheta_NoHit[locSystem][locIsTimeBased] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DThetaBins, dMinTheta, locMaxTheta, dNum2DPhiBins, dMinPhi, dMaxPhi);
761
762	gDirectory(TDirectory::CurrentDirectory())->cd("..");
763	}
764
765	//SC
766	CreateAndChangeTo_Directory("SC", "SC");
767	locHistName = "SCPaddleVsTheta_HasHit";
768	locHistTitle = locTrackString + string(", Has Other Match, SC Has Hit;#theta#circ;Projected SC Paddle");
769	dHistMap_SCPaddleVsTheta_HasHit[locIsTimeBased] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DThetaBins, dMinTheta, dMaxTheta, 30, 0.5, 30.5);
770
771	locHistName = "SCPaddleVsTheta_NoHit";
772	locHistTitle = locTrackString + string(", Has Other Match, SC No Hit;#theta#circ;Projected SC Paddle");
773	dHistMap_SCPaddleVsTheta_NoHit[locIsTimeBased] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DThetaBins, dMinTheta, dMaxTheta, 30, 0.5, 30.5);
774
775	locHistName = "SCPaddle_BarrelRegion_HasHit";
776	locHistTitle = locTrackString + string(", Has Other Match, SC Barrel Region Has Hit;Projected SC Paddle");
777	dHistMap_SCPaddle_BarrelRegion_HasHit[locIsTimeBased] = GetOrCreate_Histogram<TH1I>(locHistName, locHistTitle, 30, 0.5, 30.5);
778
779	locHistName = "SCPaddle_BarrelRegion_NoHit";
780	locHistTitle = locTrackString + string(", Has Other Match, SC Barrel Region No Hit;Projected SC Paddle");
781	dHistMap_SCPaddle_BarrelRegion_NoHit[locIsTimeBased] = GetOrCreate_Histogram<TH1I>(locHistName, locHistTitle, 30, 0.5, 30.5);
782
783	locHistName = "SCPaddle_NoseRegion_HasHit";
784	locHistTitle = locTrackString + string(", Has Other Match, SC Front Region Has Hit;Projected SC Paddle");
785	dHistMap_SCPaddle_NoseRegion_HasHit[locIsTimeBased] = GetOrCreate_Histogram<TH1I>(locHistName, locHistTitle, 30, 0.5, 30.5);
786
787	locHistName = "SCPaddle_NoseRegion_NoHit";
788	locHistTitle = locTrackString + string(", Has Other Match, SC Front Region No Hit;Projected SC Paddle");
789	dHistMap_SCPaddle_NoseRegion_NoHit[locIsTimeBased] = GetOrCreate_Histogram<TH1I>(locHistName, locHistTitle, 30, 0.5, 30.5);
790
791	locHistName = "SCTrackDeltaPhiVsP";
792	locHistTitle = locTrackString + string(";p (GeV/c);SC / Track #Delta#phi#circ");
793	dHistMap_SCTrackDeltaPhiVsP[locIsTimeBased] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNum2DDeltaPhiBins, dSCMatchMinDeltaPhi, dSCMatchMaxDeltaPhi);
794
795	locHistName = "SCTrackDeltaPhiVsTheta";
796	locHistTitle = locTrackString + string(";#theta#circ;SC / Track #Delta#phi#circ");
797	dHistMap_SCTrackDeltaPhiVsTheta[locIsTimeBased] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DThetaBins, dMinTheta, dMaxTheta, dNum2DDeltaPhiBins, dSCMatchMinDeltaPhi, dSCMatchMaxDeltaPhi);
798	gDirectory(TDirectory::CurrentDirectory())->cd("..");
799
800	//TOFPaddle
801	CreateAndChangeTo_Directory("TOFPaddle", "TOFPaddle");
802	locHistName = "VerticalPaddleTrackDeltaX";
803	locHistTitle = locTrackString + string(";Vertical TOF Paddle / Track \|#DeltaX\| (cm)");
804	dHistMap_TOFPaddleTrackDeltaX[locIsTimeBased] = GetOrCreate_Histogram<TH1I>(locHistName, locHistTitle, dNumTrackDOCABins/2, 0.0, dMaxTrackMatchDOCA);
805
806	locHistName = "HorizontalPaddleTrackDeltaY";
807	locHistTitle = locTrackString + string(";Horizontal TOF Paddle / Track \|#DeltaY\| (cm)");
808	dHistMap_TOFPaddleTrackDeltaY[locIsTimeBased] = GetOrCreate_Histogram<TH1I>(locHistName, locHistTitle, dNumTrackDOCABins/2, 0.0, dMaxTrackMatchDOCA);
809
810	locHistName = "TrackYVsVerticalPaddle_HasHit";
811	locHistTitle = locTrackString + string(", Has Other Match, TOF Paddle Has Hit;Projected Vertical Paddle;Projected TOF Hit Y (cm)");
812	dHistMap_TOFPaddleTrackYVsVerticalPaddle_HasHit[locIsTimeBased] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, 44, 0.5, 44.5, dNumFCALTOFXYBins, -130.0, 130.0);
813
814	locHistName = "TrackYVsVerticalPaddle_NoHit";
815	locHistTitle = locTrackString + string(", Has Other Match, TOF Paddle No Hit;Projected Vertical Paddle;Projected TOF Hit Y (cm)");
816	dHistMap_TOFPaddleTrackYVsVerticalPaddle_NoHit[locIsTimeBased] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, 44, 0.5, 44.5, dNumFCALTOFXYBins, -130.0, 130.0);
817
818	locHistName = "HorizontalPaddleVsTrackX_HasHit";
819	locHistTitle = locTrackString + string(", Has Other Match, TOF Paddle Has Hit;Projected TOF Hit X (cm);Projected Horizontal Paddle");
820	dHistMap_TOFPaddleHorizontalPaddleVsTrackX_HasHit[locIsTimeBased] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNumFCALTOFXYBins, -130.0, 130.0, 44, 0.5, 44.5);
821
822	locHistName = "HorizontalPaddleVsTrackX_NoHit";
823	locHistTitle = locTrackString + string(", Has Other Match, TOF Paddle No Hit;Projected TOF Hit X (cm);Projected Horizontal Paddle");
824	dHistMap_TOFPaddleHorizontalPaddleVsTrackX_NoHit[locIsTimeBased] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNumFCALTOFXYBins, -130.0, 130.0, 44, 0.5, 44.5);
825	gDirectory(TDirectory::CurrentDirectory())->cd("..");
826
827	//TOFPoint
828	CreateAndChangeTo_Directory("TOFPoint", "TOFPoint");
829	locHistName = "TrackTOFYVsX_HasHit";
830	locHistTitle = locTrackString + string(", Has Other Match, TOF Has Hit;Projected TOF Hit X (cm);Projected TOF Hit Y (cm)");
831	dHistMap_TrackTOFYVsX_HasHit[locIsTimeBased] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNumFCALTOFXYBins, -130.0, 130.0, dNumFCALTOFXYBins, -130.0, 130.0);
832
833	locHistName = "TrackTOFYVsX_NoHit";
834	locHistTitle = locTrackString + string(", Has Other Match, TOF No Hit;Projected TOF Hit X (cm);Projected TOF Hit Y (cm)");
835	dHistMap_TrackTOFYVsX_NoHit[locIsTimeBased] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNumFCALTOFXYBins, -130.0, 130.0, dNumFCALTOFXYBins, -130.0, 130.0);
836
837	locHistName = "TrackTOF2DPaddles_HasHit";
838	locHistTitle = locTrackString + string(", Has Other Match, TOF Has Hit;Projected Vertical TOF Paddle;Projected Horizontal TOF Paddle");
839	dHistMap_TrackTOF2DPaddles_HasHit[locIsTimeBased] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, 44, 0.5, 44.5, 44, 0.5, 44.5);
840
841	locHistName = "TrackTOF2DPaddles_NoHit";
842	locHistTitle = locTrackString + string(", Has Other Match, TOF No Hit;Projected Vertical TOF Paddle;Projected Horizontal TOF Paddle");
843	dHistMap_TrackTOF2DPaddles_NoHit[locIsTimeBased] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, 44, 0.5, 44.5, 44, 0.5, 44.5);
844
845	locHistName = "TOFTrackDistanceVsP";
846	locHistTitle = locTrackString + string(";p (GeV/c);TOF / Track Distance (cm)");
847	dHistMap_TOFPointTrackDistanceVsP[locIsTimeBased] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNum2DTrackDOCABins, dMinTrackDOCA, dMaxTrackMatchDOCA);
848
849	locHistName = "TOFTrackDistanceVsTheta";
850	locHistTitle = locTrackString + string(";#theta#circ;TOF / Track Distance (cm)");
851	dHistMap_TOFPointTrackDistanceVsTheta[locIsTimeBased] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DThetaBins, dMinTheta, 20.0, dNum2DTrackDOCABins, dMinTrackDOCA, dMaxTrackMatchDOCA);
852
853	locHistName = "TOFTrackDeltaXVsHorizontalPaddle";
854	locHistTitle = locTrackString + string(";TOF Horizontal Paddle;TOF / Track #DeltaX (cm)");
855	dHistMap_TOFPointTrackDeltaXVsHorizontalPaddle[locIsTimeBased] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, 44, 0.5, 44.5, dNum2DTrackDOCABins, -1.0*dMaxTrackMatchDOCA, dMaxTrackMatchDOCA);
856
857	locHistName = "TOFTrackDeltaXVsVerticalPaddle";
858	locHistTitle = locTrackString + string(";TOF Vertical Paddle;TOF / Track #DeltaX (cm)");
859	dHistMap_TOFPointTrackDeltaXVsVerticalPaddle[locIsTimeBased] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, 44, 0.5, 44.5, dNum2DTrackDOCABins, -1.0*dMaxTrackMatchDOCA, dMaxTrackMatchDOCA);
860
861	locHistName = "TOFTrackDeltaYVsHorizontalPaddle";
862	locHistTitle = locTrackString + string(";TOF Horizontal Paddle;TOF / Track #DeltaY (cm)");
863	dHistMap_TOFPointTrackDeltaYVsHorizontalPaddle[locIsTimeBased] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, 44, 0.5, 44.5, dNum2DTrackDOCABins, -1.0*dMaxTrackMatchDOCA, dMaxTrackMatchDOCA);
864
865	locHistName = "TOFTrackDeltaYVsVerticalPaddle";
866	locHistTitle = locTrackString + string(";TOF Vertical Paddle;TOF / Track #DeltaY (cm)");
867	dHistMap_TOFPointTrackDeltaYVsVerticalPaddle[locIsTimeBased] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, 44, 0.5, 44.5, dNum2DTrackDOCABins, -1.0*dMaxTrackMatchDOCA, dMaxTrackMatchDOCA);
868
869	locHistName = "TOFTrackDistance_BothPlanes";
870	locHistTitle = locTrackString + string("TOF Hit in Both Planes;TOF / Track Distance (cm)");
871	dHistMap_TOFPointTrackDistance_BothPlanes[locIsTimeBased] = GetOrCreate_Histogram<TH1I>(locHistName, locHistTitle, dNumTrackDOCABins, dMinTrackDOCA, dMaxTrackMatchDOCA);
872
873	locHistName = "TOFTrackDistance_OnePlane";
874	locHistTitle = locTrackString + string("TOF Hit in One Plane;TOF / Track Distance (cm)");
875	dHistMap_TOFPointTrackDistance_OnePlane[locIsTimeBased] = GetOrCreate_Histogram<TH1I>(locHistName, locHistTitle, dNumTrackDOCABins, dMinTrackDOCA, dMaxTrackMatchDOCA);
876	gDirectory(TDirectory::CurrentDirectory())->cd("..");
877
878	//FCAL
879	CreateAndChangeTo_Directory("FCAL", "FCAL");
880	locHistName = "TrackFCALYVsX_HasHit";
881	locHistTitle = locTrackString + string(", Has Other Match, FCAL Has Hit;Projected FCAL Hit X (cm);Projected FCAL Hit Y (cm)");
882	dHistMap_TrackFCALYVsX_HasHit[locIsTimeBased] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNumFCALTOFXYBins, -130.0, 130.0, dNumFCALTOFXYBins, -130.0, 130.0);
883
884	locHistName = "TrackFCALYVsX_NoHit";
885	locHistTitle = locTrackString + string(", Has Other Match, FCAL No Hit;Projected FCAL Hit X (cm);Projected FCAL Hit Y (cm)");
886	dHistMap_TrackFCALYVsX_NoHit[locIsTimeBased] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNumFCALTOFXYBins, -130.0, 130.0, dNumFCALTOFXYBins, -130.0, 130.0);
887
888	locHistName = "TrackFCALRowVsColumn_HasHit";
889	locHistTitle = locTrackString + string(", Has Other Match, FCAL Has Hit;Projected FCAL Hit Column;Projected FCAL Hit Row");
890	dHistMap_TrackFCALRowVsColumn_HasHit[locIsTimeBased] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, 59, -0.5, 58.5, 59, -0.5, 58.5);
891
892	locHistName = "TrackFCALRowVsColumn_NoHit";
893	locHistTitle = locTrackString + string(", Has Other Match, FCAL No Hit;Projected FCAL Hit Column;Projected FCAL Hit Row");
894	dHistMap_TrackFCALRowVsColumn_NoHit[locIsTimeBased] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, 59, -0.5, 58.5, 59, -0.5, 58.5);
895
896	locHistName = "FCALTrackDistanceVsP";
897	locHistTitle = locTrackString + string(";p (GeV/c);FCAL / Track Distance (cm)");
898	dHistMap_FCALTrackDistanceVsP[locIsTimeBased] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNum2DTrackDOCABins, dMinTrackDOCA, dMaxTrackMatchDOCA);
899
900	locHistName = "FCALTrackDistanceVsTheta";
901	locHistTitle = locTrackString + string(";#theta#circ;FCAL / Track Distance (cm)");
902	dHistMap_FCALTrackDistanceVsTheta[locIsTimeBased] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DThetaBins, dMinTheta, 20.0, dNum2DTrackDOCABins, dMinTrackDOCA, dMaxTrackMatchDOCA);
903	gDirectory(TDirectory::CurrentDirectory())->cd("..");
904
905	//BCAL
906	CreateAndChangeTo_Directory("BCAL", "BCAL");
907	locHistName = "TrackBCALModuleVsZ_HasHit";
908	locHistTitle = locTrackString + string(", Has Other Match, BCAL Has Hit;Projected BCAL Hit Z (cm);Projected BCAL Hit Module");
909	dHistMap_TrackBCALModuleVsZ_HasHit[locIsTimeBased] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DBCALZBins, 0.0, 450.0, 48, 0.5, 48.5);
910
911	locHistName = "TrackBCALModuleVsZ_NoHit";
912	locHistTitle = locTrackString + string(", Has Other Match, BCAL No Hit;Projected BCAL Hit Z (cm);Projected BCAL Hit Module");
913	dHistMap_TrackBCALModuleVsZ_NoHit[locIsTimeBased] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DBCALZBins, 0.0, 450.0, 48, 0.5, 48.5);
914
915	locHistName = "TrackBCALPhiVsZ_HasHit";
916	locHistTitle = locTrackString + string(", Has Other Match, BCAL Has Hit;Projected BCAL Hit Z (cm);Projected BCAL Hit #phi#circ");
917	dHistMap_TrackBCALPhiVsZ_HasHit[locIsTimeBased] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DBCALZBins, 0.0, 450.0, dNum2DPhiBins, dMinPhi, dMaxPhi);
918
919	locHistName = "TrackBCALPhiVsZ_NoHit";
920	locHistTitle = locTrackString + string(", Has Other Match, BCAL No Hit;Projected BCAL Hit Z (cm);Projected BCAL Hit #phi#circ");
921	dHistMap_TrackBCALPhiVsZ_NoHit[locIsTimeBased] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DBCALZBins, 0.0, 450.0, dNum2DPhiBins, dMinPhi, dMaxPhi);
922
923	locHistName = "BCALDeltaPhiVsP";
924	locHistTitle = locTrackString + string(";p (GeV/c);BCAL / Track #Delta#phi#circ");
925	dHistMap_BCALDeltaPhiVsP[locIsTimeBased] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, 4.0, dNum2DDeltaPhiBins, dMinDeltaPhi, dMaxDeltaPhi);
926
927	locHistName = "BCALDeltaZVsTheta";
928	locHistTitle = locTrackString + string(";#theta#circ;BCAL / Track #Deltaz (cm)");
929	dHistMap_BCALDeltaZVsTheta[locIsTimeBased] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DThetaBins, dMinTheta, dMaxTheta, dNum2DDeltaZBins, dMinDeltaZ, dMaxDeltaZ);
930	gDirectory(TDirectory::CurrentDirectory())->cd("..");
931
932	//TRACKING
933	locHistName = "PVsTheta_NoHitMatch";
934	locHistTitle = locTrackString + string(", No Hit Match;#theta#circ;p (GeV/c)");
935	dHistMap_TrackPVsTheta_NoHitMatch[locIsTimeBased] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DThetaBins, dMinTheta, dMaxTheta, dNum2DPBins, dMinP, dMaxP);
936
937	locHistName = "PVsTheta_HitMatch";
938	locHistTitle = locTrackString + string(", Hit Match;#theta#circ;p (GeV/c)");
939	dHistMap_TrackPVsTheta_HitMatch[locIsTimeBased] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DThetaBins, dMinTheta, dMaxTheta, dNum2DPBins, dMinP, dMaxP);
940
941	gDirectory(TDirectory::CurrentDirectory())->cd("..");
942	}
943
944	//Return to the base directory
945	ChangeTo_BaseDirectory();
946	}
947	japp->RootUnLock(); //RELEASE ROOT LOCK!!
948	}
949
950	bool DHistogramAction_DetectorMatching::Perform_Action(JEventLoop* locEventLoop, const DParticleCombo* locParticleCombo)
951	{
952	//Expect locParticleCombo to be NULL since this is a reaction-independent action.
953
954	//Optional: Quit the action if it has already been executed this event (else may result in double-counting when filling histograms)
955	if(Get_NumPreviousParticleCombos() != 0)
956	return true;
957
958	bool locIsRESTEvent = (string(locEventLoop->GetJEvent().GetJEventSource()->className()) == string("DEventSourceREST"));
959
960	Fill_MatchingHists(locEventLoop, true); //Time-based tracks
961	if(!locIsRESTEvent)
962	Fill_MatchingHists(locEventLoop, false); //Wire-based tracks
963
964	return true; //return false if you want to use this action to apply a cut (and it fails the cut!)
965	}
966
967	void DHistogramAction_DetectorMatching::Fill_MatchingHists(JEventLoop* locEventLoop, bool locIsTimeBased)
968	{
969	const DParticleID* locParticleID = NULL__null;
970	locEventLoop->GetSingle(locParticleID);
971
972	//can't make this a class member: may cause race condition
973	DCutAction_TrackHitPattern locCutAction_TrackHitPattern(NULL__null, dMinHitRingsPerCDCSuperlayer, dMinHitPlanesPerFDCPackage);
974	locCutAction_TrackHitPattern.Initialize(locEventLoop);
975
976	//get the best tracks for each candidate id, based on good hit pattern & tracking FOM
977	map<JObject::oid_t, const DKinematicData*> locBestTrackMap; //lowest tracking FOM for each candidate id
978	if(locIsTimeBased)
979	{
980	vector<const DTrackTimeBased*> locTrackTimeBasedVector;
981	locEventLoop->Get(locTrackTimeBasedVector);
982
983	//select the best DTrackTimeBased for each track: of tracks with good hit pattern, use best tracking FOM
984	for(size_t loc_i = 0; loc_i < locTrackTimeBasedVector.size(); ++loc_i)
985	{
986	if(locTrackTimeBasedVector[loc_i]->FOM < dMinTrackingFOM)
987	continue;
988	if(!locCutAction_TrackHitPattern.Cut_TrackHitPattern(locTrackTimeBasedVector[loc_i]))
989	continue;
990	JObject::oid_t locCandidateID = locTrackTimeBasedVector[loc_i]->candidateid;
991	if(locBestTrackMap.find(locCandidateID) == locBestTrackMap.end())
992	locBestTrackMap[locCandidateID] = locTrackTimeBasedVector[loc_i];
993	else if(locTrackTimeBasedVector[loc_i]->FOM > (dynamic_cast<const DTrackTimeBased*>(locBestTrackMap[locCandidateID]))->FOM)
994	locBestTrackMap[locCandidateID] = locTrackTimeBasedVector[loc_i];
995	}
996	}
997	else
998	{
999	vector<const DTrackWireBased*> locTrackWireBasedVector;
1000	locEventLoop->Get(locTrackWireBasedVector);
1001
1002	//select the best DTrackWireBased for each track: of tracks with good hit pattern, use best tracking FOM
1003	for(size_t loc_i = 0; loc_i < locTrackWireBasedVector.size(); ++loc_i)
1004	{
1005	if(locTrackWireBasedVector[loc_i]->FOM < dMinTrackingFOM)
1006	continue;
1007	if(!locCutAction_TrackHitPattern.Cut_TrackHitPattern(locTrackWireBasedVector[loc_i]))
1008	continue;
1009	JObject::oid_t locCandidateID = locTrackWireBasedVector[loc_i]->candidateid;
1010	if(locBestTrackMap.find(locCandidateID) == locBestTrackMap.end())
1011	locBestTrackMap[locCandidateID] = locTrackWireBasedVector[loc_i];
1012	else if(locTrackWireBasedVector[loc_i]->FOM > (dynamic_cast<const DTrackWireBased*>(locBestTrackMap[locCandidateID]))->FOM)
1013	locBestTrackMap[locCandidateID] = locTrackWireBasedVector[loc_i];
1014	}
1015	}
1016
1017	vector<const DBCALShower*> locBCALShowers;
1018	locEventLoop->Get(locBCALShowers);
1019
1020	vector<const DFCALShower*> locFCALShowers;
1021	locEventLoop->Get(locFCALShowers);
1022
1023	vector<const DTOFPoint*> locTOFPoints;
1024	locEventLoop->Get(locTOFPoints);
1025
1026	vector<const DTOFPaddleHit*> locTOFPaddleHits;
1027	locEventLoop->Get(locTOFPaddleHits);
1028
1029	vector<const DSCHit*> locSCHits;
1030	locEventLoop->Get(locSCHits);
1031
1032	string locDetectorMatchesTag = locIsTimeBased ? "" : "WireBased";
1033	const DDetectorMatches* locDetectorMatches = NULL__null;
1034	locEventLoop->GetSingle(locDetectorMatches, locDetectorMatchesTag.c_str());
1035
1036	map<JObject::oid_t, const DKinematicData*>::iterator locTrackIterator;
1037
1038	//TRACK / BCAL CLOSEST MATCHES
1039	map<const DKinematicData*, pair<double, double> > locBCALTrackDistanceMap; //first double is delta-phi, second delta-z
1040	for(locTrackIterator = locBestTrackMap.begin(); locTrackIterator != locBestTrackMap.end(); ++locTrackIterator)
1041	{
1042	double locBestMatchDeltaPhi = 0.0, locBestMatchDeltaZ = 0.0;
1043	if(locParticleID->Get_ClosestToTrack_BCAL(locTrackIterator->second, locBCALShowers, locBestMatchDeltaPhi, locBestMatchDeltaZ) != NULL__null)
1044	locBCALTrackDistanceMap[locTrackIterator->second] = pair<double, double>(locBestMatchDeltaPhi, locBestMatchDeltaZ);
1045	}
1046
1047	//TRACK / FCAL CLOSEST MATCHES
1048	map<const DKinematicData*, double> locFCALTrackDistanceMap;
1049	for(locTrackIterator = locBestTrackMap.begin(); locTrackIterator != locBestTrackMap.end(); ++locTrackIterator)
1050	{
1051	double locBestDistance = 999.0;
1052	if(locParticleID->Get_ClosestToTrack_FCAL(locTrackIterator->second, locFCALShowers, locBestDistance) != NULL__null)
1053	locFCALTrackDistanceMap[locTrackIterator->second] = locBestDistance;
1054	}
1055
1056	//TRACK / TOF POINT CLOSEST MATCHES
1057	map<const DKinematicData, pair<const DTOFPoint, pair<double, double> > > locTOFPointTrackDistanceMap; //doubles: delta-x, delta-y
1058	for(locTrackIterator = locBestTrackMap.begin(); locTrackIterator != locBestTrackMap.end(); ++locTrackIterator)
1059	{
1060	double locBestDeltaX, locBestDeltaY;
1061	const DTOFPoint* locClosestTOFPoint = locParticleID->Get_ClosestToTrack_TOFPoint(locTrackIterator->second, locTOFPoints, locBestDeltaX, locBestDeltaY);
1062	if(locClosestTOFPoint == NULL__null)
1063	continue;
1064	pair<double, double> locDeltas(locBestDeltaX, locBestDeltaY);
1065	locTOFPointTrackDistanceMap[locTrackIterator->second] = pair<const DTOFPoint*, pair<double, double> >(locClosestTOFPoint, locDeltas);
1066	}
1067
1068	//TRACK / TOF PADDLE CLOSEST MATCHES
1069	map<const DKinematicData, pair<const DTOFPaddleHit, double> > locHorizontalTOFPaddleTrackDistanceMap; //double: delta-y
1070	map<const DKinematicData, pair<const DTOFPaddleHit, double> > locVerticalTOFPaddleTrackDistanceMap; //double: delta-x
1071	for(locTrackIterator = locBestTrackMap.begin(); locTrackIterator != locBestTrackMap.end(); ++locTrackIterator)
1072	{
1073	double locBestDeltaX = 999.9, locBestDeltaY = 999.9;
1074	pair<const DTOFPaddleHit, const DTOFPaddleHit> locClosestTOFPaddleHits = locParticleID->Get_ClosestToTrack_TOFPaddles(locTrackIterator->second, locTOFPaddleHits, locBestDeltaX, locBestDeltaY);
1075	if(locClosestTOFPaddleHits.first != NULL__null)
1076	locVerticalTOFPaddleTrackDistanceMap[locTrackIterator->second] = pair<const DTOFPaddleHit*, double>(locClosestTOFPaddleHits.first, locBestDeltaX);
1077	if(locClosestTOFPaddleHits.second != NULL__null)
1078	locHorizontalTOFPaddleTrackDistanceMap[locTrackIterator->second] = pair<const DTOFPaddleHit*, double>(locClosestTOFPaddleHits.first, locBestDeltaY);
1079	}
1080
1081	//TRACK / SC CLOSEST MATCHES
1082	map<const DKinematicData*, double> locSCTrackDistanceMap;
1083	for(locTrackIterator = locBestTrackMap.begin(); locTrackIterator != locBestTrackMap.end(); ++locTrackIterator)
1084	{
1085	double locBestDeltaPhi = 999.0;
1086	if(locParticleID->Get_ClosestToTrack_SC(locTrackIterator->second, locSCHits, locBestDeltaPhi) != NULL__null)
1087	locSCTrackDistanceMap[locTrackIterator->second] = locBestDeltaPhi;
1088	}
1089
1090	//PROJECTED HIT POSITIONS
1091	map<const DKinematicData*, pair<int, bool> > dProjectedSCPaddleMap; //pair: paddle, hit-barrel-flag (false if bend/nose)
1092	map<const DKinematicData*, pair<int, int> > dProjectedTOF2DPaddlesMap; //pair: vertical, horizontal
1093	map<const DKinematicData*, pair<float, float> > dProjectedTOFXYMap; //pair: x, y
1094	map<const DKinematicData*, pair<int, int> > dProjectedFCALRowColumnMap; //pair: column, row
1095	map<const DKinematicData*, pair<float, float> > dProjectedFCALXYMap; //pair: x, y
1096	map<const DKinematicData*, pair<float, int> > dProjectedBCALModuleSectorMap; //pair: z, module
1097	map<const DKinematicData*, pair<float, float> > dProjectedBCALPhiZMap; //pair: z, phi
1098	for(locTrackIterator = locBestTrackMap.begin(); locTrackIterator != locBestTrackMap.end(); ++locTrackIterator)
1099	{
1100	const DKinematicData* locTrack = locTrackIterator->second;
1101	const DReferenceTrajectory* locReferenceTrajectory = NULL__null;
1102	if(locIsTimeBased)
1103	locReferenceTrajectory = (static_cast<const DTrackTimeBased*>(locTrack))->rt;
1104	else
1105	locReferenceTrajectory = (static_cast<const DTrackWireBased*>(locTrack))->rt;
1106	if(locReferenceTrajectory == NULL__null)
1107	break; //e.g. REST data: no trajectory
1108
1109	//SC
1110	DVector3 locSCIntersection;
1111	bool locProjBarrelFlag = false;
1112	unsigned int locProjectedSCPaddle = locParticleID->PredictSCSector(locReferenceTrajectory, 999.9, &locSCIntersection, &locProjBarrelFlag);
1113	if(locProjectedSCPaddle != 0)
1114	dProjectedSCPaddleMap[locTrack] = pair<int, bool>(locProjectedSCPaddle, locProjBarrelFlag);
1115
1116	//TOF
1117	DVector3 locTOFIntersection;
1118	unsigned int locHorizontalBar = 0, locVerticalBar = 0;
1119	if(locParticleID->PredictTOFPaddles(locReferenceTrajectory, locHorizontalBar, locVerticalBar, &locTOFIntersection))
1120	{
1121	dProjectedTOF2DPaddlesMap[locTrack] = pair<int, int>(locVerticalBar, locHorizontalBar);
1122	dProjectedTOFXYMap[locTrack] = pair<float, float>(locTOFIntersection.X(), locTOFIntersection.Y());
1123	}
1124
1125	//FCAL
1126	DVector3 locFCALIntersection;
1127	unsigned int locRow = 0, locColumn = 0;
1128	if(locParticleID->PredictFCALHit(locReferenceTrajectory, locRow, locColumn, &locFCALIntersection))
1129	{
1130	dProjectedFCALRowColumnMap[locTrack] = pair<int, int>(locColumn, locRow);
1131	dProjectedFCALXYMap[locTrack] = pair<float, float>(locFCALIntersection.X(), locFCALIntersection.Y());
1132	}
1133
1134	//BCAL
1135	DVector3 locBCALIntersection;
1136	unsigned int locModule = 0, locSector = 0;
1137	if(locParticleID->PredictBCALWedge(locReferenceTrajectory, locModule, locSector, &locBCALIntersection))
1138	{
1139	dProjectedBCALModuleSectorMap[locTrack] = pair<float, int>(locBCALIntersection.Z(), locModule);
1140	dProjectedBCALPhiZMap[locTrack] = pair<float, float>(locBCALIntersection.Z(), locBCALIntersection.Phi()*180.0/TMath::Pi());
1141	}
1142	}
1143
1144	//Fill Histograms
1145	japp->RootWriteLock(); //ACQUIRE ROOT LOCK!!
1146	{
1147
1148	/******************************************************** MATCHING DISTANCE ********************************************************/
1149
1150	//BCAL
1151	map<const DKinematicData*, pair<double, double> >::iterator locBCALIterator = locBCALTrackDistanceMap.begin();
1152	for(; locBCALIterator != locBCALTrackDistanceMap.end(); ++locBCALIterator)
1153	{
1154	const DKinematicData* locTrack = locBCALIterator->first;
1155	dHistMap_BCALDeltaPhiVsP[locIsTimeBased]->Fill(locTrack->momentum().Mag(), locBCALIterator->second.first*180.0/TMath::Pi());
1156	dHistMap_BCALDeltaZVsTheta[locIsTimeBased]->Fill(locTrack->momentum().Theta()*180.0/TMath::Pi(), locBCALIterator->second.second);
1157	}
1158
1159	//FCAL
1160	map<const DKinematicData*, double>::iterator locFCALIterator = locFCALTrackDistanceMap.begin();
1161	for(; locFCALIterator != locFCALTrackDistanceMap.end(); ++locFCALIterator)
1162	{
1163	const DKinematicData* locTrack = locFCALIterator->first;
1164	dHistMap_FCALTrackDistanceVsP[locIsTimeBased]->Fill(locTrack->momentum().Mag(), locFCALIterator->second);
1165	dHistMap_FCALTrackDistanceVsTheta[locIsTimeBased]->Fill(locTrack->momentum().Theta()*180.0/TMath::Pi(), locFCALIterator->second);
1166	}
1167
1168	//TOF Paddle
1169	//Horizontal
1170	map<const DKinematicData, pair<const DTOFPaddleHit, double> >::iterator locTOFPaddleIterator = locHorizontalTOFPaddleTrackDistanceMap.begin();
1171	for(; locTOFPaddleIterator != locHorizontalTOFPaddleTrackDistanceMap.end(); ++locTOFPaddleIterator)
1172	{
1173	double locDistance = locTOFPaddleIterator->second.second;
1174	dHistMap_TOFPaddleTrackDeltaY[locIsTimeBased]->Fill(locDistance);
1175	}
1176	//Vertical
1177	locTOFPaddleIterator = locVerticalTOFPaddleTrackDistanceMap.begin();
1178	for(; locTOFPaddleIterator != locVerticalTOFPaddleTrackDistanceMap.end(); ++locTOFPaddleIterator)
1179	{
1180	double locDistance = locTOFPaddleIterator->second.second;
1181	dHistMap_TOFPaddleTrackDeltaX[locIsTimeBased]->Fill(locDistance);
1182	}
1183
1184	//TOF Point
1185	map<const DKinematicData, pair<const DTOFPoint, pair<double, double> > >::iterator locTOFPointIterator = locTOFPointTrackDistanceMap.begin();
1186	for(; locTOFPointIterator != locTOFPointTrackDistanceMap.end(); ++locTOFPointIterator)
1187	{
1188	const DKinematicData* locTrack = locTOFPointIterator->first;
1189	const DTOFPoint* locTOFPoint = locTOFPointIterator->second.first;
1190	double locDeltaX = locTOFPointIterator->second.second.first;
1191	double locDeltaY = locTOFPointIterator->second.second.second;
1192
1193	double locDistance = sqrt(locDeltaXlocDeltaX + locDeltaYlocDeltaY);
1194	if((locDeltaX < 500.0) && (locDeltaY < 500.0)) //else position not well-defined
1195	{
1196	dHistMap_TOFPointTrackDistanceVsP[locIsTimeBased]->Fill(locTrack->momentum().Mag(), locDistance);
1197	dHistMap_TOFPointTrackDistanceVsTheta[locIsTimeBased]->Fill(locTrack->momentum().Theta()*180.0/TMath::Pi(), locDistance);
1198	if((locTOFPoint->dHorizontalBar != 0) && (locTOFPoint->dVerticalBar != 0))
1199	dHistMap_TOFPointTrackDistance_BothPlanes[locIsTimeBased]->Fill(locDistance);
1200	else
1201	dHistMap_TOFPointTrackDistance_OnePlane[locIsTimeBased]->Fill(locDistance);
1202	}
1203
1204	dHistMap_TOFPointTrackDeltaXVsHorizontalPaddle[locIsTimeBased]->Fill(locTOFPoint->dHorizontalBar, locDeltaX);
1205	dHistMap_TOFPointTrackDeltaXVsVerticalPaddle[locIsTimeBased]->Fill(locTOFPoint->dVerticalBar, locDeltaX);
1206
1207	dHistMap_TOFPointTrackDeltaYVsHorizontalPaddle[locIsTimeBased]->Fill(locTOFPoint->dHorizontalBar, locDeltaY);
1208	dHistMap_TOFPointTrackDeltaYVsVerticalPaddle[locIsTimeBased]->Fill(locTOFPoint->dVerticalBar, locDeltaY);
1209	}
1210
1211	//SC
1212	map<const DKinematicData*, double>::iterator locSCIterator = locSCTrackDistanceMap.begin();
1213	if(locSCHits.size() <= 4) //don't fill if every paddle fired!
1214	{
1215	for(; locSCIterator != locSCTrackDistanceMap.end(); ++locSCIterator)
1216	{
1217	const DKinematicData* locTrack = locSCIterator->first;
1218	double locDeltaPhi = locSCIterator->second*180.0/TMath::Pi();
1219	dHistMap_SCTrackDeltaPhiVsP[locIsTimeBased]->Fill(locTrack->momentum().Mag(), locDeltaPhi);
1220	dHistMap_SCTrackDeltaPhiVsTheta[locIsTimeBased]->Fill(locTrack->momentum().Theta()*180.0/TMath::Pi(), locDeltaPhi);
1221	}
1222	}
1223
1224	/******************************************************* MATCHING EFFICINECY *******************************************************/
1225
1226	//Does-it-match, by detector
1227	for(locTrackIterator = locBestTrackMap.begin(); locTrackIterator != locBestTrackMap.end(); ++locTrackIterator)
1228	{
1229	const DKinematicData* locTrack = locTrackIterator->second;
1230	double locTheta = locTrack->momentum().Theta()*180.0/TMath::Pi();
1231	double locPhi = locTrack->momentum().Phi()*180.0/TMath::Pi();
1232	double locP = locTrack->momentum().Mag();
1233
1234	//BCAL
1235	if(locDetectorMatches->Get_IsMatchedToDetector(locTrack, SYS_START))
1236	{
1237	if(locDetectorMatches->Get_IsMatchedToDetector(locTrack, SYS_BCAL))
1238	{
1239	dHistMap_PVsTheta_HasHit[SYS_BCAL][locIsTimeBased]->Fill(locTheta, locP);
1240	dHistMap_PhiVsTheta_HasHit[SYS_BCAL][locIsTimeBased]->Fill(locTheta, locPhi);
1241	if(dProjectedBCALModuleSectorMap.find(locTrack) != dProjectedBCALModuleSectorMap.end())
1242	{
1243	pair<float, float>& locPositionPair = dProjectedBCALPhiZMap[locTrack];
1244	dHistMap_TrackBCALPhiVsZ_HasHit[locIsTimeBased]->Fill(locPositionPair.first, locPositionPair.second);
1245	pair<float, int>& locElementPair = dProjectedBCALModuleSectorMap[locTrack];
1246	dHistMap_TrackBCALModuleVsZ_HasHit[locIsTimeBased]->Fill(locElementPair.first, locElementPair.second);
1247	}
1248	}
1249	else
1250	{
1251	dHistMap_PVsTheta_NoHit[SYS_BCAL][locIsTimeBased]->Fill(locTheta, locP);
1252	dHistMap_PhiVsTheta_NoHit[SYS_BCAL][locIsTimeBased]->Fill(locTheta, locPhi);
1253	if(dProjectedBCALModuleSectorMap.find(locTrack) != dProjectedBCALModuleSectorMap.end())
1254	{
1255	pair<float, float>& locPositionPair = dProjectedBCALPhiZMap[locTrack];
1256	dHistMap_TrackBCALPhiVsZ_NoHit[locIsTimeBased]->Fill(locPositionPair.first, locPositionPair.second);
1257	pair<float, int>& locElementPair = dProjectedBCALModuleSectorMap[locTrack];
1258	dHistMap_TrackBCALModuleVsZ_NoHit[locIsTimeBased]->Fill(locElementPair.first, locElementPair.second);
1259	}
1260	}
1261	}
1262
1263	//FCAL
1264	if(locDetectorMatches->Get_IsMatchedToDetector(locTrack, SYS_TOF))
1265	{
1266	if(locDetectorMatches->Get_IsMatchedToDetector(locTrack, SYS_FCAL))
1267	{
1268	dHistMap_PVsTheta_HasHit[SYS_FCAL][locIsTimeBased]->Fill(locTheta, locP);
1269	dHistMap_PhiVsTheta_HasHit[SYS_FCAL][locIsTimeBased]->Fill(locTheta, locPhi);
1270	if(dProjectedFCALRowColumnMap.find(locTrack) != dProjectedFCALRowColumnMap.end())
1271	{
1272	pair<float, float>& locPositionPair = dProjectedFCALXYMap[locTrack];
1273	dHistMap_TrackFCALYVsX_HasHit[locIsTimeBased]->Fill(locPositionPair.first, locPositionPair.second);
1274	pair<int, int>& locElementPair = dProjectedFCALRowColumnMap[locTrack];
1275	dHistMap_TrackFCALRowVsColumn_HasHit[locIsTimeBased]->Fill(locElementPair.first, locElementPair.second);
1276	}
1277	}
1278	else
1279	{
1280	dHistMap_PVsTheta_NoHit[SYS_FCAL][locIsTimeBased]->Fill(locTheta, locP);
1281	dHistMap_PhiVsTheta_NoHit[SYS_FCAL][locIsTimeBased]->Fill(locTheta, locPhi);
1282	if(dProjectedFCALRowColumnMap.find(locTrack) != dProjectedFCALRowColumnMap.end())
1283	{
1284	pair<float, float>& locPositionPair = dProjectedFCALXYMap[locTrack];
1285	dHistMap_TrackFCALYVsX_NoHit[locIsTimeBased]->Fill(locPositionPair.first, locPositionPair.second);
1286	pair<int, int>& locElementPair = dProjectedFCALRowColumnMap[locTrack];
1287	dHistMap_TrackFCALRowVsColumn_NoHit[locIsTimeBased]->Fill(locElementPair.first, locElementPair.second);
1288	}
1289	}
1290	}
1291
1292	//TOF Paddle
1293	if((dProjectedTOFXYMap.find(locTrack) != dProjectedTOFXYMap.end()) && locDetectorMatches->Get_IsMatchedToDetector(locTrack, SYS_FCAL))
1294	{
1295	pair<float, float>& locPositionPair = dProjectedTOFXYMap[locTrack];
1296	pair<int, int>& locPaddlePair = dProjectedTOF2DPaddlesMap[locTrack]; //vertical, horizontal
1297
1298	//Horizontal
1299	if(locHorizontalTOFPaddleTrackDistanceMap.find(locTrack) != locHorizontalTOFPaddleTrackDistanceMap.end())
1300	{
1301	double locDistance = locHorizontalTOFPaddleTrackDistanceMap[locTrack].second;
1302	if(locDistance <= dMinTOFPaddleMatchDistance) //match
1303	dHistMap_TOFPaddleHorizontalPaddleVsTrackX_HasHit[locIsTimeBased]->Fill(locPositionPair.first, locPaddlePair.second);
1304	else //no match
1305	dHistMap_TOFPaddleHorizontalPaddleVsTrackX_NoHit[locIsTimeBased]->Fill(locPositionPair.first, locPaddlePair.second);
1306	}
1307	else // no match
1308	dHistMap_TOFPaddleHorizontalPaddleVsTrackX_NoHit[locIsTimeBased]->Fill(locPositionPair.first, locPaddlePair.second);
1309
1310	//Vertical
1311	if(locVerticalTOFPaddleTrackDistanceMap.find(locTrack) != locVerticalTOFPaddleTrackDistanceMap.end())
1312	{
1313	double locDistance = locVerticalTOFPaddleTrackDistanceMap[locTrack].second;
1314	if(locDistance <= dMinTOFPaddleMatchDistance) //match
1315	dHistMap_TOFPaddleTrackYVsVerticalPaddle_HasHit[locIsTimeBased]->Fill(locPaddlePair.first, locPositionPair.second);
1316	else //no match
1317	dHistMap_TOFPaddleTrackYVsVerticalPaddle_NoHit[locIsTimeBased]->Fill(locPaddlePair.first, locPositionPair.second);
1318	}
1319	else // no match
1320	dHistMap_TOFPaddleTrackYVsVerticalPaddle_NoHit[locIsTimeBased]->Fill(locPaddlePair.first, locPositionPair.second);
1321	}
1322
1323	//TOF Point
1324	if(locDetectorMatches->Get_IsMatchedToDetector(locTrack, SYS_FCAL))
1325	{
1326	if(locDetectorMatches->Get_IsMatchedToDetector(locTrack, SYS_TOF))
1327	{
1328	dHistMap_PVsTheta_HasHit[SYS_TOF][locIsTimeBased]->Fill(locTheta, locP);
1329	dHistMap_PhiVsTheta_HasHit[SYS_TOF][locIsTimeBased]->Fill(locTheta, locPhi);
1330	if(dProjectedTOFXYMap.find(locTrack) != dProjectedTOFXYMap.end())
1331	{
1332	pair<float, float>& locPositionPair = dProjectedTOFXYMap[locTrack];
1333	dHistMap_TrackTOFYVsX_HasHit[locIsTimeBased]->Fill(locPositionPair.first, locPositionPair.second);
1334	pair<int, int>& locElementPair = dProjectedTOF2DPaddlesMap[locTrack];
1335	dHistMap_TrackTOF2DPaddles_HasHit[locIsTimeBased]->Fill(locElementPair.first, locElementPair.second);
1336	}
1337	}
1338	else
1339	{
1340	dHistMap_PVsTheta_NoHit[SYS_TOF][locIsTimeBased]->Fill(locTheta, locP);
1341	dHistMap_PhiVsTheta_NoHit[SYS_TOF][locIsTimeBased]->Fill(locTheta, locPhi);
1342	if(dProjectedTOFXYMap.find(locTrack) != dProjectedTOFXYMap.end())
1343	{
1344	pair<float, float>& locPositionPair = dProjectedTOFXYMap[locTrack];
1345	dHistMap_TrackTOFYVsX_NoHit[locIsTimeBased]->Fill(locPositionPair.first, locPositionPair.second);
1346	pair<int, int>& locElementPair = dProjectedTOF2DPaddlesMap[locTrack];
1347	dHistMap_TrackTOF2DPaddles_NoHit[locIsTimeBased]->Fill(locElementPair.first, locElementPair.second);
1348	}
1349	}
1350	}
1351
1352	//SC
1353	if(locDetectorMatches->Get_IsMatchedToDetector(locTrack, SYS_BCAL) \|\| locDetectorMatches->Get_IsMatchedToDetector(locTrack, SYS_FCAL) \|\| locDetectorMatches->Get_IsMatchedToDetector(locTrack, SYS_TOF))
1354	{
1355	if(locDetectorMatches->Get_IsMatchedToDetector(locTrack, SYS_START))
1356	{
1357	dHistMap_PVsTheta_HasHit[SYS_START][locIsTimeBased]->Fill(locTheta, locP);
1358	dHistMap_PhiVsTheta_HasHit[SYS_START][locIsTimeBased]->Fill(locTheta, locPhi);
1359	if(dProjectedSCPaddleMap.find(locTrack) != dProjectedSCPaddleMap.end())
1360	{
1361	dHistMap_SCPaddleVsTheta_HasHit[locIsTimeBased]->Fill(locTheta, dProjectedSCPaddleMap[locTrack].first);
1362	if(dProjectedSCPaddleMap[locTrack].second)
1363	dHistMap_SCPaddle_BarrelRegion_HasHit[locIsTimeBased]->Fill(dProjectedSCPaddleMap[locTrack].first);
1364	else
1365	dHistMap_SCPaddle_NoseRegion_HasHit[locIsTimeBased]->Fill(dProjectedSCPaddleMap[locTrack].first);
1366	}
1367	}
1368	else
1369	{
1370	dHistMap_PVsTheta_NoHit[SYS_START][locIsTimeBased]->Fill(locTheta, locP);
1371	dHistMap_PhiVsTheta_NoHit[SYS_START][locIsTimeBased]->Fill(locTheta, locPhi);
1372	if(dProjectedSCPaddleMap.find(locTrack) != dProjectedSCPaddleMap.end())
1373	{
1374	dHistMap_SCPaddleVsTheta_NoHit[locIsTimeBased]->Fill(locTheta, dProjectedSCPaddleMap[locTrack].first);
1375	if(dProjectedSCPaddleMap[locTrack].second)
1376	dHistMap_SCPaddle_BarrelRegion_NoHit[locIsTimeBased]->Fill(dProjectedSCPaddleMap[locTrack].first);
1377	else
1378	dHistMap_SCPaddle_NoseRegion_NoHit[locIsTimeBased]->Fill(dProjectedSCPaddleMap[locTrack].first);
1379	}
1380	}
1381	}
1382	}
1383	//Is-Matched to Something
1384	for(locTrackIterator = locBestTrackMap.begin(); locTrackIterator != locBestTrackMap.end(); ++locTrackIterator)
1385	{
1386	const DKinematicData* locTrack = locTrackIterator->second;
1387	double locTheta = locTrack->momentum().Theta()*180.0/TMath::Pi();
1388	double locP = locTrack->momentum().Mag();
1389	if(locDetectorMatches->Get_IsMatchedToHit(locTrack))
1390	dHistMap_TrackPVsTheta_HitMatch[locIsTimeBased]->Fill(locTheta, locP);
1391	else
1392	dHistMap_TrackPVsTheta_NoHitMatch[locIsTimeBased]->Fill(locTheta, locP);
1393	}
1394	}
1395	japp->RootUnLock(); //RELEASE ROOT LOCK!!
1396	}
1397
1398	void DHistogramAction_DetectorPID::Initialize(JEventLoop* locEventLoop)
1399	{
1400	//Create any histograms/trees/etc. within a ROOT lock.
1401	//This is so that when running multithreaded, only one thread is writing to the ROOT file at a time.
1402
1403	//When creating a reaction-independent action, only modify member variables within a ROOT lock.
1404	//Objects created within a plugin (such as reaction-independent actions) can be accessed by many threads simultaneously.
1405
1406	string locHistName, locHistTitle, locParticleROOTName;
1407
1408	string locTrackSelectionTag = "NotATag", locShowerSelectionTag = "NotATag";
1409	if(gPARMS->Exists("COMBO:TRACK_SELECT_TAG"))
1410	gPARMS->GetParameter("COMBO:TRACK_SELECT_TAG", locTrackSelectionTag);
1411	if(gPARMS->Exists("COMBO:SHOWER_SELECT_TAG"))
1412	gPARMS->GetParameter("COMBO:SHOWER_SELECT_TAG", locShowerSelectionTag);
1413
1414	japp->RootWriteLock(); //ACQUIRE ROOT LOCK!!
1415	{
1416	//So: Default tag is "", User can set it to something else
1417	//In here, if tag is "", get from gparms, if not, leave it alone
1418	//If gparms value does not exist, set it to (and use) "PreSelect"
1419	if(dTrackSelectionTag == "NotATag")
1420	dTrackSelectionTag = (locTrackSelectionTag == "NotATag") ? "PreSelect" : locTrackSelectionTag;
1421	if(dShowerSelectionTag == "NotATag")
1422	dShowerSelectionTag = (locShowerSelectionTag == "NotATag") ? "PreSelect" : locShowerSelectionTag;
1423
1424	//Required: Create a folder in the ROOT output file that will contain all of the output ROOT objects (if any) for this action.
1425	//If another thread has already created the folder, it just changes to it.
1426	CreateAndChangeTo_ActionDirectory();
1427
1428	//q = 0
1429	locParticleROOTName = ParticleName_ROOT(Gamma);
1430	//BCAL
1431	CreateAndChangeTo_Directory("BCAL", "BCAL");
1432
1433	locHistName = "BetaVsP_q0";
1434	locHistTitle = "BCAL q^{0};Shower Energy (GeV);#beta";
1435	dHistMap_BetaVsP[SYS_BCAL][0] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxBCALP, dNum2DBetaBins, dMinBeta, dMaxBeta);
1436
1437	locHistName = "DeltaTVsShowerE_Photon";
1438	locHistTitle = string("BCAL q^{0}") + locParticleROOTName + string(";Shower Energy (GeV);#Deltat_{BCAL - RF}");
1439	dHistMap_DeltaTVsP[SYS_BCAL][Gamma] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNum2DDeltaTBins, dMinDeltaT, dMaxDeltaT);
1440
1441	/*
1442	//Uncomment when ready!
1443	locHistName = "TimePullVsShowerE_Photon";
1444	locHistTitle = string("BCAL ") + locParticleROOTName + string(";Shower Energy (GeV);#Deltat/#sigma_{#Deltat}");
1445	dHistMap_TimePullVsP[SYS_BCAL][Gamma] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNum2DPullBins, dMinPull, dMaxPull);
1446
1447	locHistName = "TimeFOMVsShowerE_Photon";
1448	locHistTitle = string("BCAL ") + locParticleROOTName + string(";Shower Energy (GeV);Timing PID Confidence Level");
1449	dHistMap_TimeFOMVsP[SYS_BCAL][Gamma] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNum2DFOMBins, 0.0, 1.0);
1450	*/
1451	gDirectory(TDirectory::CurrentDirectory())->cd("..");
1452
1453	//FCAL
1454	CreateAndChangeTo_Directory("FCAL", "FCAL");
1455
1456	locHistName = "BetaVsP_q0";
1457	locHistTitle = "FCAL q^{0};Shower Energy (GeV);#beta";
1458	dHistMap_BetaVsP[SYS_FCAL][0] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNum2DBetaBins, dMinBeta, dMaxBeta);
1459
1460	locHistName = "DeltaTVsShowerE_Photon";
1461	locHistTitle = string("FCAL ") + locParticleROOTName + string(";Shower Energy (GeV);#Deltat_{FCAL - RF}");
1462	dHistMap_DeltaTVsP[SYS_FCAL][Gamma] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNum2DDeltaTBins, dMinDeltaT, dMaxDeltaT);
1463
1464	/*
1465	//Uncomment when ready!
1466	locHistName = "TimePullVsShowerE_Photon";
1467	locHistTitle = string("FCAL ") + locParticleROOTName + string(";Shower Energy (GeV);#Deltat/#sigma_{#Deltat}");
1468	dHistMap_TimePullVsP[SYS_FCAL][Gamma] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNum2DPullBins, dMinPull, dMaxPull);
1469
1470	locHistName = "TimeFOMVsShowerE_Photon";
1471	locHistTitle = string("FCAL ") + locParticleROOTName + string(";Shower Energy (GeV);Timing PID Confidence Level");
1472	dHistMap_TimeFOMVsP[SYS_FCAL][Gamma] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNum2DFOMBins, 0.0, 1.0);
1473	*/
1474	gDirectory(TDirectory::CurrentDirectory())->cd("..");
1475
1476	//q +/-
1477	for(int locCharge = -1; locCharge <= 1; locCharge += 2)
1478	{
1479	string locParticleName = (locCharge == -1) ? "q-" : "q+";
1480	string locParticleROOTName = (locCharge == -1) ? "q^{-}" : "q^{+}";
1481
1482	//SC
1483	CreateAndChangeTo_Directory("SC", "SC");
1484
1485	locHistName = string("dEdXVsP_") + locParticleName;
1486	locHistTitle = locParticleROOTName + ";p (GeV/c);SC dE/dX (MeV/cm)";
1487	dHistMap_dEdXVsP[SYS_START][locCharge] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNum2DdEdxBins, dMindEdX, dMaxdEdX);
1488
1489	locHistName = string("BetaVsP_") + locParticleName;
1490	locHistTitle = string("SC ") + locParticleROOTName + string(";p (GeV/c);#beta");
1491	dHistMap_BetaVsP[SYS_START][locCharge] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNum2DBetaBins, dMinBeta, dMaxBeta);
1492
1493	gDirectory(TDirectory::CurrentDirectory())->cd("..");
1494
1495	//TOF
1496	CreateAndChangeTo_Directory("TOF", "TOF");
1497
1498	locHistName = string("dEdXVsP_") + locParticleName;
1499	locHistTitle = locParticleROOTName + ";p (GeV/c);TOF dE/dX (MeV/cm)";
1500	dHistMap_dEdXVsP[SYS_TOF][locCharge] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNum2DdEdxBins, dMindEdX, dMaxdEdX);
1501
1502	locHistName = string("BetaVsP_") + locParticleName;
1503	locHistTitle = string("TOF ") + locParticleROOTName + string(";p (GeV/c);#beta");
1504	dHistMap_BetaVsP[SYS_TOF][locCharge] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNum2DBetaBins, dMinBeta, dMaxBeta);
1505
1506	gDirectory(TDirectory::CurrentDirectory())->cd("..");
1507
1508	//BCAL
1509	CreateAndChangeTo_Directory("BCAL", "BCAL");
1510
1511	locHistName = string("BetaVsP_") + locParticleName;
1512	locHistTitle = string("BCAL ") + locParticleROOTName + string(";p (GeV/c);#beta");
1513	dHistMap_BetaVsP[SYS_BCAL][locCharge] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxBCALP, dNum2DBetaBins, dMinBeta, dMaxBeta);
1514
1515	locHistName = string("EOverPVsP_") + locParticleName;
1516	locHistTitle = string("BCAL ") + locParticleROOTName + string(";p (GeV/c);E_{Shower}/p_{Track} (c);");
1517	dHistMap_BCALEOverPVsP[locCharge] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxBCALP, dNum2DEOverPBins, dMinEOverP, dMaxEOverP);
1518
1519	locHistName = string("EOverPVsTheta_") + locParticleName;
1520	locHistTitle = string("BCAL ") + locParticleROOTName + string(";#theta#circ;E_{Shower}/p_{Track} (c);");
1521	dHistMap_BCALEOverPVsTheta[locCharge] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DBCALThetaBins, dMinBCALTheta, dMaxBCALTheta, dNum2DEOverPBins, dMinEOverP, dMaxEOverP);
1522
1523	gDirectory(TDirectory::CurrentDirectory())->cd("..");
1524
1525	//FCAL
1526	CreateAndChangeTo_Directory("FCAL", "FCAL");
1527
1528	locHistName = string("BetaVsP_") + locParticleName;
1529	locHistTitle = string("FCAL ") + locParticleROOTName + string(";p (GeV/c);#beta");
1530	dHistMap_BetaVsP[SYS_FCAL][locCharge] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNum2DBetaBins, dMinBeta, dMaxBeta);
1531
1532	locHistName = string("EOverPVsP_") + locParticleName;
1533	locHistTitle = string("FCAL ") + locParticleROOTName + string(";p (GeV/c);E_{Shower}/p_{Track} (c);");
1534	dHistMap_FCALEOverPVsP[locCharge] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNum2DEOverPBins, dMinEOverP, dMaxEOverP);
1535
1536	locHistName = string("EOverPVsTheta_") + locParticleName;
1537	locHistTitle = string("FCAL ") + locParticleROOTName + string(";#theta#circ;E_{Shower}/p_{Track} (c);");
1538	dHistMap_FCALEOverPVsTheta[locCharge] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DFCALThetaBins, dMinFCALTheta, dMaxFCALTheta, dNum2DEOverPBins, dMinEOverP, dMaxEOverP);
1539
1540	gDirectory(TDirectory::CurrentDirectory())->cd("..");
1541
1542	//CDC
1543	CreateAndChangeTo_Directory("CDC", "CDC");
1544
1545	locHistName = string("dEdXVsP_") + locParticleName;
1546	locHistTitle = locParticleROOTName + string(";p (GeV/c);CDC dE/dx (keV/cm)");
1547	dHistMap_dEdXVsP[SYS_CDC][locCharge] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNum2DdEdxBins, dMindEdX, dMaxdEdX);
1548
1549	gDirectory(TDirectory::CurrentDirectory())->cd("..");
1550
1551	//FDC
1552	CreateAndChangeTo_Directory("FDC", "FDC");
1553
1554	locHistName = string("dEdXVsP_") + locParticleName;
1555	locHistTitle = locParticleROOTName + string(";p (GeV/c);FDC dE/dx (keV/cm)");
1556	dHistMap_dEdXVsP[SYS_FDC][locCharge] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNum2DdEdxBins, dMindEdX, dMaxdEdX);
1557
1558	gDirectory(TDirectory::CurrentDirectory())->cd("..");
1559	}
1560
1561	//delta's by PID
1562	for(size_t loc_i = 0; loc_i < dFinalStatePIDs.size(); ++loc_i)
1563	{
1564	Particle_t locPID = dFinalStatePIDs[loc_i];
1565	string locParticleName = ParticleType(locPID);
1566	string locParticleROOTName = ParticleName_ROOT(locPID);
1567
1568	//SC
1569	CreateAndChangeTo_Directory("SC", "SC");
1570
1571	locHistName = string("DeltadEdXVsP_") + locParticleName;
1572	locHistTitle = locParticleROOTName + " Candidates;p (GeV/c);SC #Delta(dE/dX) (MeV/cm)";
1573	dHistMap_DeltadEdXVsP[SYS_START][locPID] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNum2DDeltadEdxBins, dMinDeltadEdx, dMaxDeltadEdx);
1574
1575	locHistName = string("DeltaBetaVsP_") + locParticleName;
1576	locHistTitle = string("SC ") + locParticleROOTName + string(" Candidates;p (GeV/c);#Delta#beta");
1577	dHistMap_DeltaBetaVsP[SYS_START][locPID] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNum2DDeltaBetaBins, dMinDeltaBeta, dMaxDeltaBeta);
1578
1579	locHistName = string("DeltaTVsP_") + locParticleName;
1580	locHistTitle = string("SC ") + locParticleROOTName + string(" Candidates;p (GeV/c);#Deltat_{SC - RF}");
1581	dHistMap_DeltaTVsP[SYS_START][locPID] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNum2DDeltaTBins, dMinDeltaT, dMaxDeltaT);
1582
1583	/*
1584	//Uncomment when ready!
1585	locHistName = string("TimePullVsP_") + locParticleName;
1586	locHistTitle = string("SC ") + locParticleROOTName + string(" Candidates;p (GeV/c);#Deltat/#sigma_{#Deltat}");
1587	dHistMap_TimePullVsP[SYS_START][locPID] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNum2DPullBins, dMinPull, dMaxPull);
1588
1589	locHistName = string("TimeFOMVsP_") + locParticleName;
1590	locHistTitle = string("SC ") + locParticleROOTName + string(" Candidates;p (GeV/c);Timing PID Confidence Level");
1591	dHistMap_TimeFOMVsP[SYS_START][locPID] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNum2DFOMBins, 0.0, 1.0);
1592	*/
1593	gDirectory(TDirectory::CurrentDirectory())->cd("..");
1594
1595	//TOF
1596	CreateAndChangeTo_Directory("TOF", "TOF");
1597
1598	locHistName = string("DeltadEdXVsP_") + locParticleName;
1599	locHistTitle = locParticleROOTName + " Candidates;p (GeV/c);TOF #Delta(dE/dX) (MeV/cm)";
1600	dHistMap_DeltadEdXVsP[SYS_TOF][locPID] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNum2DDeltadEdxBins, dMinDeltadEdx, dMaxDeltadEdx);
1601
1602	locHistName = string("DeltaBetaVsP_") + locParticleName;
1603	locHistTitle = string("TOF ") + locParticleROOTName + string(" Candidates;p (GeV/c);#Delta#beta");
1604	dHistMap_DeltaBetaVsP[SYS_TOF][locPID] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNum2DDeltaBetaBins, dMinDeltaBeta, dMaxDeltaBeta);
1605
1606	locHistName = string("DeltaTVsP_") + locParticleName;
1607	locHistTitle = string("TOF ") + locParticleROOTName + string(" Candidates;p (GeV/c);#Deltat_{TOF - RF}");
1608	dHistMap_DeltaTVsP[SYS_TOF][locPID] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNum2DDeltaTBins, dMinDeltaT, dMaxDeltaT);
1609
1610	/*
1611	//Uncomment when ready!
1612	locHistName = string("TimePullVsP_") + locParticleName;
1613	locHistTitle = string("TOF ") + locParticleROOTName + string(" Candidates;p (GeV/c);#Deltat/#sigma_{#Deltat}");
1614	dHistMap_TimePullVsP[SYS_TOF][locPID] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNum2DPullBins, dMinPull, dMaxPull);
1615
1616	locHistName = string("TimeFOMVsP_") + locParticleName;
1617	locHistTitle = string("TOF ") + locParticleROOTName + string(" Candidates;p (GeV/c);Timing PID Confidence Level");
1618	dHistMap_TimeFOMVsP[SYS_TOF][locPID] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNum2DFOMBins, 0.0, 1.0);
1619	*/
1620	gDirectory(TDirectory::CurrentDirectory())->cd("..");
1621
1622	//BCAL
1623	CreateAndChangeTo_Directory("BCAL", "BCAL");
1624
1625	locHistName = string("DeltaBetaVsP_") + locParticleName;
1626	locHistTitle = string("BCAL ") + locParticleROOTName + string(" Candidates;p (GeV/c);#Delta#beta");
1627	dHistMap_DeltaBetaVsP[SYS_BCAL][locPID] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxBCALP, dNum2DDeltaBetaBins, dMinDeltaBeta, dMaxDeltaBeta);
1628
1629	locHistName = string("DeltaTVsP_") + locParticleName;
1630	locHistTitle = string("BCAL ") + locParticleROOTName + string(" Candidates;p (GeV/c);#Deltat_{BCAL - RF}");
1631	dHistMap_DeltaTVsP[SYS_BCAL][locPID] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNum2DDeltaTBins, dMinDeltaT, dMaxDeltaT);
1632
1633	/*
1634	//Uncomment when ready!
1635	locHistName = string("TimePullVsP_") + locParticleName;
1636	locHistTitle = string("BCAL ") + locParticleROOTName + string(" Candidates;p (GeV/c);#Deltat/#sigma_{#Deltat}");
1637	dHistMap_TimePullVsP[SYS_BCAL][locPID] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNum2DPullBins, dMinPull, dMaxPull);
1638
1639	locHistName = string("TimeFOMVsP_") + locParticleName;
1640	locHistTitle = string("BCAL ") + locParticleROOTName + string(" Candidates;p (GeV/c);Timing PID Confidence Level");
1641	dHistMap_TimeFOMVsP[SYS_BCAL][locPID] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNum2DFOMBins, 0.0, 1.0);
1642	*/
1643	gDirectory(TDirectory::CurrentDirectory())->cd("..");
1644
1645	//FCAL
1646	CreateAndChangeTo_Directory("FCAL", "FCAL");
1647
1648	locHistName = string("DeltaBetaVsP_") + locParticleName;
1649	locHistTitle = string("FCAL ") + locParticleROOTName + string(" Candidates;p (GeV/c);#Delta#beta");
1650	dHistMap_DeltaBetaVsP[SYS_FCAL][locPID] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNum2DDeltaBetaBins, dMinDeltaBeta, dMaxDeltaBeta);
1651
1652	locHistName = string("DeltaTVsP_") + locParticleName;
1653	locHistTitle = string("FCAL ") + locParticleROOTName + string(" Candidates;p (GeV/c);#Deltat_{FCAL - RF}");
1654	dHistMap_DeltaTVsP[SYS_FCAL][locPID] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNum2DDeltaTBins, dMinDeltaT, dMaxDeltaT);
1655
1656	/*
1657	//Uncomment when ready!
1658	locHistName = string("TimePullVsP_") + locParticleName;
1659	locHistTitle = string("FCAL ") + locParticleROOTName + string(" Candidates;p (GeV/c);#Deltat/#sigma_{#Deltat}");
1660	dHistMap_TimePullVsP[SYS_FCAL][locPID] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNum2DPullBins, dMinPull, dMaxPull);
1661
1662	locHistName = string("TimeFOMVsP_") + locParticleName;
1663	locHistTitle = string("FCAL ") + locParticleROOTName + string(" Candidates;p (GeV/c);Timing PID Confidence Level");
1664	dHistMap_TimeFOMVsP[SYS_FCAL][locPID] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNum2DFOMBins, 0.0, 1.0);
1665	*/
1666	gDirectory(TDirectory::CurrentDirectory())->cd("..");
1667
1668	//CDC
1669	CreateAndChangeTo_Directory("CDC", "CDC");
1670
1671	locHistName = string("DeltadEdXVsP_") + locParticleName;
1672	locHistTitle = locParticleROOTName + string(" Candidates;p (GeV/c);CDC #Delta(dE/dX) (keV/cm)");
1673	dHistMap_DeltadEdXVsP[SYS_CDC][locPID] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNum2DDeltadEdxBins, dMinDeltadEdx, dMaxDeltadEdx);
1674
1675	/*
1676	//Uncomment when ready!
1677	locHistName = string("dEdXPullVsP_") + locParticleName;
1678	locHistTitle = locParticleROOTName + string(" Candidates;p (GeV/c);CDC #Delta(dE/dX)/#sigma_{#Delta(dE/dX)}");
1679	dHistMap_dEdXPullVsP[SYS_CDC][locPID] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNum2DPullBins, dMinPull, dMaxPull);
1680
1681	locHistName = string("dEdXFOMVsP_") + locParticleName;
1682	locHistTitle = locParticleROOTName + string(" Candidates;p (GeV/c);CDC dE/dx PID Confidence Level");
1683	dHistMap_dEdXFOMVsP[SYS_CDC][locPID] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNum2DFOMBins, 0.0, 1.0);
1684	*/
1685	gDirectory(TDirectory::CurrentDirectory())->cd("..");
1686
1687	//FDC
1688	CreateAndChangeTo_Directory("FDC", "FDC");
1689
1690	locHistName = string("DeltadEdXVsP_") + locParticleName;
1691	locHistTitle = locParticleROOTName + string(" Candidates;p (GeV/c);FDC #Delta(dE/dX) (keV/cm)");
1692	dHistMap_DeltadEdXVsP[SYS_FDC][locPID] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNum2DDeltadEdxBins, dMinDeltadEdx, dMaxDeltadEdx);
1693
1694	/*
1695	//Uncomment when ready!
1696	locHistName = string("dEdXPullVsP_") + locParticleName;
1697	locHistTitle = locParticleROOTName + string(" Candidates;p (GeV/c);FDC #Delta(dE/dX)/#sigma_{#Delta(dE/dX)}");
1698	dHistMap_dEdXPullVsP[SYS_FDC][locPID] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNum2DPullBins, dMinPull, dMaxPull);
1699
1700	locHistName = string("dEdXFOMVsP_") + locParticleName;
1701	locHistTitle = locParticleROOTName + string(" Candidates;p (GeV/c);FDC dE/dx PID Confidence Level");
1702	dHistMap_dEdXFOMVsP[SYS_FDC][locPID] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNum2DFOMBins, 0.0, 1.0);
1703	*/
1704	gDirectory(TDirectory::CurrentDirectory())->cd("..");
1705	}
1706
1707	//Return to the base directory
1708	ChangeTo_BaseDirectory();
1709	}
1710	japp->RootUnLock(); //RELEASE ROOT LOCK!!
1711	}
1712
1713	bool DHistogramAction_DetectorPID::Perform_Action(JEventLoop* locEventLoop, const DParticleCombo* locParticleCombo)
1714	{
1715	//Expect locParticleCombo to be NULL since this is a reaction-independent action.
1716
1717	//Optional: Quit the action if it has already been executed this event (else may result in double-counting when filling histograms)
1718	if(Get_NumPreviousParticleCombos() != 0)
1719	return true;
1720
1721	vector<const DChargedTrack*> locChargedTracks;
1722	locEventLoop->Get(locChargedTracks, dTrackSelectionTag.c_str());
1723
1724	vector<const DNeutralParticle*> locNeutralParticles;
1725	locEventLoop->Get(locNeutralParticles, dShowerSelectionTag.c_str());
1726
1727	const DDetectorMatches* locDetectorMatches = NULL__null;
1728	locEventLoop->GetSingle(locDetectorMatches);
1729
1730	const DEventRFBunch* locEventRFBunch = NULL__null;
1731	locEventLoop->GetSingle(locEventRFBunch);
1732
1733	const DParticleID* locParticleID = NULL__null;
1734	locEventLoop->GetSingle(locParticleID);
1735
1736	//Fill Histograms
1737	japp->RootWriteLock(); //ACQUIRE ROOT LOCK!!
1738	{
1739	if(locEventRFBunch->dTimeSource != SYS_NULL) //only histogram beta for neutrals if the t0 is well known
1740	{
1741	for(size_t loc_i = 0; loc_i < locNeutralParticles.size(); ++loc_i)
1742	{
1743	//doesn't matter which hypothesis you use for beta: t0 is from DEventVertex time
1744	const DNeutralParticleHypothesis* locNeutralParticleHypothesis = locNeutralParticles[loc_i]->Get_BestFOM();
1745	double locBeta_Timing = locNeutralParticleHypothesis->measuredBeta();
1746	const DNeutralShower* locNeutralShower = locNeutralParticles[loc_i]->dNeutralShower;
1747	double locShowerEnergy = locNeutralShower->dEnergy;
1748
1749	double locDeltaT = locNeutralParticleHypothesis->time() - locEventRFBunch->dTime;
1750	if(locNeutralShower->dDetectorSystem == SYS_BCAL)
1751	{
1752	dHistMap_BetaVsP[SYS_BCAL][0]->Fill(locShowerEnergy, locBeta_Timing);
1753	dHistMap_DeltaTVsP[SYS_BCAL][Gamma]->Fill(locShowerEnergy, locDeltaT);
1754	}
1755	else
1756	{
1757	dHistMap_BetaVsP[SYS_FCAL][0]->Fill(locShowerEnergy, locBeta_Timing);
1758	dHistMap_DeltaTVsP[SYS_FCAL][Gamma]->Fill(locShowerEnergy, locDeltaT);
1759	}
1760	}
1761	}
1762
1763	for(size_t loc_i = 0; loc_i < locChargedTracks.size(); ++loc_i)
1764	{
1765	const DChargedTrackHypothesis* locChargedTrackHypothesis = locChargedTracks[loc_i]->Get_BestTrackingFOM();
1766	int locCharge = ParticleCharge(locChargedTrackHypothesis->PID());
1767	if(dHistMap_dEdXVsP[SYS_START].find(locCharge) == dHistMap_dEdXVsP[SYS_START].end())
1768	continue;
1769
1770	double locStartTime = locParticleID->Calc_PropagatedRFTime(locChargedTrackHypothesis, locEventRFBunch);
1771
1772	const DTrackTimeBased* locTrackTimeBased = NULL__null;
1773	locChargedTrackHypothesis->GetSingle(locTrackTimeBased);
1774
1775	Particle_t locPID = locChargedTrackHypothesis->PID();
1776	double locP = locTrackTimeBased->momentum().Mag();
1777	double locTheta = locTrackTimeBased->momentum().Theta()*180.0/TMath::Pi();
1778
1779	//if RF time is indeterminate, start time will be NaN
1780	const DBCALShowerMatchParams* locBCALShowerMatchParams = locChargedTrackHypothesis->Get_BCALShowerMatchParams();
1781	const DFCALShowerMatchParams* locFCALShowerMatchParams = locChargedTrackHypothesis->Get_FCALShowerMatchParams();
1782	const DTOFHitMatchParams* locTOFHitMatchParams = locChargedTrackHypothesis->Get_TOFHitMatchParams();
1783	const DSCHitMatchParams* locSCHitMatchParams = locChargedTrackHypothesis->Get_SCHitMatchParams();
1784
1785	if(locSCHitMatchParams != NULL__null)
1786	{
1787	dHistMap_dEdXVsP[SYS_START][locCharge]->Fill(locP, locSCHitMatchParams->dEdx*1.0E3);
1788	if((locEventRFBunch->dTimeSource != SYS_START) && (locEventRFBunch->dNumParticleVotes > 1))
1789	{
1790	//If SC was used for RF time, don't compute delta-beta
1791	double locBeta_Timing = locSCHitMatchParams->dPathLength/(29.9792458*(locSCHitMatchParams->dHitTime - locChargedTrackHypothesis->t0()));
1792	dHistMap_BetaVsP[SYS_START][locCharge]->Fill(locP, locBeta_Timing);
1793	if(dHistMap_DeltaBetaVsP[SYS_START].find(locPID) != dHistMap_DeltaBetaVsP[SYS_START].end())
1794	{
1795	double locDeltaBeta = locChargedTrackHypothesis->lorentzMomentum().Beta() - locBeta_Timing;
1796	dHistMap_DeltaBetaVsP[SYS_START][locPID]->Fill(locP, locDeltaBeta);
1797	double locDeltaT = locSCHitMatchParams->dHitTime - locSCHitMatchParams->dFlightTime - locStartTime;
1798	dHistMap_DeltaTVsP[SYS_START][locPID]->Fill(locP, locDeltaT);
1799	}
1800	}
1801	if(dHistMap_DeltadEdXVsP[SYS_START].find(locPID) != dHistMap_DeltadEdXVsP[SYS_START].end())
1802	{
1803	double locdx = locSCHitMatchParams->dHitEnergy/locSCHitMatchParams->dEdx;
1804	double locProbabledEdx = 0.0, locSigmadEdx = 0.0;
1805	locParticleID->GetScintMPdEandSigma(locP, locChargedTrackHypothesis->mass(), locdx, locProbabledEdx, locSigmadEdx);
1806	dHistMap_DeltadEdXVsP[SYS_START][locPID]->Fill(locP, (locSCHitMatchParams->dEdx - locProbabledEdx)*1.0E3);
1807	}
1808	}
1809	if(locTOFHitMatchParams != NULL__null)
1810	{
1811	dHistMap_dEdXVsP[SYS_TOF][locCharge]->Fill(locP, locTOFHitMatchParams->dEdx*1.0E3);
1812	if(locEventRFBunch->dNumParticleVotes > 1)
1813	{
1814	double locBeta_Timing = locTOFHitMatchParams->dPathLength/(29.9792458*(locTOFHitMatchParams->dHitTime - locChargedTrackHypothesis->t0()));
1815	dHistMap_BetaVsP[SYS_TOF][locCharge]->Fill(locP, locBeta_Timing);
1816	if(dHistMap_DeltaBetaVsP[SYS_TOF].find(locPID) != dHistMap_DeltaBetaVsP[SYS_TOF].end())
1817	{
1818	double locDeltaBeta = locChargedTrackHypothesis->lorentzMomentum().Beta() - locBeta_Timing;
1819	dHistMap_DeltaBetaVsP[SYS_TOF][locPID]->Fill(locP, locDeltaBeta);
1820	double locDeltaT = locTOFHitMatchParams->dHitTime - locTOFHitMatchParams->dFlightTime - locStartTime;
1821	dHistMap_DeltaTVsP[SYS_TOF][locPID]->Fill(locP, locDeltaT);
1822	}
1823	}
1824	if(dHistMap_DeltadEdXVsP[SYS_TOF].find(locPID) != dHistMap_DeltadEdXVsP[SYS_TOF].end())
1825	{
1826	double locdx = locTOFHitMatchParams->dHitEnergy/locTOFHitMatchParams->dEdx;
1827	double locProbabledEdx = 0.0, locSigmadEdx = 0.0;
1828	locParticleID->GetScintMPdEandSigma(locP, locChargedTrackHypothesis->mass(), locdx, locProbabledEdx, locSigmadEdx);
1829	dHistMap_DeltadEdXVsP[SYS_TOF][locPID]->Fill(locP, (locTOFHitMatchParams->dEdx - locProbabledEdx)*1.0E3);
1830	}
1831	}
1832	if(locBCALShowerMatchParams != NULL__null)
1833	{
1834	const DBCALShower* locBCALShower = locBCALShowerMatchParams->dBCALShower;
1835	double locEOverP = locBCALShower->E/locP;
1836	dHistMap_BCALEOverPVsP[locCharge]->Fill(locP, locEOverP);
1837	dHistMap_BCALEOverPVsTheta[locCharge]->Fill(locTheta, locEOverP);
1838	if(locEventRFBunch->dNumParticleVotes > 1)
1839	{
1840	double locBeta_Timing = locBCALShowerMatchParams->dPathLength/(29.9792458*(locBCALShower->t - locChargedTrackHypothesis->t0()));
1841	dHistMap_BetaVsP[SYS_BCAL][locCharge]->Fill(locP, locBeta_Timing);
1842	if(dHistMap_DeltaBetaVsP[SYS_BCAL].find(locPID) != dHistMap_DeltaBetaVsP[SYS_BCAL].end())
1843	{
1844	double locDeltaBeta = locChargedTrackHypothesis->lorentzMomentum().Beta() - locBeta_Timing;
1845	dHistMap_DeltaBetaVsP[SYS_BCAL][locPID]->Fill(locP, locDeltaBeta);
1846	double locDeltaT = locBCALShower->t - locBCALShowerMatchParams->dFlightTime - locStartTime;
1847	dHistMap_DeltaTVsP[SYS_BCAL][locPID]->Fill(locP, locDeltaT);
1848	}
1849	}
1850	}
1851	if(locFCALShowerMatchParams != NULL__null)
1852	{
1853	const DFCALShower* locFCALShower = locFCALShowerMatchParams->dFCALShower;
1854	double locEOverP = locFCALShower->getEnergy()/locP;
1855	dHistMap_FCALEOverPVsP[locCharge]->Fill(locP, locEOverP);
1856	dHistMap_FCALEOverPVsTheta[locCharge]->Fill(locTheta, locEOverP);
1857	if(locEventRFBunch->dNumParticleVotes > 1)
1858	{
1859	double locBeta_Timing = locFCALShowerMatchParams->dPathLength/(29.9792458*(locFCALShower->getTime() - locChargedTrackHypothesis->t0()));
1860	dHistMap_BetaVsP[SYS_FCAL][locCharge]->Fill(locP, locBeta_Timing);
1861	if(dHistMap_DeltaBetaVsP[SYS_FCAL].find(locPID) != dHistMap_DeltaBetaVsP[SYS_FCAL].end())
1862	{
1863	double locDeltaBeta = locChargedTrackHypothesis->lorentzMomentum().Beta() - locBeta_Timing;
1864	dHistMap_DeltaBetaVsP[SYS_FCAL][locPID]->Fill(locP, locDeltaBeta);
1865	double locDeltaT = locFCALShower->getTime() - locFCALShowerMatchParams->dFlightTime - locStartTime;
1866	dHistMap_DeltaTVsP[SYS_FCAL][locPID]->Fill(locP, locDeltaT);
1867	}
1868	}
1869	}
1870
1871	if(locTrackTimeBased->dNumHitsUsedFordEdx_CDC > 0)
1872	{
1873	dHistMap_dEdXVsP[SYS_CDC][locCharge]->Fill(locP, locTrackTimeBased->ddEdx_CDC*1.0E6);
1874	if(dHistMap_DeltadEdXVsP[SYS_CDC].find(locPID) != dHistMap_DeltadEdXVsP[SYS_CDC].end())
1875	{
1876	double locProbabledEdx = locParticleID->GetMostProbabledEdx_DC(locP, locChargedTrackHypothesis->mass(), locTrackTimeBased->ddx_CDC, true);
1877	dHistMap_DeltadEdXVsP[SYS_CDC][locPID]->Fill(locP, (locTrackTimeBased->ddEdx_CDC - locProbabledEdx)*1.0E6);
1878	}
1879	}
1880	if(locTrackTimeBased->dNumHitsUsedFordEdx_FDC > 0)
1881	{
1882	dHistMap_dEdXVsP[SYS_FDC][locCharge]->Fill(locP, locTrackTimeBased->ddEdx_FDC*1.0E6);
1883	if(dHistMap_DeltadEdXVsP[SYS_FDC].find(locPID) != dHistMap_DeltadEdXVsP[SYS_FDC].end())
1884	{
1885	double locProbabledEdx = locParticleID->GetMostProbabledEdx_DC(locP, locChargedTrackHypothesis->mass(), locTrackTimeBased->ddx_FDC, false);
1886	dHistMap_DeltadEdXVsP[SYS_FDC][locPID]->Fill(locP, (locTrackTimeBased->ddEdx_FDC - locProbabledEdx)*1.0E6);
1887	}
1888	}
1889	}
1890	}
1891	japp->RootUnLock(); //RELEASE ROOT LOCK!!
1892
1893	return true; //return false if you want to use this action to apply a cut (and it fails the cut!)
1894	}
1895
1896	void DHistogramAction_Neutrals::Initialize(JEventLoop* locEventLoop)
1897	{
1898	//Create any histograms/trees/etc. within a ROOT lock.
1899	//This is so that when running multithreaded, only one thread is writing to the ROOT file at a time.
1900
1901	//When creating a reaction-independent action, only modify member variables within a ROOT lock.
1902	//Objects created within a plugin (such as reaction-independent actions) can be accessed by many threads simultaneously.
1903
1904	string locHistName;
1905
1906	DApplication* locApplication = dynamic_cast<DApplication*>(locEventLoop->GetJApplication());
1907	DGeometry* locGeometry = locApplication->GetDGeometry(locEventLoop->GetJEvent().GetRunNumber());
1908	double locTargetZCenter = 0.0;
1909	locGeometry->GetTargetZ(locTargetZCenter);
1910
1911	japp->RootWriteLock(); //ACQUIRE ROOT LOCK!!
1912	{
1913	//Required: Create a folder in the ROOT output file that will contain all of the output ROOT objects (if any) for this action.
1914	//If another thread has already created the folder, it just changes to it.
1915	CreateAndChangeTo_ActionDirectory();
1916
1917	if(dTargetCenter.Z() < -9.8E9)
1918	dTargetCenter.SetXYZ(0.0, 0.0, locTargetZCenter);
1919
1920	//BCAL
1921	locHistName = "BCALTrackDOCA";
1922	dHist_BCALTrackDOCA = GetOrCreate_Histogram<TH1I>(locHistName, ";BCAL Shower Distance to Nearest Track (cm)", dNumTrackDOCABins, dMinTrackDOCA, dMaxTrackDOCA);
1923	locHistName = "BCALTrackDeltaPhi";
1924	dHist_BCALTrackDeltaPhi = GetOrCreate_Histogram<TH1I>(locHistName, ";BCAL Shower #Delta#phi#circ to Nearest Track", dNumDeltaPhiBins, dMinDeltaPhi, dMaxDeltaPhi);
1925	locHistName = "BCALTrackDeltaZ";
1926	dHist_BCALTrackDeltaZ = GetOrCreate_Histogram<TH1I>(locHistName, ";BCAL Shower #DeltaZ to Nearest Track (cm)", dNumTrackDOCABins, dMinTrackDOCA, dMaxTrackDOCA);
1927	locHistName = "BCALNeutralShowerEnergy";
1928	dHist_BCALNeutralShowerEnergy = GetOrCreate_Histogram<TH1I>(locHistName, ";BCAL Neutral Shower Energy (GeV)", dNumShowerEnergyBins, dMinShowerEnergy, dMaxBCALP);
1929	locHistName = "BCALNeutralShowerDeltaT";
1930	dHist_BCALNeutralShowerDeltaT = GetOrCreate_Histogram<TH1I>(locHistName, ";BCAL Neutral Shower #Deltat (Propagated - RF) (ns)", dNumDeltaTBins, dMinDeltaT, dMaxDeltaT);
1931	locHistName = "BCALNeutralShowerDeltaTVsE";
1932	dHist_BCALNeutralShowerDeltaTVsE = GetOrCreate_Histogram<TH2I>(locHistName, ";BCAL Neutral Shower Energy (GeV);BCAL Neutral Shower #Deltat (ns)", dNum2DShowerEnergyBins, dMinShowerEnergy, dMaxBCALP, dNum2DDeltaTBins, dMinDeltaT, dMaxDeltaT);
1933	locHistName = "BCALNeutralShowerDeltaTVsZ";
1934	dHist_BCALNeutralShowerDeltaTVsZ = GetOrCreate_Histogram<TH2I>(locHistName, ";BCAL Neutral Shower Z (cm);BCAL Neutral Shower #Deltat (ns)", dNum2DBCALZBins, 0.0, 450.0, dNum2DDeltaTBins, dMinDeltaT, dMaxDeltaT);
1935
1936	//FCAL
1937	locHistName = "FCALTrackDOCA";
1938	dHist_FCALTrackDOCA = GetOrCreate_Histogram<TH1I>(locHistName, ";FCAL Shower Distance to Nearest Track (cm)", dNumTrackDOCABins, dMinTrackDOCA, dMaxTrackDOCA);
1939	locHistName = "FCALNeutralShowerEnergy";
1940	dHist_FCALNeutralShowerEnergy = GetOrCreate_Histogram<TH1I>(locHistName, ";FCAL Neutral Shower Energy (GeV)", dNumShowerEnergyBins, dMinShowerEnergy, dMaxShowerEnergy);
1941	locHistName = "FCALNeutralShowerDeltaT";
1942	dHist_FCALNeutralShowerDeltaT = GetOrCreate_Histogram<TH1I>(locHistName, ";FCAL Neutral Shower #Deltat (Propagated - RF) (ns)", dNumDeltaTBins, dMinDeltaT, dMaxDeltaT);
1943	locHistName = "FCALNeutralShowerDeltaTVsE";
1944	dHist_FCALNeutralShowerDeltaTVsE = GetOrCreate_Histogram<TH2I>(locHistName, ";FCAL Neutral Shower Energy (GeV);FCAL Neutral Shower #Deltat (ns)", dNum2DShowerEnergyBins, dMinShowerEnergy, dMaxShowerEnergy, dNum2DDeltaTBins, dMinDeltaT, dMaxDeltaT);
1945
1946	//Return to the base directory
1947	ChangeTo_BaseDirectory();
1948	}
1949	japp->RootUnLock(); //RELEASE ROOT LOCK!!
1950	}
1951
1952	bool DHistogramAction_Neutrals::Perform_Action(JEventLoop* locEventLoop, const DParticleCombo* locParticleCombo)
1953	{
1954	//Expect locParticleCombo to be NULL since this is a reaction-independent action.
1955
1956	//Optional: Quit the action if it has already been executed this event (else may result in double-counting when filling histograms)
1957	if(Get_NumPreviousParticleCombos() != 0)
1958	return true;
1959
1960	vector<const DNeutralShower*> locNeutralShowers;
1961	locEventLoop->Get(locNeutralShowers);
1962
1963	vector<const DTrackTimeBased*> locTrackTimeBasedVector;
1964	locEventLoop->Get(locTrackTimeBasedVector);
1965
1966	const DDetectorMatches* locDetectorMatches = NULL__null;
1967	locEventLoop->GetSingle(locDetectorMatches);
1968
1969	vector<const DEventRFBunch*> locEventRFBunches;
1970	locEventLoop->Get(locEventRFBunches);
1971	double locStartTime = locEventRFBunches.empty() ? 0.0 : locEventRFBunches[0]->dTime;
1972
1973	//Fill Histograms
1974	japp->RootWriteLock(); //ACQUIRE ROOT LOCK!!
1975	{
1976	for(size_t loc_i = 0; loc_i < locNeutralShowers.size(); ++loc_i)
1977	{
1978	//assume is photon
1979	double locPathLength = (locNeutralShowers[loc_i]->dSpacetimeVertex.Vect() - dTargetCenter).Mag();
1980	double locDeltaT = locNeutralShowers[loc_i]->dSpacetimeVertex.T() - locPathLength/29.9792458 - locStartTime;
1981
1982	if(locNeutralShowers[loc_i]->dDetectorSystem == SYS_FCAL)
1983	{
1984	const DFCALShower* locFCALShower = NULL__null;
1985	locNeutralShowers[loc_i]->GetSingle(locFCALShower);
1986
1987	double locDistance = 9.9E9;
1988	if(locDetectorMatches->Get_DistanceToNearestTrack(locFCALShower, locDistance))
1989	dHist_FCALTrackDOCA->Fill(locDistance);
1990
1991	dHist_FCALNeutralShowerEnergy->Fill(locNeutralShowers[loc_i]->dEnergy);
1992	dHist_FCALNeutralShowerDeltaT->Fill(locDeltaT);
1993	dHist_FCALNeutralShowerDeltaTVsE->Fill(locNeutralShowers[loc_i]->dEnergy, locDeltaT);
1994	}
1995	else
1996	{
1997	const DBCALShower* locBCALShower = NULL__null;
1998	locNeutralShowers[loc_i]->GetSingle(locBCALShower);
1999
2000	double locDistance = 9.9E9, locDeltaPhi = 9.9E9, locDeltaZ = 9.9E9;
2001	if(locDetectorMatches->Get_DistanceToNearestTrack(locBCALShower, locDistance))
2002	dHist_BCALTrackDOCA->Fill(locDistance);
2003	if(locDetectorMatches->Get_DistanceToNearestTrack(locBCALShower, locDeltaPhi, locDeltaZ))
2004	{
2005	dHist_BCALTrackDeltaPhi->Fill(180.0*locDeltaPhi/TMath::Pi());
2006	dHist_BCALTrackDeltaZ->Fill(locDeltaZ);
2007	}
2008
2009	dHist_BCALNeutralShowerEnergy->Fill(locNeutralShowers[loc_i]->dEnergy);
2010	dHist_BCALNeutralShowerDeltaT->Fill(locDeltaT);
2011	dHist_BCALNeutralShowerDeltaTVsE->Fill(locNeutralShowers[loc_i]->dEnergy, locDeltaT);
2012	dHist_BCALNeutralShowerDeltaTVsZ->Fill(locNeutralShowers[loc_i]->dSpacetimeVertex.Z(), locDeltaT);
2013	}
2014	}
2015	}
2016	japp->RootUnLock(); //RELEASE ROOT LOCK!!
2017
2018	return true; //return false if you want to use this action to apply a cut (and it fails the cut!)
2019	}
2020
2021
2022	void DHistogramAction_DetectorMatchParams::Initialize(JEventLoop* locEventLoop)
2023	{
2024	//Create any histograms/trees/etc. within a ROOT lock.
2025	//This is so that when running multithreaded, only one thread is writing to the ROOT file at a time.
2026
2027	//When creating a reaction-independent action, only modify member variables within a ROOT lock.
2028	//Objects created within a plugin (such as reaction-independent actions) can be accessed by many threads simultaneously.
2029
2030	string locHistName, locHistTitle;
2031
2032	vector<const DMCThrown*> locMCThrowns;
2033	locEventLoop->Get(locMCThrowns);
2034
2035	DApplication* locApplication = dynamic_cast<DApplication*>(locEventLoop->GetJApplication());
2036	DGeometry* locGeometry = locApplication->GetDGeometry(locEventLoop->GetJEvent().GetRunNumber());
2037	double locTargetZCenter = 0.0;
2038	locGeometry->GetTargetZ(locTargetZCenter);
2039
2040	string locTrackSelectionTag = "NotATag";
2041	if(gPARMS->Exists("COMBO:TRACK_SELECT_TAG"))
2042	gPARMS->GetParameter("COMBO:TRACK_SELECT_TAG", locTrackSelectionTag);
2043
2044	japp->RootWriteLock(); //ACQUIRE ROOT LOCK!!
2045	{
2046	//So: Default tag is "", User can set it to something else
2047	//In here, if tag is "", get from gparms, if not, leave it alone
2048	//If gparms value does not exist, set it to (and use) "PreSelect"
2049	if(dTrackSelectionTag == "NotATag")
2050	dTrackSelectionTag = (locTrackSelectionTag == "NotATag") ? "PreSelect" : locTrackSelectionTag;
2051
2052	//Required: Create a folder in the ROOT output file that will contain all of the output ROOT objects (if any) for this action.
2053	//If another thread has already created the folder, it just changes to it.
2054	CreateAndChangeTo_ActionDirectory();
2055
2056	if(dTargetCenterZ < -9.8E9)
2057	dTargetCenterZ = locTargetZCenter; //only set if not already set
2058
2059	//Track Matched to Hit
2060	for(int locTruePIDFlag = 0; locTruePIDFlag < 2; ++locTruePIDFlag)
2061	{
2062	if(locMCThrowns.empty() && (locTruePIDFlag == 1))
2063	continue; //not a simulated event: don't histogram thrown info!
2064
2065	string locDirName = (locTruePIDFlag == 1) ? "TruePID" : "ReconstructedPID";
2066	CreateAndChangeTo_Directory(locDirName.c_str(), locDirName.c_str());
2067
2068	//By PID
2069	for(int loc_i = -2; loc_i < int(dTrackingPIDs.size()); ++loc_i) //-2 = q-, -1 = q+
2070	{
2071	string locParticleName, locParticleROOTName;
2072	int locPID = loc_i;
2073	if(loc_i == -2)
2074	{
2075	locParticleName = "q-";
2076	locParticleROOTName = "#it{q}^{-}";
2077	}
2078	else if(loc_i == -1)
2079	{
2080	locParticleName = "q+";
2081	locParticleROOTName = "#it{q}^{+}";
2082	}
2083	else
2084	{
2085	locParticleName = ParticleType(dTrackingPIDs[loc_i]);
2086	locParticleROOTName = ParticleName_ROOT(dTrackingPIDs[loc_i]);
2087	locPID = int(dTrackingPIDs[loc_i]);
2088	}
2089	CreateAndChangeTo_Directory(locParticleName, locParticleName);
2090	pair<int, bool> locPIDPair(locPID, bool(locTruePIDFlag));
2091
2092	//BCAL
2093	locHistName = "BCALShowerEnergy";
2094	locHistTitle = locParticleROOTName + ";BCAL Shower Energy (GeV)";
2095	dHistMap_BCALShowerEnergy[locPIDPair] = GetOrCreate_Histogram<TH1I>(locHistName, locHistTitle, dNumShowerEnergyBins, dMinShowerEnergy, dMaxBCALP);
2096
2097	locHistName = "BCALShowerTrackDepth";
2098	locHistTitle = locParticleROOTName + ";BCAL Shower Track Depth (cm)";
2099	dHistMap_BCALShowerTrackDepth[locPIDPair] = GetOrCreate_Histogram<TH1I>(locHistName, locHistTitle, dNumShowerDepthBins, dMinShowerDepth, dMaxShowerDepth);
2100
2101	locHistName = "BCALShowerTrackDepthVsP";
2102	locHistTitle = locParticleROOTName + ";p (GeV/c);BCAL Shower Track Depth (cm)";
2103	dHistMap_BCALShowerTrackDepthVsP[locPIDPair] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxBCALP, dNumShowerDepthBins, dMinShowerDepth, dMaxShowerDepth);
2104
2105	//FCAL
2106	locHistName = "FCALShowerEnergy";
2107	locHistTitle = locParticleROOTName + ";FCAL Shower Energy (GeV)";
2108	dHistMap_FCALShowerEnergy[locPIDPair] = GetOrCreate_Histogram<TH1I>(locHistName, locHistTitle, dNumShowerEnergyBins, dMinShowerEnergy, dMaxShowerEnergy);
2109
2110	locHistName = "FCALShowerTrackDepth";
2111	locHistTitle = locParticleROOTName + ";FCAL Shower Track Depth (cm)";
2112	dHistMap_FCALShowerTrackDepth[locPIDPair] = GetOrCreate_Histogram<TH1I>(locHistName, locHistTitle, dNumShowerDepthBins, dMinShowerDepth, dMaxShowerDepth);
2113
2114	locHistName = "FCALShowerTrackDepthVsP";
2115	locHistTitle = locParticleROOTName + ";p (GeV/c);FCAL Shower Track Depth (cm)";
2116	dHistMap_FCALShowerTrackDepthVsP[locPIDPair] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNumShowerDepthBins, dMinShowerDepth, dMaxShowerDepth);
2117
2118	//SC
2119	locHistName = "SCEnergyVsTheta";
2120	locHistTitle = locParticleROOTName + ";#theta#circ;SC Point Energy (MeV)";
2121	dHistMap_SCEnergyVsTheta[locPIDPair] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DThetaBins, dMinTheta, dMaxTheta, dNum2DHitEnergyBins, dMinHitEnergy, dMaxHitEnergy);
2122
2123	locHistName = "SCPhiVsTheta";
2124	locHistTitle = locParticleROOTName + ";#theta#circ;#phi#circ";
2125	dHistMap_SCPhiVsTheta[locPIDPair] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DThetaBins, dMinTheta, dMaxTheta, dNum2DPhiBins, dMinPhi, dMaxPhi);
2126
2127	gDirectory(TDirectory::CurrentDirectory())->cd(".."); //end of PID
2128	}
2129	gDirectory(TDirectory::CurrentDirectory())->cd(".."); //end of true/recon PID
2130	}
2131
2132	//Return to the base directory
2133	ChangeTo_BaseDirectory();
2134	}
2135	japp->RootUnLock(); //RELEASE ROOT LOCK!!
2136	}
2137
2138	bool DHistogramAction_DetectorMatchParams::Perform_Action(JEventLoop* locEventLoop, const DParticleCombo* locParticleCombo)
2139	{
2140	//Expect locParticleCombo to be NULL since this is a reaction-independent action.
2141
2142	//Optional: Quit the action if it has already been executed this event (else may result in double-counting when filling histograms)
2143	if(Get_NumPreviousParticleCombos() != 0)
2144	return true;
2145
2146	vector<const DMCThrown*> locMCThrowns;
2147	locEventLoop->Get(locMCThrowns);
2148
2149	Fill_Hists(locEventLoop, false);
2150	if(!locMCThrowns.empty())
2151	Fill_Hists(locEventLoop, true);
2152
2153	return true; //return false if you want to use this action to apply a cut (and it fails the cut!)
2154	}
2155
2156	void DHistogramAction_DetectorMatchParams::Fill_Hists(JEventLoop* locEventLoop, bool locUseTruePIDFlag)
2157	{
2158	vector<const DChargedTrack*> locChargedTracks;
2159	locEventLoop->Get(locChargedTracks, dTrackSelectionTag.c_str());
2160
2161	vector<const DMCThrownMatching*> locMCThrownMatchingVector;
2162	locEventLoop->Get(locMCThrownMatchingVector);
2163
2164	const DEventRFBunch* locEventRFBunch = NULL__null;
2165	locEventLoop->GetSingle(locEventRFBunch);
2166
2167	//Fill Histograms
2168	japp->RootWriteLock(); //ACQUIRE ROOT LOCK!!
2169	{
2170	for(size_t loc_i = 0; loc_i < locChargedTracks.size(); ++loc_i)
2171	{
2172	const DChargedTrackHypothesis* locChargedTrackHypothesis = locChargedTracks[loc_i]->Get_BestFOM();
2173
2174	if(locUseTruePIDFlag && (!locMCThrownMatchingVector.empty()))
2175	{
2176	double locMatchFOM = 0.0;
2177	const DMCThrown* locMCThrown = locMCThrownMatchingVector[0]->Get_MatchingMCThrown(locChargedTrackHypothesis, locMatchFOM);
2178	if(locMCThrown == NULL__null)
2179	continue;
2180	//OK, have the thrown. Now, grab the best charged track hypothesis to get the best matching
2181	locChargedTrackHypothesis = locMCThrownMatchingVector[0]->Get_MatchingChargedHypothesis(locMCThrown, locMatchFOM);
2182	}
2183
2184	pair<int, bool> locPIDPair(int(locChargedTrackHypothesis->PID()), locUseTruePIDFlag);
2185	bool locDisregardPIDFlag = (dHistMap_BCALShowerEnergy.find(locPIDPair) == dHistMap_BCALShowerEnergy.end());
2186	int locQIndex = (locChargedTrackHypothesis->charge() > 0.0) ? -1 : -2;
2187	pair<int, bool> locChargePair(locQIndex, locUseTruePIDFlag);
2188
2189	DVector3 locMomentum = locChargedTrackHypothesis->momentum();
2190	const DFCALShowerMatchParams* locFCALShowerMatchParams = locChargedTrackHypothesis->Get_FCALShowerMatchParams();
2191	const DSCHitMatchParams* locSCHitMatchParams = locChargedTrackHypothesis->Get_SCHitMatchParams();
2192	const DBCALShowerMatchParams* locBCALShowerMatchParams = locChargedTrackHypothesis->Get_BCALShowerMatchParams();
2193
2194	//BCAL
2195	if(locBCALShowerMatchParams != NULL__null)
2196	{
2197	const DBCALShower* locBCALShower = locBCALShowerMatchParams->dBCALShower;
2198	dHistMap_BCALShowerEnergy[locChargePair]->Fill(locBCALShower->E);
2199	dHistMap_BCALShowerTrackDepth[locChargePair]->Fill(locBCALShowerMatchParams->dx);
2200	dHistMap_BCALShowerTrackDepthVsP[locChargePair]->Fill(locMomentum.Mag(), locBCALShowerMatchParams->dx);
2201
2202	if(!locDisregardPIDFlag)
2203	{
2204	dHistMap_BCALShowerEnergy[locPIDPair]->Fill(locBCALShower->E);
2205	dHistMap_BCALShowerTrackDepth[locPIDPair]->Fill(locBCALShowerMatchParams->dx);
2206	dHistMap_BCALShowerTrackDepthVsP[locPIDPair]->Fill(locMomentum.Mag(), locBCALShowerMatchParams->dx);
2207	}
2208	}
2209
2210	//FCAL
2211	if(locFCALShowerMatchParams != NULL__null)
2212	{
2213	const DFCALShower* locFCALShower = locFCALShowerMatchParams->dFCALShower;
2214	dHistMap_FCALShowerEnergy[locChargePair]->Fill(locFCALShower->getEnergy());
2215	dHistMap_FCALShowerTrackDepth[locChargePair]->Fill(locFCALShowerMatchParams->dx);
2216	dHistMap_FCALShowerTrackDepthVsP[locChargePair]->Fill(locMomentum.Mag(), locFCALShowerMatchParams->dx);
2217
2218	if(!locDisregardPIDFlag)
2219	{
2220	dHistMap_FCALShowerEnergy[locPIDPair]->Fill(locFCALShower->getEnergy());
2221	dHistMap_FCALShowerTrackDepth[locPIDPair]->Fill(locFCALShowerMatchParams->dx);
2222	dHistMap_FCALShowerTrackDepthVsP[locPIDPair]->Fill(locMomentum.Mag(), locFCALShowerMatchParams->dx);
2223	}
2224	}
2225
2226	//SC
2227	if(locSCHitMatchParams != NULL__null)
2228	{
2229	dHistMap_SCEnergyVsTheta[locChargePair]->Fill(locMomentum.Theta()180.0/TMath::Pi(), locSCHitMatchParams->dHitEnergy1.0E3);
2230	dHistMap_SCPhiVsTheta[locChargePair]->Fill(locMomentum.Theta()180.0/TMath::Pi(), locMomentum.Phi()180.0/TMath::Pi());
2231
2232	if(!locDisregardPIDFlag)
2233	{
2234	dHistMap_SCEnergyVsTheta[locPIDPair]->Fill(locMomentum.Theta()180.0/TMath::Pi(), locSCHitMatchParams->dHitEnergy1.0E3);
2235	dHistMap_SCPhiVsTheta[locPIDPair]->Fill(locMomentum.Theta()180.0/TMath::Pi(), locMomentum.Phi()180.0/TMath::Pi());
2236	}
2237	}
2238	}
2239	}
2240	japp->RootUnLock(); //RELEASE ROOT LOCK!!
2241	}
2242
2243	void DHistogramAction_EventVertex::Initialize(JEventLoop* locEventLoop)
2244	{
2245	string locHistName, locHistTitle;
2246
2247	string locTrackSelectionTag = "NotATag";
2248	if(gPARMS->Exists("COMBO:TRACK_SELECT_TAG"))
2249	gPARMS->GetParameter("COMBO:TRACK_SELECT_TAG", locTrackSelectionTag);
2250
2251	japp->RootWriteLock(); //ACQUIRE ROOT LOCK!!
2252	{
2253	//So: Default tag is "", User can set it to something else
2254	//In here, if tag is "", get from gparms, if not, leave it alone
2255	//If gparms value does not exist, set it to (and use) "PreSelect"
2256	if(dTrackSelectionTag == "NotATag")
2257	dTrackSelectionTag = (locTrackSelectionTag == "NotATag") ? "PreSelect" : locTrackSelectionTag;
2258
2259	CreateAndChangeTo_ActionDirectory();
2260
2261	// Event RF Bunch Time
2262	locHistName = "RFTrackDeltaT";
2263	locHistTitle = ";#Deltat_{RF - Track} (ns)";
2264	dRFTrackDeltaT = GetOrCreate_Histogram<TH1I>(locHistName, locHistTitle, dNumRFTBins, -3.0, 3.0);
2265
2266	// ALL EVENTS
2267	CreateAndChangeTo_Directory("AllEvents", "AllEvents");
2268
2269	// Event Vertex-Z
2270	locHistName = "EventVertexZ";
2271	locHistTitle = ";Event Vertex-Z (cm)";
2272	dEventVertexZ_AllEvents = GetOrCreate_Histogram<TH1I>(locHistName, locHistTitle, dNumVertexZBins, dMinVertexZ, dMaxVertexZ);
2273
2274	// Event Vertex-Y Vs Vertex-X
2275	locHistName = "EventVertexYVsX";
2276	locHistTitle = ";Event Vertex-X (cm);Event Vertex-Y (cm)";
2277	dEventVertexYVsX_AllEvents = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNumVertexXYBins, dMinVertexXY, dMaxVertexXY, dNumVertexXYBins, dMinVertexXY, dMaxVertexXY);
2278
2279	// Event Vertex-T
2280	locHistName = "EventVertexT";
2281	locHistTitle = ";Event Vertex Time (ns)";
2282	dEventVertexT_AllEvents = GetOrCreate_Histogram<TH1I>(locHistName, locHistTitle, dNumTBins, dMinT, dMaxT);
2283
2284	gDirectory(TDirectory::CurrentDirectory())->cd("..");
2285
2286
2287	// 2+ Good Tracks
2288	CreateAndChangeTo_Directory("2+GoodTracks", "2+GoodTracks");
2289
2290	// Event Vertex-Z
2291	locHistName = "EventVertexZ";
2292	locHistTitle = ";Event Vertex-Z (cm)";
2293	dEventVertexZ_2OrMoreGoodTracks = GetOrCreate_Histogram<TH1I>(locHistName, locHistTitle, dNumVertexZBins, dMinVertexZ, dMaxVertexZ);
2294
2295	// Event Vertex-Y Vs Vertex-X
2296	locHistName = "EventVertexYVsX";
2297	locHistTitle = ";Event Vertex-X (cm);Event Vertex-Y (cm)";
2298	dEventVertexYVsX_2OrMoreGoodTracks = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNumVertexXYBins, dMinVertexXY, dMaxVertexXY, dNumVertexXYBins, dMinVertexXY, dMaxVertexXY);
2299
2300	// Event Vertex-T
2301	locHistName = "EventVertexT";
2302	locHistTitle = ";Event Vertex Time (ns)";
2303	dEventVertexT_2OrMoreGoodTracks = GetOrCreate_Histogram<TH1I>(locHistName, locHistTitle, dNumTBins, dMinT, dMaxT);
2304
2305	// Confidence Level
2306	locHistName = "ConfidenceLevel";
2307	dHist_KinFitConfidenceLevel = GetOrCreate_Histogram<TH1I>(locHistName, "Event Vertex Kinematic Fit;Confidence Level;# Events", dNumConfidenceLevelBins, 0.0, 1.0);
2308
2309	//final particle pulls
2310	for(size_t loc_i = 0; loc_i < dFinalStatePIDs.size(); ++loc_i)
2311	{
2312	Particle_t locPID = dFinalStatePIDs[loc_i];
2313	string locParticleDirName = string("Pulls_") + string(ParticleType(locPID));
2314	string locParticleROOTName = ParticleName_ROOT(locPID);
2315	CreateAndChangeTo_Directory(locParticleDirName, locParticleDirName);
2316
2317	//Px Pull
2318	locHistName = "Pull_Px";
2319	locHistTitle = locParticleROOTName + string(", Vertex Fit;p_{x} Pull;# Events");
2320	dHistMap_KinFitPulls[locPID][d_PxPull] = GetOrCreate_Histogram<TH1I>(locHistName, locHistTitle, dNumPullBins, dMinPull, dMaxPull);
2321
2322	//Py Pull
2323	locHistName = "Pull_Py";
2324	locHistTitle = locParticleROOTName + string(", Vertex Fit;p_{y} Pull;# Events");
2325	dHistMap_KinFitPulls[locPID][d_PyPull] = GetOrCreate_Histogram<TH1I>(locHistName, locHistTitle, dNumPullBins, dMinPull, dMaxPull);
2326
2327	//Pz Pull
2328	locHistName = "Pull_Pz";
2329	locHistTitle = locParticleROOTName + string(", Vertex Fit;p_{z} Pull;# Events");
2330	dHistMap_KinFitPulls[locPID][d_PzPull] = GetOrCreate_Histogram<TH1I>(locHistName, locHistTitle, dNumPullBins, dMinPull, dMaxPull);
2331
2332	//Xx Pull
2333	locHistName = "Pull_Xx";
2334	locHistTitle = locParticleROOTName + string(", Vertex Fit;x_{x} Pull;# Events");
2335	dHistMap_KinFitPulls[locPID][d_XxPull] = GetOrCreate_Histogram<TH1I>(locHistName, locHistTitle, dNumPullBins, dMinPull, dMaxPull);
2336
2337	//Xy Pull
2338	locHistName = "Pull_Xy";
2339	locHistTitle = locParticleROOTName + string(", Vertex Fit;x_{y} Pull;# Events");
2340	dHistMap_KinFitPulls[locPID][d_XyPull] = GetOrCreate_Histogram<TH1I>(locHistName, locHistTitle, dNumPullBins, dMinPull, dMaxPull);
2341
2342	//Xz Pull
2343	locHistName = "Pull_Xz";
2344	locHistTitle = locParticleROOTName + string(", Vertex Fit;x_{z} Pull;# Events");
2345	dHistMap_KinFitPulls[locPID][d_XzPull] = GetOrCreate_Histogram<TH1I>(locHistName, locHistTitle, dNumPullBins, dMinPull, dMaxPull);
2346
2347	gDirectory(TDirectory::CurrentDirectory())->cd("..");
2348	}
2349
2350	gDirectory(TDirectory::CurrentDirectory())->cd("..");
2351
2352	//Return to the base directory
2353	ChangeTo_BaseDirectory();
2354	}
2355	japp->RootUnLock(); //RELEASE ROOT LOCK!!
2356	}
2357
2358	bool DHistogramAction_EventVertex::Perform_Action(JEventLoop* locEventLoop, const DParticleCombo* locParticleCombo)
2359	{
2360	if(Get_NumPreviousParticleCombos() != 0)
2361	return true; //else double-counting!
2362
2363	const DVertex* locVertex = NULL__null;
2364	locEventLoop->GetSingle(locVertex);
2365
2366	vector<const DChargedTrack*> locChargedTracks;
2367	locEventLoop->Get(locChargedTracks, dTrackSelectionTag.c_str());
2368
2369	const DDetectorMatches* locDetectorMatches = NULL__null;
2370	locEventLoop->GetSingle(locDetectorMatches);
2371
2372	const DParticleID* locParticleID = NULL__null;
2373	locEventLoop->GetSingle(locParticleID);
2374
2375	const DEventRFBunch* locEventRFBunch = NULL__null;
2376	locEventLoop->GetSingle(locEventRFBunch);
2377
2378	//Make sure that brun() is called (to get rf period) before using.
2379	//Cannot call JEventLoop->Get() because object may be in datastream (REST), bypassing factory brun() call.
2380	//Must do here rather than in Initialize() function because this object is shared by all threads (which each have their own factory)
2381	DRFTime_factory* locRFTimeFactory = static_cast<DRFTime_factory*>(locEventLoop->GetFactory("DRFTime"));
2382	if(!locRFTimeFactory->brun_was_called())
2383	{
2384	locRFTimeFactory->brun(locEventLoop, locEventLoop->GetJEvent().GetRunNumber());
2385	locRFTimeFactory->Set_brun_called();
2386	}
2387
2388	//Get time-based tracks: use best PID FOM
2389	//Note that these may not be the PIDs that were used in the fit!!!
2390	//e.g. for a DTrackTimeBased the proton hypothesis has the highest tracking FOM, so it is used in the fit, but the pi+ PID has the highest PID FOM
2391	vector<const DTrackTimeBased*> locTrackTimeBasedVector;
2392	for(size_t loc_i = 0; loc_i < locChargedTracks.size(); ++loc_i)
2393	{
2394	const DChargedTrackHypothesis* locChargedTrackHypothesis = locChargedTracks[loc_i]->Get_BestFOM();
2395	const DTrackTimeBased* locTrackTimeBased = NULL__null;
2396	locChargedTrackHypothesis->GetSingle(locTrackTimeBased);
2397	if(locTrackTimeBased != NULL__null)
2398	locTrackTimeBasedVector.push_back(locTrackTimeBased);
2399	}
2400
2401	//Event Vertex
2402	japp->RootWriteLock();
2403	{
2404	for(size_t loc_i = 0; loc_i < locChargedTracks.size(); ++loc_i)
2405	{
2406	const DChargedTrackHypothesis* locChargedTrackHypothesis = locChargedTracks[loc_i]->Get_BestFOM();
2407	double locPropagatedRFTime = locParticleID->Calc_PropagatedRFTime(locChargedTrackHypothesis, locEventRFBunch);
2408	double locShiftedRFTime = locRFTimeFactory->Step_TimeToNearInputTime(locPropagatedRFTime, locChargedTrackHypothesis->time());
2409	double locDeltaT = locShiftedRFTime - locChargedTrackHypothesis->time();
2410	dRFTrackDeltaT->Fill(locDeltaT);
2411	}
2412	dEventVertexZ_AllEvents->Fill(locVertex->dSpacetimeVertex.Z());
2413	dEventVertexYVsX_AllEvents->Fill(locVertex->dSpacetimeVertex.X(), locVertex->dSpacetimeVertex.Y());
2414	dEventVertexT_AllEvents->Fill(locVertex->dSpacetimeVertex.T());
2415
2416	if(locChargedTracks.size() >= 2)
2417	{
2418	dEventVertexZ_2OrMoreGoodTracks->Fill(locVertex->dSpacetimeVertex.Z());
2419	dEventVertexYVsX_2OrMoreGoodTracks->Fill(locVertex->dSpacetimeVertex.X(), locVertex->dSpacetimeVertex.Y());
2420	dEventVertexT_2OrMoreGoodTracks->Fill(locVertex->dSpacetimeVertex.T());
2421	}
2422	}
2423	japp->RootUnLock();
2424
2425	if(locVertex->dKinFitNDF == 0)
2426	return true; //kin fit not performed or didn't converge: no results to histogram
2427
2428	double locConfidenceLevel = TMath::Prob(locVertex->dKinFitChiSq, locVertex->dKinFitNDF);
2429
2430	japp->RootWriteLock();
2431	{
2432	dHist_KinFitConfidenceLevel->Fill(locConfidenceLevel);
2433
2434	//pulls
2435	if(locConfidenceLevel > dPullHistConfidenceLevelCut)
2436	{
2437	for(size_t loc_i = 0; loc_i < locTrackTimeBasedVector.size(); ++loc_i)
2438	{
2439	const DKinematicData* locKinematicData = static_cast<const DKinematicData*>(locTrackTimeBasedVector[loc_i]);
2440	Particle_t locPID = locKinematicData->PID();
2441	if(dHistMap_KinFitPulls.find(locPID) == dHistMap_KinFitPulls.end())
2442	continue; //PID not histogrammed
2443
2444	map<const JObject*, map<DKinFitPullType, double> >::const_iterator locParticleIterator = locVertex->dKinFitPulls.find(locKinematicData);
2445	if(locParticleIterator == locVertex->dKinFitPulls.end())
2446	continue;
2447
2448	const map<DKinFitPullType, double>& locPullMap = locParticleIterator->second;
2449	map<DKinFitPullType, double>::const_iterator locPullIterator = locPullMap.begin();
2450	for(; locPullIterator != locPullMap.end(); ++locPullIterator)
2451	dHistMap_KinFitPulls[locPID][locPullIterator->first]->Fill(locPullIterator->second);
2452	}
2453	}
2454	}
2455	japp->RootUnLock();
2456
2457	return true;
2458	}
2459
2460	void DHistogramAction_DetectedParticleKinematics::Initialize(JEventLoop* locEventLoop)
2461	{
2462	string locHistName, locHistTitle, locParticleName, locParticleROOTName;
2463	Particle_t locPID;
2464
2465	string locTrackSelectionTag = "NotATag", locShowerSelectionTag = "NotATag";
2466	if(gPARMS->Exists("COMBO:TRACK_SELECT_TAG"))
2467	gPARMS->GetParameter("COMBO:TRACK_SELECT_TAG", locTrackSelectionTag);
2468	if(gPARMS->Exists("COMBO:SHOWER_SELECT_TAG"))
2469	gPARMS->GetParameter("COMBO:SHOWER_SELECT_TAG", locShowerSelectionTag);
2470
2471	japp->RootWriteLock(); //ACQUIRE ROOT LOCK!!
2472	{
2473	//So: Default tag is "", User can set it to something else
2474	//In here, if tag is "", get from gparms, if not, leave it alone
2475	//If gparms value does not exist, set it to (and use) "PreSelect"
2476	if(dTrackSelectionTag == "NotATag")
2477	dTrackSelectionTag = (locTrackSelectionTag == "NotATag") ? "PreSelect" : locTrackSelectionTag;
2478	if(dShowerSelectionTag == "NotATag")
2479	dShowerSelectionTag = (locShowerSelectionTag == "NotATag") ? "PreSelect" : locShowerSelectionTag;
2480
2481	CreateAndChangeTo_ActionDirectory();
2482
2483	// Beam Particle
2484	locPID = Gamma;
2485	locParticleName = string("Beam_") + ParticleType(locPID);
2486	locParticleROOTName = ParticleName_ROOT(locPID);
2487	CreateAndChangeTo_Directory(locParticleName, locParticleName);
2488	locHistName = "Momentum";
2489	locHistTitle = string("Beam ") + locParticleROOTName + string(";p (GeV/c)");
2490	dBeamParticle_P = GetOrCreate_Histogram<TH1I>(locHistName, locHistTitle, dNumPBins, dMinP, dMaxP);
2491	gDirectory(TDirectory::CurrentDirectory())->cd("..");
2492
2493	//PID
2494	CreateAndChangeTo_Directory("PID", "PID");
2495	{
2496	//beta vs p
2497	locHistName = "BetaVsP_Q+";
2498	locHistTitle = "q^{+};p (GeV/c);#beta";
2499	dHistMap_QBetaVsP[1] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNumBetaBins, dMinBeta, dMaxBeta);
2500
2501	locHistName = "BetaVsP_Q-";
2502	locHistTitle = "q^{-};p (GeV/c);#beta";
2503	dHistMap_QBetaVsP[-1] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNumBetaBins, dMinBeta, dMaxBeta);
2504	}
2505	gDirectory(TDirectory::CurrentDirectory())->cd("..");
2506
2507	for(size_t loc_i = 0; loc_i < dFinalStatePIDs.size(); ++loc_i)
2508	{
2509	locPID = dFinalStatePIDs[loc_i];
2510	locParticleName = ParticleType(locPID);
2511	locParticleROOTName = ParticleName_ROOT(locPID);
2512	CreateAndChangeTo_Directory(locParticleName, locParticleName);
2513
2514	// Momentum
2515	locHistName = "Momentum";
2516	locHistTitle = locParticleROOTName + string(";p (GeV/c)");
2517	dHistMap_P[locPID] = GetOrCreate_Histogram<TH1I>(locHistName, locHistTitle, dNumPBins, dMinP, dMaxP);
2518
2519	// Theta
2520	locHistName = "Theta";
2521	locHistTitle = locParticleROOTName + string(";#theta#circ");
2522	dHistMap_Theta[locPID] = GetOrCreate_Histogram<TH1I>(locHistName, locHistTitle, dNumThetaBins, dMinTheta, dMaxTheta);
2523
2524	// Phi
2525	locHistName = "Phi";
2526	locHistTitle = locParticleROOTName + string(";#phi#circ");
2527	dHistMap_Phi[locPID] = GetOrCreate_Histogram<TH1I>(locHistName, locHistTitle, dNumPhiBins, dMinPhi, dMaxPhi);
2528
2529	// P Vs Theta
2530	locHistName = "PVsTheta";
2531	locHistTitle = locParticleROOTName + string(";#theta#circ;p (GeV/c)");
2532	dHistMap_PVsTheta[locPID] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DThetaBins, dMinTheta, dMaxTheta, dNum2DPBins, dMinP, dMaxP);
2533
2534	// Phi Vs Theta
2535	locHistName = "PhiVsTheta";
2536	locHistTitle = locParticleROOTName + string(";#theta#circ;#phi#circ");
2537	dHistMap_PhiVsTheta[locPID] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DThetaBins, dMinTheta, dMaxTheta, dNum2DPhiBins, dMinPhi, dMaxPhi);
2538
2539	//beta vs p
2540	locHistName = "BetaVsP";
2541	locHistTitle = locParticleROOTName + string(";p (GeV/c);#beta");
2542	dHistMap_BetaVsP[locPID] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNumBetaBins, dMinBeta, dMaxBeta);
2543
2544	//delta-beta vs p
2545	locHistName = "DeltaBetaVsP";
2546	locHistTitle = locParticleROOTName + string(";p (GeV/c);#Delta#beta");
2547	dHistMap_DeltaBetaVsP[locPID] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNum2DPBins, dMinP, dMaxP, dNum2DDeltaBetaBins, dMinDeltaBeta, dMaxDeltaBeta);
2548
2549	// Vertex-Z
2550	locHistName = "VertexZ";
2551	locHistTitle = locParticleROOTName + string(";Vertex-Z (cm)");
2552	dHistMap_VertexZ[locPID] = GetOrCreate_Histogram<TH1I>(locHistName, locHistTitle, dNumVertexZBins, dMinVertexZ, dMaxVertexZ);
2553
2554	// Vertex-Y Vs Vertex-X
2555	locHistName = "VertexYVsX";
2556	locHistTitle = locParticleROOTName + string(";Vertex-X (cm);Vertex-Y (cm)");
2557	dHistMap_VertexYVsX[locPID] = GetOrCreate_Histogram<TH2I>(locHistName, locHistTitle, dNumVertexXYBins, dMinVertexXY, dMaxVertexXY, dNumVertexXYBins, dMinVertexXY, dMaxVertexXY);
2558
2559	// Vertex-T
2560	locHistName = "VertexT";
2561	locHistTitle = locParticleROOTName + string(";Vertex-T (ns)");
2562	dHistMap_VertexT[locPID] = GetOrCreate_Histogram<TH1I>(locHistName, locHistTitle, dNumTBins, dMinT, dMaxT);
2563
2564	gDirectory(TDirectory::CurrentDirectory())->cd("..");
2565	}
2566
2567	//Return to the base directory
2568	ChangeTo_BaseDirectory();
2569	}
2570	japp->RootUnLock(); //RELEASE ROOT LOCK!!
2571	}
2572
2573	bool DHistogramAction_DetectedParticleKinematics::Perform_Action(JEventLoop* locEventLoop, const DParticleCombo* locParticleCombo)
2574	{
2575	if(Get_NumPreviousParticleCombos() != 0)
2576	return true; //else double-counting!
2577
2578	vector<const DBeamPhoton*> locBeamPhotons;
2579	locEventLoop->Get(locBeamPhotons);
2580	japp->RootWriteLock();
2581	{
2582	for(size_t loc_i = 0; loc_i < locBeamPhotons.size(); ++loc_i)
2583	dBeamParticle_P->Fill(locBeamPhotons[loc_i]->energy());
2584	}
2585	japp->RootUnLock();
2586
2587	vector<const DChargedTrack*> locPreSelectChargedTracks;
2588	locEventLoop->Get(locPreSelectChargedTracks, dTrackSelectionTag.c_str());
2589
2590	for(size_t loc_i = 0; loc_i < locPreSelectChargedTracks.size(); ++loc_i)
2591	{
2592	const DChargedTrackHypothesis* locChargedTrackHypothesis = locPreSelectChargedTracks[loc_i]->Get_BestTrackingFOM();
2593	int locCharge = ParticleCharge(locChargedTrackHypothesis->PID());
2594
2595	DVector3 locMomentum = locChargedTrackHypothesis->momentum();
2596	double locP = locMomentum.Mag();
2597	double locBeta_Timing = locChargedTrackHypothesis->measuredBeta();
2598
2599	if(dHistMap_QBetaVsP.find(locCharge) == dHistMap_QBetaVsP.end())
2600	continue;
2601
2602	japp->RootWriteLock();
2603	{
2604	//Extremely inefficient, I know ...
2605	dHistMap_QBetaVsP[locCharge]->Fill(locP, locBeta_Timing);
2606	}
2607	japp->RootUnLock();
2608	}
2609
2610	for(size_t loc_i = 0; loc_i < locPreSelectChargedTracks.size(); ++loc_i)
2611	{
2612	const DChargedTrackHypothesis* locChargedTrackHypothesis = locPreSelectChargedTracks[loc_i]->Get_BestFOM();
2613
2614	DVector3 locMomentum = locChargedTrackHypothesis->momentum();
2615	double locPhi = locMomentum.Phi()*180.0/TMath::Pi();
2616	double locTheta = locMomentum.Theta()*180.0/TMath::Pi();
2617	double locP = locMomentum.Mag();
2618	double locBeta_Timing = locChargedTrackHypothesis->measuredBeta();
2619	double locDeltaBeta = locChargedTrackHypothesis->deltaBeta();
2620
2621	Particle_t locPID = (locChargedTrackHypothesis->dFOM < dMinPIDFOM) ? Unknown : locChargedTrackHypothesis->PID();
2622	if(dHistMap_P.find(locPID) == dHistMap_P.end())
2623	continue; //not interested in histogramming
2624
2625	japp->RootWriteLock();
2626	{
2627	dHistMap_P[locPID]->Fill(locP);
2628	dHistMap_Phi[locPID]->Fill(locPhi);
2629	dHistMap_Theta[locPID]->Fill(locTheta);
2630	dHistMap_PVsTheta[locPID]->Fill(locTheta, locP);
2631	dHistMap_PhiVsTheta[locPID]->Fill(locTheta, locPhi);
2632	dHistMap_BetaVsP[locPID]->Fill(locP, locBeta_Timing);
2633	dHistMap_DeltaBetaVsP[locPID]->Fill(locP, locDeltaBeta);
2634	dHistMap_VertexZ[locPID]->Fill(locChargedTrackHypothesis->position().Z());
2635	dHistMap_VertexYVsX[locPID]->Fill(locChargedTrackHypothesis->position().X(), locChargedTrackHypothesis->position().Y());
2636	dHistMap_VertexT[locPID]->Fill(locChargedTrackHypothesis->time());
2637	}
2638	japp->RootUnLock();
2639	}
2640
2641	vector<const DNeutralParticle*> locNeutralParticles;
2642	locEventLoop->Get(locNeutralParticles);
2643
2644	vector<const DNeutralShower*> locNeutralShowers;
2645	locEventLoop->Get(locNeutralShowers, dShowerSelectionTag.c_str());
2646
2647	for(size_t loc_i = 0; loc_i < locNeutralShowers.size(); ++loc_i)
2648	{
2649	const DNeutralParticle* locNeutralParticle = NULL__null;
2650	for(size_t loc_j = 0; loc_j < locNeutralParticles.size(); ++loc_j)
2651	{
2652	if(locNeutralParticles[loc_j]->dNeutralShower != locNeutralShowers[loc_i])
2653	continue;
2654	locNeutralParticle = locNeutralParticles[loc_j];
2655	break;
2656	}
2657	if(locNeutralParticle == NULL__null)
2658	continue;
2659	const DNeutralParticleHypothesis* locNeutralParticleHypothesis = locNeutralParticle->Get_Hypothesis(Gamma);
2660	if(locNeutralParticleHypothesis->dFOM < dMinPIDFOM)
2661	continue;
2662
2663	Particle_t locPID = locNeutralParticleHypothesis->PID();
2664
2665	if(dHistMap_P.find(locPID) == dHistMap_P.end())
2666	continue; //e.g. a decaying particle, or not interested in histogramming
2667
2668	DVector3 locMomentum = locNeutralParticleHypothesis->momentum();
2669	double locPhi = locMomentum.Phi()*180.0/TMath::Pi();
2670	double locTheta = locMomentum.Theta()*180.0/TMath::Pi();
2671	double locP = locMomentum.Mag();
2672
2673	double locBeta_Timing = locNeutralParticleHypothesis->measuredBeta();
2674	double locDeltaBeta = locNeutralParticleHypothesis->deltaBeta();
2675
2676	japp->RootWriteLock();
2677	{
2678	dHistMap_P[locPID]->Fill(locP);
2679	dHistMap_Phi[locPID]->Fill(locPhi);
2680	dHistMap_Theta[locPID]->Fill(locTheta);
2681	dHistMap_PVsTheta[locPID]->Fill(locTheta, locP);
2682	dHistMap_PhiVsTheta[locPID]->Fill(locTheta, locPhi);
2683	dHistMap_BetaVsP[locPID]->Fill(locP, locBeta_Timing);
2684	dHistMap_DeltaBetaVsP[locPID]->Fill(locP, locDeltaBeta);
2685	dHistMap_VertexZ[locPID]->Fill(locNeutralParticleHypothesis->position().Z());
2686	dHistMap_VertexYVsX[locPID]->Fill(locNeutralParticleHypothesis->position().X(), locNeutralParticleHypothesis->position().Y());
2687	dHistMap_VertexT[locPID]->Fill(locNeutralParticleHypothesis->time());
2688	}
2689	japp->RootUnLock();
2690	}
2691	return true;
2692	}
2693
2694	void DHistogramAction_NumReconstructedObjects::Initialize(JEventLoop* locEventLoop)
2695	{
2696	string locHistName;
2697
2698	bool locIsRESTEvent = (string(locEventLoop->GetJEvent().GetJEventSource()->className()) == string("DEventSourceREST"));
2699
2700	//CREATE THE HISTOGRAMS
2701	japp->RootWriteLock(); //ACQUIRE ROOT LOCK!!
2702	{
2703	CreateAndChangeTo_ActionDirectory();
2704
2705	//2D Summary
2706	locHistName = "NumHighLevelObjects";
2707	if(gDirectory(TDirectory::CurrentDirectory())->Get(locHistName.c_str()) != NULL__null) //already created by another thread, or directory name is duplicate (e.g. two identical steps)
2708	dHist_NumHighLevelObjects = static_cast<TH2D*>(gDirectory(TDirectory::CurrentDirectory())->Get(locHistName.c_str()));
2709	else
2710	{
2711	dHist_NumHighLevelObjects = new TH2D(locHistName.c_str(), ";;# Objects / Event", 13, 0.5, 13.5, dMaxNumObjects + 1, -0.5, (float)dMaxNumObjects + 0.5);
2712	dHist_NumHighLevelObjects->GetXaxis()->SetBinLabel(1, "DRFTime");
2713	dHist_NumHighLevelObjects->GetXaxis()->SetBinLabel(2, "DSCHit");
2714	dHist_NumHighLevelObjects->GetXaxis()->SetBinLabel(3, "DTOFPoint");
2715	dHist_NumHighLevelObjects->GetXaxis()->SetBinLabel(4, "DBCALShower");
2716	dHist_NumHighLevelObjects->GetXaxis()->SetBinLabel(5, "DFCALShower");
2717	dHist_NumHighLevelObjects->GetXaxis()->SetBinLabel(6, "DTimeBasedTrack");
2718	dHist_NumHighLevelObjects->GetXaxis()->SetBinLabel(7, "TrackSCMatches");
2719	dHist_NumHighLevelObjects->GetXaxis()->SetBinLabel(8, "TrackTOFMatches");
2720	dHist_NumHighLevelObjects->GetXaxis()->SetBinLabel(9, "TrackBCALMatches");
2721	dHist_NumHighLevelObjects->GetXaxis()->SetBinLabel(10, "TrackFCALMatches");
2722	dHist_NumHighLevelObjects->GetXaxis()->SetBinLabel(11, "DBeamPhoton");
2723	dHist_NumHighLevelObjects->GetXaxis()->SetBinLabel(12, "DChargedTrack");
2724	dHist_NumHighLevelObjects->GetXaxis()->SetBinLabel(13, "DNeutralShower");
2725	}
2726
2727	//Charged
2728	locHistName = "NumChargedTracks";
2729	dHist_NumChargedTracks = GetOrCreate_Histogram<TH1D>(locHistName, ";# DChargedTrack", dMaxNumObjects + 1, -0.5, (float)dMaxNumObjects + 0.5);
2730	locHistName = "NumPosChargedTracks";
2731	dHist_NumPosChargedTracks = GetOrCreate_Histogram<TH1D>(locHistName, ";# #it{q}^{+} DChargedTrack", dMaxNumObjects + 1, -0.5, (float)dMaxNumObjects + 0.5);
2732	locHistName = "NumNegChargedTracks";
2733	dHist_NumNegChargedTracks = GetOrCreate_Histogram<TH1D>(locHistName, ";# #it{q}^{-} DChargedTrack", dMaxNumObjects + 1, -0.5, (float)dMaxNumObjects + 0.5);
2734
2735	//TBT
2736	locHistName = "NumTimeBasedTracks";
2737	dHist_NumTimeBasedTracks = GetOrCreate_Histogram<TH1D>(locHistName, ";# Time-Based Tracks", dMaxNumObjects + 1, -0.5, (float)dMaxNumObjects + 0.5);
2738	locHistName = "NumPosTimeBasedTracks";
2739	dHist_NumPosTimeBasedTracks = GetOrCreate_Histogram<TH1D>(locHistName, ";# #it{q}^{+} Time-Based Tracks", dMaxNumObjects + 1, -0.5, (float)dMaxNumObjects + 0.5);
2740	locHistName = "NumNegTimeBasedTracks";
2741	dHist_NumNegTimeBasedTracks = GetOrCreate_Histogram<TH1D>(locHistName, ";# #it{q}^{-} Time-Based Tracks", dMaxNumObjects + 1, -0.5, (float)dMaxNumObjects + 0.5);
2742
2743	if(!locIsRESTEvent)
2744	{
2745	//WBT
2746	locHistName = "NumWireBasedTracks";
2747	dHist_NumWireBasedTracks = GetOrCreate_Histogram<TH1D>(locHistName, ";# Wire-Based Tracks", dMaxNumObjects + 1, -0.5, (float)dMaxNumObjects + 0.5);
2748	locHistName = "NumPosWireBasedTracks";
2749	dHist_NumPosWireBasedTracks = GetOrCreate_Histogram<TH1D>(locHistName, ";# #it{q}^{-} Wire-Based Tracks", dMaxNumObjects + 1, -0.5, (float)dMaxNumObjects + 0.5);
2750	locHistName = "NumNegWireBasedTracks";
2751	dHist_NumNegWireBasedTracks = GetOrCreate_Histogram<TH1D>(locHistName, ";# #it{q}^{-} Wire-Based Tracks", dMaxNumObjects + 1, -0.5, (float)dMaxNumObjects + 0.5);
2752
2753	//Track Candidates
2754	locHistName = "NumTrackCandidates";
2755	dHist_NumTrackCandidates = GetOrCreate_Histogram<TH1D>(locHistName, ";# Track Candidates", dMaxNumObjects + 1, -0.5, (float)dMaxNumObjects + 0.5);
2756	locHistName = "NumPosTrackCandidates";
2757	dHist_NumPosTrackCandidates = GetOrCreate_Histogram<TH1D>(locHistName, ";# #it{q}^{+} Track Candidates", dMaxNumObjects + 1, -0.5, (float)dMaxNumObjects + 0.5);
2758	locHistName = "NumNegTrackCandidates";
2759	dHist_NumNegTrackCandidates = GetOrCreate_Histogram<TH1D>(locHistName, ";# #it{q}^{-} Track Candidates", dMaxNumObjects + 1, -0.5, (float)dMaxNumObjects + 0.5);
2760
2761	//CDC Track Candidates
2762	locHistName = "NumPosTrackCandidates_CDC";
2763	dHist_NumPosTrackCandidates_CDC = GetOrCreate_Histogram<TH1D>(locHistName, ";# #it{q}^{+} CDC Track Candidates", dMaxNumObjects + 1, -0.5, (float)dMaxNumObjects + 0.5);
2764	locHistName = "NumNegTrackCandidates_CDC";
2765	dHist_NumNegTrackCandidates_CDC = GetOrCreate_Histogram<TH1D>(locHistName, ";# #it{q}^{-} CDC Track Candidates", dMaxNumObjects + 1, -0.5, (float)dMaxNumObjects + 0.5);
2766
2767	//FDC Track Candidates
2768	locHistName = "NumPosTrackCandidates_FDC";
2769	dHist_NumPosTrackCandidates_FDC = GetOrCreate_Histogram<TH1D>(locHistName, ";# #it{q}^{+} FDC Track Candidates", dMaxNumObjects + 1, -0.5, (float)dMaxNumObjects + 0.5);
2770	locHistName = "NumNegTrackCandidates_FDC";
2771	dHist_NumNegTrackCandidates_FDC = GetOrCreate_Histogram<TH1D>(locHistName, ";# #it{q}^{-} FDC Track Candidates", dMaxNumObjects + 1, -0.5, (float)dMaxNumObjects + 0.5);
2772	}
2773
2774	//Beam Photons
2775	locHistName = "NumBeamPhotons";
2776	dHist_NumBeamPhotons = GetOrCreate_Histogram<TH1D>(locHistName, ";# DBeamPhoton", dMaxNumObjects + 1, -0.5, (float)dMaxNumObjects + 0.5);
2777
2778	//Showers / Neutrals / TOF / SC
2779	locHistName = "NumFCALShowers";
2780	dHist_NumFCALShowers = GetOrCreate_Histogram<TH1D>(locHistName, ";# DFCALShower", dMaxNumObjects + 1, -0.5, (float)dMaxNumObjects + 0.5);
2781	locHistName = "NumBCALShowers";
2782	dHist_NumBCALShowers = GetOrCreate_Histogram<TH1D>(locHistName, ";# DBCALShower", dMaxNumObjects + 1, -0.5, (float)dMaxNumObjects + 0.5);
2783	locHistName = "NumNeutralShowers";
2784	dHist_NumNeutralShowers = GetOrCreate_Histogram<TH1D>(locHistName, ";# DNeutralShower", dMaxNumObjects + 1, -0.5, (float)dMaxNumObjects + 0.5);
2785	locHistName = "NumTOFPoints";
2786	dHist_NumTOFPoints = GetOrCreate_Histogram<TH1D>(locHistName, ";# DTOFPoint", dMaxNumObjects + 1, -0.5, (float)dMaxNumObjects + 0.5);
2787	locHistName = "NumSCHits";
2788	dHist_NumSCHits = GetOrCreate_Histogram<TH1D>(locHistName, ";# DSCHit", dMaxNumObjects + 1, -0.5, (float)dMaxNumObjects + 0.5);
2789
2790	if(!locIsRESTEvent)
2791	{
2792	locHistName = "NumTAGMHits";
2793	dHist_NumTAGMHits = GetOrCreate_Histogram<TH1D>(locHistName, ";# DTAGMHit", dMaxNumObjects + 1, -0.5, (float)dMaxNumObjects + 0.5);
2794	locHistName = "NumTAGHHits";
2795	dHist_NumTAGHHits = GetOrCreate_Histogram<TH1D>(locHistName, ";# DTAGHHit", dMaxNumObjects + 1, -0.5, (float)dMaxNumObjects + 0.5);
2796	}
2797
2798	//Matches
2799	locHistName = "NumTrackBCALMatches";
2800	dHist_NumTrackBCALMatches = GetOrCreate_Histogram<TH1D>(locHistName, ";# Track-BCAL Matches", dMaxNumMatchObjects + 1, -0.5, (float)dMaxNumMatchObjects + 0.5);
2801	locHistName = "NumTrackFCALMatches";
2802	dHist_NumTrackFCALMatches = GetOrCreate_Histogram<TH1D>(locHistName, ";# Track-FCAL Matches", dMaxNumMatchObjects + 1, -0.5, (float)dMaxNumMatchObjects + 0.5);
2803	locHistName = "NumTrackTOFMatches";
2804	dHist_NumTrackTOFMatches = GetOrCreate_Histogram<TH1D>(locHistName, ";# Track-TOF Matches", dMaxNumMatchObjects + 1, -0.5, (float)dMaxNumMatchObjects + 0.5);
2805	locHistName = "NumTrackSCMatches";
2806	dHist_NumTrackSCMatches = GetOrCreate_Histogram<TH1D>(locHistName, ";# Track-SC Matches", dMaxNumMatchObjects + 1, -0.5, (float)dMaxNumMatchObjects + 0.5);
2807
2808	if(!locIsRESTEvent)
2809	{
2810	//Hits
2811	locHistName = "NumCDCHits";
2812	dHist_NumCDCHits = GetOrCreate_Histogram<TH1I>(locHistName, ";# DCDCHit", dMaxNumCDCHits + 1, -0.5, (float)dMaxNumCDCHits + 0.5);
2813	locHistName = "NumFDCWireHits";
2814	dHist_NumFDCWireHits = GetOrCreate_Histogram<TH1I>(locHistName, ";# Wire DFDCHit", dMaxNumFDCHits/2 + 1, -0.5, (float)dMaxNumFDCHits - 0.5 + 2.0);
2815	locHistName = "NumFDCCathodeHits";
2816	dHist_NumFDCCathodeHits = GetOrCreate_Histogram<TH1I>(locHistName, ";# Cathode DFDCHit", dMaxNumFDCHits/2 + 1, -0.5, (float)dMaxNumFDCHits - 0.5 + 2.0);
2817
2818	locHistName = "NumTOFHits";
2819	dHist_NumTOFHits = GetOrCreate_Histogram<TH1I>(locHistName, ";# DTOFHit", dMaxNumTOFCalorimeterHits + 1, -0.5, (float)dMaxNumTOFCalorimeterHits + 0.5);
2820	locHistName = "NumBCALHits";
2821	dHist_NumBCALHits = GetOrCreate_Histogram<TH1I>(locHistName, ";# DBCALHit", dMaxNumTOFCalorimeterHits + 1, -0.5, (float)dMaxNumTOFCalorimeterHits + 0.5);
2822	locHistName = "NumFCALHits";
2823	dHist_NumFCALHits = GetOrCreate_Histogram<TH1I>(locHistName, ";# DFCALHit", dMaxNumTOFCalorimeterHits + 1, -0.5, (float)dMaxNumTOFCalorimeterHits + 0.5);
2824
2825	locHistName = "NumRFSignals";
2826	dHist_NumRFSignals = GetOrCreate_Histogram<TH1I>(locHistName, ";# DRFDigiTime + # DRFTDCDigiTime", dMaxNumObjects + 1, -0.5, (float)dMaxNumObjects + 0.5);
2827	}
2828
2829	//Return to the base directory
2830	ChangeTo_BaseDirectory();
2831	}
2832	japp->RootUnLock(); //RELEASE ROOT LOCK!!
2833	}
2834
2835	bool DHistogramAction_NumReconstructedObjects::Perform_Action(JEventLoop* locEventLoop, const DParticleCombo* locParticleCombo)
2836	{
2837	if(Get_NumPreviousParticleCombos() != 0)
2838	return true; //else double-counting!
2839
2840	bool locIsRESTEvent = (string(locEventLoop->GetJEvent().GetJEventSource()->className()) == string("DEventSourceREST"));
2841
2842	vector<const DTrackTimeBased*> locTrackTimeBasedVector;
2843	locEventLoop->Get(locTrackTimeBasedVector);
2844
2845	vector<const DBeamPhoton*> locBeamPhotons;
2846	locEventLoop->Get(locBeamPhotons);
2847
2848	vector<const DFCALShower*> locFCALShowers;
2849	locEventLoop->Get(locFCALShowers);
2850
2851	vector<const DChargedTrack*> locChargedTracks;
2852	locEventLoop->Get(locChargedTracks);
2853
2854	vector<const DBCALShower*> locBCALShowers;
2855	locEventLoop->Get(locBCALShowers);
2856
2857	vector<const DNeutralShower*> locNeutralShowers;
2858	locEventLoop->Get(locNeutralShowers);
2859
2860	vector<const DTOFPoint*> locTOFPoints;
2861	locEventLoop->Get(locTOFPoints);
2862
2863	vector<const DSCHit*> locSCHits;
2864	locEventLoop->Get(locSCHits);
2865
2866	vector<const DRFTime*> locRFTimes;
2867	locEventLoop->Get(locRFTimes);
2868
2869	const DDetectorMatches* locDetectorMatches = NULL__null;
2870	locEventLoop->GetSingle(locDetectorMatches);
2871
2872	//if not REST
2873	vector<const DTrackWireBased*> locTrackWireBasedVector;
2874	vector<const DTrackCandidate*> locTrackCandidates;
2875	vector<const DTrackCandidate*> locTrackCandidates_CDC;
2876	vector<const DTrackCandidate*> locTrackCandidates_FDC;
2877	vector<const DCDCHit*> locCDCHits;
2878	vector<const DFDCHit*> locFDCHits;
2879	vector<const DTOFHit*> locTOFHits;
2880	vector<const DBCALHit*> locBCALHits;
2881	vector<const DFCALHit*> locFCALHits;
2882	vector<const DTAGMHit*> locTAGMHits;
2883	vector<const DTAGHHit*> locTAGHHits;
2884	vector<const DRFDigiTime*> locRFDigiTimes;
2885	vector<const DRFTDCDigiTime*> locRFTDCDigiTimes;
2886
2887	size_t locNumFDCWireHits = 0, locNumFDCCathodeHits = 0;
2888	if(!locIsRESTEvent)
2889	{
2890	locEventLoop->Get(locTrackWireBasedVector);
2891	locEventLoop->Get(locTrackCandidates);
2892	locEventLoop->Get(locTrackCandidates_CDC, "CDC");
2893	locEventLoop->Get(locTrackCandidates_FDC, "FDCCathodes");
2894	locEventLoop->Get(locCDCHits);
2895	locEventLoop->Get(locFDCHits);
2896	locEventLoop->Get(locTOFHits);
2897	locEventLoop->Get(locBCALHits);
2898	locEventLoop->Get(locFCALHits);
2899	locEventLoop->Get(locTAGHHits);
2900	locEventLoop->Get(locTAGMHits);
2901	locEventLoop->Get(locRFDigiTimes);
2902	locEventLoop->Get(locRFTDCDigiTimes);
2903
2904	for(size_t loc_i = 0; loc_i < locFDCHits.size(); ++loc_i)
2905	{
2906	if(locFDCHits[loc_i]->type == DFDCHit::AnodeWire)
2907	++locNumFDCWireHits;
2908	else
2909	++locNumFDCCathodeHits;
2910	}
2911	}
2912
2913	japp->RootWriteLock(); //ACQUIRE ROOT LOCK!!
2914	{
2915	//# High-Level Objects
2916	dHist_NumHighLevelObjects->Fill(1, (Double_t)locRFTimes.size());
2917	dHist_NumHighLevelObjects->Fill(2, (Double_t)locSCHits.size());
2918	dHist_NumHighLevelObjects->Fill(3, (Double_t)locTOFPoints.size());
2919	dHist_NumHighLevelObjects->Fill(4, (Double_t)locBCALShowers.size());
2920	dHist_NumHighLevelObjects->Fill(5, (Double_t)locFCALShowers.size());
2921	dHist_NumHighLevelObjects->Fill(6, (Double_t)locTrackTimeBasedVector.size());
2922	dHist_NumHighLevelObjects->Fill(7, (Double_t)locDetectorMatches->Get_NumTrackSCMatches());
2923	dHist_NumHighLevelObjects->Fill(8, (Double_t)locDetectorMatches->Get_NumTrackTOFMatches());
2924	dHist_NumHighLevelObjects->Fill(9, (Double_t)locDetectorMatches->Get_NumTrackBCALMatches());
2925	dHist_NumHighLevelObjects->Fill(10, (Double_t)locDetectorMatches->Get_NumTrackFCALMatches());
2926	dHist_NumHighLevelObjects->Fill(11, (Double_t)locBeamPhotons.size());
2927	dHist_NumHighLevelObjects->Fill(12, (Double_t)locChargedTracks.size());
2928	dHist_NumHighLevelObjects->Fill(13, (Double_t)locNeutralShowers.size());
2929
2930	//Charged
2931	unsigned int locNumPos = 0, locNumNeg = 0;
2932	for(size_t loc_i = 0; loc_i < locChargedTracks.size(); ++loc_i)
2933	{
2934	if(ParticleCharge(locChargedTracks[loc_i]->Get_BestFOM()->PID()) > 0)
2935	++locNumPos;
2936	else
2937	++locNumNeg;
2938	}
2939	dHist_NumChargedTracks->Fill(locChargedTracks.size());
2940	dHist_NumPosChargedTracks->Fill(locNumPos);
2941	dHist_NumNegChargedTracks->Fill(locNumNeg);
2942
2943	//TBT
2944	locNumPos = 0; locNumNeg = 0;
2945	for(size_t loc_i = 0; loc_i < locTrackTimeBasedVector.size(); ++loc_i)
2946	{
2947	if(ParticleCharge(locTrackTimeBasedVector[loc_i]->PID()) > 0)
2948	++locNumPos;
2949	else
2950	++locNumNeg;
2951	}
2952	dHist_NumTimeBasedTracks->Fill(locTrackTimeBasedVector.size());
2953	dHist_NumPosTimeBasedTracks->Fill(locNumPos);
2954	dHist_NumNegTimeBasedTracks->Fill(locNumNeg);
2955
2956	if(!locIsRESTEvent)
2957	{
2958	//WBT
2959	locNumPos = 0; locNumNeg = 0;
2960	for(size_t loc_i = 0; loc_i < locTrackWireBasedVector.size(); ++loc_i)
2961	{
2962	if(ParticleCharge(locTrackWireBasedVector[loc_i]->PID()) > 0)
2963	++locNumPos;
2964	else
2965	++locNumNeg;
2966	}
2967	dHist_NumWireBasedTracks->Fill(locTrackWireBasedVector.size());
2968	dHist_NumPosWireBasedTracks->Fill(locNumPos);
2969	dHist_NumNegWireBasedTracks->Fill(locNumNeg);
2970
2971	//Candidates
2972	locNumPos = 0; locNumNeg = 0;
2973	for(size_t loc_i = 0; loc_i < locTrackCandidates.size(); ++loc_i)
2974	{
2975	if(locTrackCandidates[loc_i]->charge() > 0.0)
2976	++locNumPos;
2977	else
2978	++locNumNeg;
2979	}
2980	dHist_NumTrackCandidates->Fill(locTrackCandidates.size());
2981	dHist_NumPosTrackCandidates->Fill(locNumPos);
2982	dHist_NumNegTrackCandidates->Fill(locNumNeg);
2983
2984	//CDC Candidates
2985	locNumPos = 0; locNumNeg = 0;
2986	for(size_t loc_i = 0; loc_i < locTrackCandidates_CDC.size(); ++loc_i)
2987	{
2988	if(locTrackCandidates_CDC[loc_i]->charge() > 0.0)
2989	++locNumPos;
2990	else
2991	++locNumNeg;
2992	}
2993	dHist_NumPosTrackCandidates_CDC->Fill(locNumPos);
2994	dHist_NumNegTrackCandidates_CDC->Fill(locNumNeg);
2995
2996	//FDC Candidates
2997	locNumPos = 0; locNumNeg = 0;
2998	for(size_t loc_i = 0; loc_i < locTrackCandidates_FDC.size(); ++loc_i)
2999	{
3000	if(locTrackCandidates_FDC[loc_i]->charge() > 0.0)
3001	++locNumPos;
3002	else
3003	++locNumNeg;
3004	}
3005	dHist_NumPosTrackCandidates_FDC->Fill(locNumPos);
3006	dHist_NumNegTrackCandidates_FDC->Fill(locNumNeg);
3007	}
3008
3009	//Beam Photons
3010	dHist_NumBeamPhotons->Fill((Double_t)locBeamPhotons.size());
3011
3012	//Showers
3013	dHist_NumFCALShowers->Fill((Double_t)locFCALShowers.size());
3014	dHist_NumBCALShowers->Fill((Double_t)locBCALShowers.size());
3015	dHist_NumNeutralShowers->Fill((Double_t)locNeutralShowers.size());
3016
3017	//TOF & SC
3018	dHist_NumTOFPoints->Fill((Double_t)locTOFPoints.size());
3019	dHist_NumSCHits->Fill((Double_t)locSCHits.size());
3020
3021	//TAGGER
3022	if(!locIsRESTEvent)
3023	{
3024	dHist_NumTAGMHits->Fill((Double_t)locTAGMHits.size());
3025	dHist_NumTAGHHits->Fill((Double_t)locTAGHHits.size());
3026	}
3027
3028	//Matches
3029	dHist_NumTrackBCALMatches->Fill((Double_t)locDetectorMatches->Get_NumTrackBCALMatches());
3030	dHist_NumTrackFCALMatches->Fill((Double_t)locDetectorMatches->Get_NumTrackFCALMatches());
3031	dHist_NumTrackTOFMatches->Fill((Double_t)locDetectorMatches->Get_NumTrackTOFMatches());
3032	dHist_NumTrackSCMatches->Fill((Double_t)locDetectorMatches->Get_NumTrackSCMatches());
3033
3034	//Hits
3035	if(!locIsRESTEvent)
3036	{
3037	dHist_NumCDCHits->Fill((Double_t)locCDCHits.size());
3038	dHist_NumFDCWireHits->Fill((Double_t)locNumFDCWireHits);
3039	dHist_NumFDCCathodeHits->Fill((Double_t)locNumFDCCathodeHits);
3040	dHist_NumTOFHits->Fill((Double_t)locTOFHits.size());
3041	dHist_NumBCALHits->Fill((Double_t)locBCALHits.size());
3042	dHist_NumFCALHits->Fill((Double_t)locFCALHits.size());
3043	dHist_NumRFSignals->Fill((Double_t)(locRFDigiTimes.size() + locRFTDCDigiTimes.size()));
3044	}
3045	}
3046	japp->RootUnLock(); //RELEASE ROOT LOCK!!
3047
3048	return true;
3049	}
3050
3051	void DHistogramAction_TrackMultiplicity::Initialize(JEventLoop* locEventLoop)
3052	{
3053	string locTrackSelectionTag = "NotATag", locShowerSelectionTag = "NotATag";
3054	if(gPARMS->Exists("COMBO:TRACK_SELECT_TAG"))
3055	gPARMS->GetParameter("COMBO:TRACK_SELECT_TAG", locTrackSelectionTag);
3056	if(gPARMS->Exists("COMBO:SHOWER_SELECT_TAG"))
3057	gPARMS->GetParameter("COMBO:SHOWER_SELECT_TAG", locShowerSelectionTag);
3058
3059	japp->RootWriteLock(); //ACQUIRE ROOT LOCK!!
3060	{
3061	//So: Default tag is "", User can set it to something else
3062	//In here, if tag is "", get from gparms, if not, leave it alone
3063	//If gparms value does not exist, set it to (and use) "PreSelect"
3064	if(dTrackSelectionTag == "NotATag")
3065	dTrackSelectionTag = (locTrackSelectionTag == "NotATag") ? "PreSelect" : locTrackSelectionTag;
3066	if(dShowerSelectionTag == "NotATag")
3067	dShowerSelectionTag = (locShowerSelectionTag == "NotATag") ? "PreSelect" : locShowerSelectionTag;
3068
3069	CreateAndChangeTo_ActionDirectory();
3070
3071	string locHistName("NumReconstructedParticles");
3072	if(gDirectory(TDirectory::CurrentDirectory())->Get(locHistName.c_str()) != NULL__null) //already created by another thread, or directory name is duplicate (e.g. two identical steps)
3073	dHist_NumReconstructedParticles = static_cast<TH2D*>(gDirectory(TDirectory::CurrentDirectory())->Get(locHistName.c_str()));
3074	else
3075	{
3076	dHist_NumReconstructedParticles = new TH2D("NumReconstructedParticles", ";Particle Type;Num Particles / Event", 5 + dFinalStatePIDs.size(), -0.5, 4.5 + dFinalStatePIDs.size(), dMaxNumTracks + 1, -0.5, (float)dMaxNumTracks + 0.5);
3077	dHist_NumReconstructedParticles->GetXaxis()->SetBinLabel(1, "# Total");
3078	dHist_NumReconstructedParticles->GetXaxis()->SetBinLabel(2, "# q != 0");
3079	dHist_NumReconstructedParticles->GetXaxis()->SetBinLabel(3, "# q = 0");
3080	dHist_NumReconstructedParticles->GetXaxis()->SetBinLabel(4, "# q = +");
3081	dHist_NumReconstructedParticles->GetXaxis()->SetBinLabel(5, "# q = -");
3082	for(size_t loc_i = 0; loc_i < dFinalStatePIDs.size(); ++loc_i)
3083	{
3084	string locLabelName = string("# ") + string(ParticleName_ROOT(dFinalStatePIDs[loc_i]));
3085	dHist_NumReconstructedParticles->GetXaxis()->SetBinLabel(6 + loc_i, locLabelName.c_str());
3086	}
3087	}
3088
3089	locHistName = "NumGoodReconstructedParticles";
3090	if(gDirectory(TDirectory::CurrentDirectory())->Get(locHistName.c_str()) != NULL__null) //already created by another thread, or directory name is duplicate (e.g. two identical steps)
3091	dHist_NumGoodReconstructedParticles = static_cast<TH2D*>(gDirectory(TDirectory::CurrentDirectory())->Get(locHistName.c_str()));
3092	else
3093	{
3094	dHist_NumGoodReconstructedParticles = new TH2D("NumGoodReconstructedParticles", ";Particle Type;Num Particles / Event", 5 + dFinalStatePIDs.size(), -0.5, 4.5 + dFinalStatePIDs.size(), dMaxNumTracks + 1, -0.5, (float)dMaxNumTracks + 0.5);
3095	dHist_NumGoodReconstructedParticles->GetXaxis()->SetBinLabel(1, "# Total");
3096	dHist_NumGoodReconstructedParticles->GetXaxis()->SetBinLabel(2, "# q != 0");
3097	dHist_NumGoodReconstructedParticles->GetXaxis()->SetBinLabel(3, "# q = 0");
3098	dHist_NumGoodReconstructedParticles->GetXaxis()->SetBinLabel(4, "# q = +");
3099	dHist_NumGoodReconstructedParticles->GetXaxis()->SetBinLabel(5, "# q = -");
3100	for(size_t loc_i = 0; loc_i < dFinalStatePIDs.size(); ++loc_i)
3101	{
3102	string locLabelName = string("# ") + string(ParticleName_ROOT(dFinalStatePIDs[loc_i]));
3103	dHist_NumGoodReconstructedParticles->GetXaxis()->SetBinLabel(6 + loc_i, locLabelName.c_str());
3104	}
3105	}
3106
3107	//Return to the base directory
3108	ChangeTo_BaseDirectory();
3109	}
3110	japp->RootUnLock(); //RELEASE ROOT LOCK!!
3111	}
3112
3113	bool DHistogramAction_TrackMultiplicity::Perform_Action(JEventLoop* locEventLoop, const DParticleCombo* locParticleCombo)
3114	{
3115	if(Get_NumPreviousParticleCombos() != 0)
3116	return true; //else double-counting!
3117
3118	vector<const DChargedTrack*> locChargedTracks;
3119	locEventLoop->Get(locChargedTracks);
3120
3121	vector<const DChargedTrack*> locGoodChargedTracks;
3122	locEventLoop->Get(locGoodChargedTracks, dTrackSelectionTag.c_str());
3123
3124	// get #tracks by PID/q type
3125	size_t locNumPositiveTracks = 0, locNumNegativeTracks = 0;
3126	map<Particle_t, size_t> locNumTracksByPID;
3127	for(size_t loc_i = 0; loc_i < locChargedTracks.size(); ++loc_i)
3128	{
3129	const DChargedTrackHypothesis* locChargedTrackHypothesis = locChargedTracks[loc_i]->Get_BestFOM();
3130	Particle_t locPID = locChargedTrackHypothesis->PID();
3131
3132	if(locChargedTrackHypothesis->charge() > 0.0)
3133	++locNumPositiveTracks;
3134	else
3135	++locNumNegativeTracks;
3136
3137	if(locNumTracksByPID.find(locPID) != locNumTracksByPID.end())
3138	++locNumTracksByPID[locPID];
3139	else
3140	locNumTracksByPID[locPID] = 1;
3141	}
3142
3143	// get # good tracks by PID/q type
3144	size_t locNumGoodPositiveTracks = 0, locNumGoodNegativeTracks = 0;
3145	map<Particle_t, size_t> locNumGoodTracksByPID;
3146	for(size_t loc_i = 0; loc_i < locGoodChargedTracks.size(); ++loc_i)
3147	{
3148	const DChargedTrackHypothesis* locChargedTrackHypothesis = locGoodChargedTracks[loc_i]->Get_BestFOM();
3149	Particle_t locPID = locChargedTrackHypothesis->PID();
3150
3151	double locPIDFOM = locChargedTrackHypothesis->dFOM;
3152
3153	if(locChargedTrackHypothesis->charge() > 0.0)
3154	++locNumGoodPositiveTracks;
3155	else
3156	++locNumGoodNegativeTracks;
3157
3158	if(locPIDFOM < dMinPIDFOM)
3159	continue;
3160
3161	if(locNumGoodTracksByPID.find(locPID) != locNumGoodTracksByPID.end())
3162	++locNumGoodTracksByPID[locPID];
3163	else
3164	locNumGoodTracksByPID[locPID] = 1;
3165	}
3166
3167	vector<const DNeutralParticle*> locNeutralParticles;
3168	locEventLoop->Get(locNeutralParticles);
3169
3170	vector<const DNeutralShower*> locGoodNeutralShowers;
3171	locEventLoop->Get(locGoodNeutralShowers, dShowerSelectionTag.c_str());
3172
3173	// neutrals by pid
3174	for(size_t loc_i = 0; loc_i < locNeutralParticles.size(); ++loc_i)
3175	{
3176	const DNeutralParticleHypothesis* locNeutralParticleHypothesis = locNeutralParticles[loc_i]->Get_BestFOM();
3177	if(locNeutralParticleHypothesis->dFOM < dMinPIDFOM)
3178	continue;
3179
3180	Particle_t locPID = locNeutralParticleHypothesis->PID();
3181	if(locNumTracksByPID.find(locPID) != locNumTracksByPID.end())
3182	++locNumTracksByPID[locPID];
3183	else
3184	locNumTracksByPID[locPID] = 1;
3185	}
3186
3187	// good neutrals
3188	size_t locNumGoodNeutrals = 0;
3189	for(size_t loc_i = 0; loc_i < locGoodNeutralShowers.size(); ++loc_i)
3190	{
3191	const DNeutralParticle* locNeutralParticle = NULL__null;
3192	for(size_t loc_j = 0; loc_j < locNeutralParticles.size(); ++loc_j)
3193	{
3194	if(locNeutralParticles[loc_j]->dNeutralShower != locGoodNeutralShowers[loc_i])
3195	continue;
3196	locNeutralParticle = locNeutralParticles[loc_j];
3197	break;
3198	}
3199	if(locNeutralParticle == NULL__null)
3200	continue;
3201	++locNumGoodNeutrals;
3202
3203	const DNeutralParticleHypothesis* locNeutralParticleHypothesis = locNeutralParticles[loc_i]->Get_BestFOM();
3204	if(locNeutralParticleHypothesis->dFOM < dMinPIDFOM)
3205	continue;
3206
3207	Particle_t locPID = locNeutralParticleHypothesis->PID();
3208	if(locNumGoodTracksByPID.find(locPID) != locNumGoodTracksByPID.end())
3209	++locNumGoodTracksByPID[locPID];
3210	else
3211	locNumGoodTracksByPID[locPID] = 1;
3212	}
3213
3214	size_t locNumGoodTracks = locNumGoodPositiveTracks + locNumGoodNegativeTracks;
3215	japp->RootWriteLock();
3216	{
3217	dHist_NumReconstructedParticles->Fill(0.0, (Double_t)(locChargedTracks.size() + locNeutralParticles.size()));
3218	dHist_NumReconstructedParticles->Fill(1.0, (Double_t)locChargedTracks.size());
3219	dHist_NumReconstructedParticles->Fill(2.0, (Double_t)locNeutralParticles.size());
3220	dHist_NumReconstructedParticles->Fill(3.0, (Double_t)locNumPositiveTracks);
3221	dHist_NumReconstructedParticles->Fill(4.0, (Double_t)locNumNegativeTracks);
3222	for(size_t loc_i = 0; loc_i < dFinalStatePIDs.size(); ++loc_i)
3223	dHist_NumReconstructedParticles->Fill(5.0 + (Double_t)loc_i, (Double_t)locNumTracksByPID[dFinalStatePIDs[loc_i]]);
3224
3225	dHist_NumGoodReconstructedParticles->Fill(0.0, (Double_t)(locNumGoodTracks + locNumGoodNeutrals));
3226	dHist_NumGoodReconstructedParticles->Fill(1.0, (Double_t)locNumGoodTracks);
3227	dHist_NumGoodReconstructedParticles->Fill(2.0, (Double_t)locNumGoodNeutrals);
3228	dHist_NumGoodReconstructedParticles->Fill(3.0, (Double_t)locNumGoodPositiveTracks);
3229	dHist_NumGoodReconstructedParticles->Fill(4.0, (Double_t)locNumGoodNegativeTracks);
3230	for(size_t loc_i = 0; loc_i < dFinalStatePIDs.size(); ++loc_i)
3231	dHist_NumGoodReconstructedParticles->Fill(5.0 + (Double_t)loc_i, (Double_t)locNumGoodTracksByPID[dFinalStatePIDs[loc_i]]);
3232	}
3233	japp->RootUnLock();
3234
3235	return true;
3236	}
3237