33 #include "Utilities/AssociationUtil.h" 100 produces<std::vector<brb::WCTrack> >();
101 produces<art::Assns<brb::WCTrack, rawdata::RawBeamlineWC> >();
102 produces<std::vector<rb::Vertex> > ();
103 produces<art::Assns<rb::Vertex, brb::WCTrack> > ();
120 fWCDist = tfs->
make<TH1F>(
"WCCond",
"WC Conditions",7,0,7);
121 fReco_P = tfs->
make<TH1F>(
"Reco_P",
"Reconstructed momentum", 180, 0, 1800);
122 fY_Kink = tfs->
make<TH1F>(
"Y_Kink",
"Angle between US/DS tracks in Y direction (degrees)",200,-10*3.1415926/180,10*3.141592654/180);
123 fX_Dist = tfs->
make<TH1F>(
"X_Dist",
"X distance between US/DS tracks at midplane (mm)",1200,-60,1260);
124 fY_Dist = tfs->
make<TH1F>(
"Y_Dist",
"Y distance between US/DS tracks at midplane (mm)",1200,-600,600);
125 fZ_Dist = tfs->
make<TH1F>(
"Z_Dist",
"Z distance between US/DS tracks at midplane (mm)",1200,-60,1260);
126 fX_Mag_Int = tfs->
make<TH1F>(
"X_Mag_Int",
"X intercept of WC track with magnet front face (cm)",100,-150,-100);
127 fX_Mag_Int->GetXaxis()->SetTitle(
"X intercept [cm]");
128 fX_Mag_Int->GetYaxis()->SetTitle(
"Tracks per 2 cm");
129 fY_Mag_Int = tfs->
make<TH1F>(
"Y_Mag_Int",
"Y intercept of WC track with magnet front face (cm)",100,-10,10);
130 fY_Mag_Int->GetXaxis()->SetTitle(
"Y intercept [cm]");
131 fY_Mag_Int->GetYaxis()->SetTitle(
"Tracks per 0.22 cm");
132 fZ_Mag_Int = tfs->
make<TH1F>(
"Z_Mag_Int",
"Z intercept of WC track with magnet front face (cm)",100,415,425);
133 fZ_Mag_Int->GetXaxis()->SetTitle(
"Z intercept [cm]");
134 fZ_Mag_Int->GetYaxis()->SetTitle(
"Tracks per 0.1 cm");
135 fDist_to_Mag_Axis = tfs->
make<TH1F>(
"Dist_to_Mag_Axis",
"Transverse distance to central axis of magnet [cm]",80,-20,20);
136 fDist_to_Mag_Axis->GetXaxis()->SetTitle(
"Transverse distance from WC track to central axis of magnet [cm]");
138 fX_Face_Dist = tfs->
make<TH1F>(
"X_Face",
"X Location of Track's Detector Entry (mm)",1600,-200,1400);
139 fY_Face_Dist = tfs->
make<TH1F>(
"Y_Face",
"Y Location of Track's Detector Entry (mm)",800,-400,400);
140 fTheta_Dist = tfs->
make<TH1F>(
"Theta",
"Track Theta (w.r.t. Detector Z axis), (radians),",400,-.4,0.4);
141 fPhi_Dist = tfs->
make<TH1F>(
"Phi",
"Track Phi (w.r.t. Detector X axis), (radians)",2000,-6.28318,6.28318);
142 fReco_P->GetXaxis()->SetTitle(
"Reconstructed momentum (MeV/c)");
143 fReco_P->GetYaxis()->SetTitle(
"Tracks per 10 MeV/c");
144 fY_Kink->GetXaxis()->SetTitle(
"Reconstructed y_kink (radians)");
145 fY_Kink->GetYaxis()->SetTitle(
"Tracks per 0.000872 radians");
146 fX_Dist->GetXaxis()->SetTitle(
"X distance between US and DS track ends");
147 fX_Dist->GetYaxis()->SetTitle(
"Tracks per 1 mm - DO NOT TRUST THESE VALUES");
148 fY_Dist->GetXaxis()->SetTitle(
"Y distance between US and DS track ends");
149 fY_Dist->GetYaxis()->SetTitle(
"Tracks per 1 mm - DO NOT TRUST THESE VALUES");
150 fZ_Dist->GetXaxis()->SetTitle(
"Z distance between US and DS track ends");
151 fZ_Dist->GetYaxis()->SetTitle(
"Tracks per 1 mm - DO NOT TRUST THESE VALUES");
154 fY_Face_Dist->GetXaxis()->SetTitle(
"Y (mm)");
155 fY_Face_Dist->GetYaxis()->SetTitle(
"Tracks per 1 mm");
156 fTheta_Dist->GetXaxis()->SetTitle(
"Theta (radians)");
157 fTheta_Dist->GetYaxis()->SetTitle(
"Tracks per .002 radians");
158 fPhi_Dist->GetXaxis()->SetTitle(
"Phi (radians)");
159 fPhi_Dist->GetYaxis()->SetTitle(
"Tracks per 0.00628 radians");
160 fTrack_Type = tfs->
make<TH1F>(
"TrackType",
"WCTrack conditions: 1=missHit,2=uniqueHits,3=lonelyHit,4=socialHits",4,0,4);
162 fTrack_Type->GetXaxis()->SetTitle(
"Track Conditions");
170 std::unique_ptr<std::vector<brb::WCTrack> > wc_reco(
new std::vector<brb::WCTrack>);
171 std::unique_ptr<art::Assns<brb::WCTrack, rawdata::RawBeamlineWC> >
173 std::unique_ptr<std::vector<rb::Vertex> > wc_vtx(
new std::vector<rb::Vertex>);
174 std::unique_ptr<art::Assns<rb::Vertex, brb::WCTrack> >
179 std::vector<art::Ptr<rawdata::RawBeamlineWC> > wcs;
192 mf::LogInfo(
"TBSpillInfo") <<
"MW: Change in reconstruction, March 27.\n" 193 <<
"The magnet information is now obtained from IFDB via the TBSpillInfo module, which need to be run first.\n" 194 <<
"Be sure to have run this module before WCTrackReco. (Alternatively, turn off the fUseIFDB flag (temporary))." 200 unsigned int beamline_run = 0;
216 std::vector<std::vector<WCHitList> > good_hits;
222 std::vector<double> reco_p_list;
223 std::vector<TVector3> mag_int_list;
224 std::vector<double> dist_to_mag_axis_list;
225 std::vector<TVector2> face_list;
226 std::vector<std::vector<TVector3> > wc_hit_pos_list;
227 std::vector<TVector3> dir_list;
228 std::vector<double> y_kink_list;
229 std::vector<TVector3> dist_list;
230 std::vector<WCHitList> final_tracks;
235 dist_to_mag_axis_list,
247 for (
unsigned int iNewTrack = 0; iNewTrack < final_tracks.size(); ++iNewTrack) {
249 WCHitList final_track = final_tracks[iNewTrack];
265 y_kink_list[iNewTrack],
266 dist_list[iNewTrack],
267 mag_int_list[iNewTrack],
268 dist_to_mag_axis_list[iNewTrack],
269 face_list[iNewTrack],
270 wc_hit_pos_list[iNewTrack],
272 dir_list[iNewTrack].
Theta(),
273 dir_list[iNewTrack].
Phi(),
276 fReco_P->Fill(reco_p_list[iNewTrack]);
277 fY_Kink->Fill(y_kink_list[iNewTrack]);
278 fX_Dist->Fill(dist_list[iNewTrack].
X());
279 fY_Dist->Fill(dist_list[iNewTrack].
Y());
280 fZ_Dist->Fill(dist_list[iNewTrack].
Z());
294 (*wc_reco).push_back(the_track);
296 (*wc_vtx).push_back(vtx);
297 util::CreateAssn(*
this, evt, *(vtx_to_wc_trk.get()), *(wc_vtx.get()), *(wc_reco.get()));
301 evt.
put(std::move(wc_reco));
302 evt.
put(std::move(wc_assn));
303 evt.
put(std::move(wc_vtx));
304 evt.
put(std::move(vtx_to_wc_trk));
art::InputTag fRawWCDataLabel
A 3D position and time representing an interaction vertex.
static bool CreateAssn(art::EDProducer const &prod, art::Event &evt, std::vector< T > &a, art::Ptr< U > b, art::Assns< T, U > &assn, size_t indx=UINT_MAX, std::string const &instance=std::string())
Create a 1 to 1 association between a new product and one already in the event.
WCTrackReco(const fhicl::ParameterSet &pset)
void createHits(const std::vector< art::Ptr< rawdata::RawBeamlineWC > > &wcs, std::vector< std::vector< WCHitList > > &good_hits, bool verbose)
MaybeLogger_< ELseverityLevel::ELsev_info, false > LogInfo
void setMagneticField(float field)
unsigned int BeamlineRun() const
WCHitFinderAlg fWCHitFinderAlg
::xsd::cxx::tree::exception< char > exception
void produce(art::Event &evt)
This module reads IFDB DB and then stores BNB spill info.
Definition of the beamline offline geometry. Also implementation of a service to obtain this informat...
DEFINE_ART_MODULE(TestTMapFile)
art::InputTag fRawBeamlineConfigLabel
art::ServiceHandle< beamlineutil::BeamlineMagneticField > fMagneticField
ProductID put(std::unique_ptr< PROD > &&product)
Utility to provide the magnetic field in the beamline for a given run number/run conditions. NB/ Currently this will use a lookup defined in a fhicl file, it should and will be replaced by an interface to a database when the field is stored in such a way.
void reconfigure(const fhicl::ParameterSet &pset)
T get(std::string const &key) const
Encapsulation of reconstructed Wire Chamber track. Part of beamline reconstruction for NOvA test beam...
fvar< T > Phi(const fvar< T > &x)
Vertex location in position and time.
T * make(ARGS...args) const
bool getByLabel(std::string const &label, std::string const &productInstanceName, Handle< PROD > &result) const
Algorithm implementation for wire chamber reconstruction. Part of the beamline reconstruction for NOv...
float FieldFromCurrent(float magnet_current)
Convert magnet current to magnetic field using field mapping information.
double mc7magnet
Tertiary beamline magnet current [Amps].
Raw data definitions for beamline data used in NOvA test beam experiment.
float MagneticField(unsigned int beamline_run)
Return the magnetic field for a given beamline data run.
void fill_ptr_vector(std::vector< Ptr< T >> &ptrs, H const &h)
art::InputTag fMCenterDataLabel
void reconstructTracks(std::vector< double > &reco_p_list, std::vector< TVector3 > &mag_int_list, std::vector< double > &dist_to_mag_axis_list, std::vector< TVector2 > &face_list, std::vector< std::vector< TVector3 > > &wc_hit_pos_list, std::vector< TVector3 > &dir_list, std::vector< WCHitList > &event_final_tracks, std::vector< std::vector< WCHitList > > &good_hits, std::vector< double > &y_kink_list, std::vector< TVector3 > &dist_list, int &WCMissed, TH1F *&WCdistribution, float &residual)
Main function called for each event.
const float Theta() const