Public Member Functions | Private Attributes | List of all members
dune::apa::APAGeometryAlg Class Reference

#include <APAGeometryAlg.h>

Public Member Functions

 APAGeometryAlg (fhicl::ParameterSet const &pset)
 
 APAGeometryAlg ()
 
virtual ~APAGeometryAlg ()
 
void reconfigure (fhicl::ParameterSet const &p)
 
void Init ()
 Initialize some chanel numbers to speed up other methods. More...
 
bool APAChannelsIntersect (uint32_t chan1, uint32_t chan2, std::vector< geo::WireIDIntersection > &IntersectVector)
 If the channels intersect, get all intersections. More...
 
bool LineSegChanIntersect (TVector3 xyzStart, TVector3 xyzEnd, uint32_t chan, std::vector< geo::WireID > &widsCrossed, bool ExtendLine)
 If a line given by start/end points intersects a channel. More...
 
std::vector< geo::WireIDChanSegsPerSide (uint32_t chan, unsigned int side)
 
std::vector< geo::WireIDChanSegsPerSide (std::vector< geo::WireID > wids, unsigned int side)
 
std::vector< double > ThreeChanPos (uint32_t u, uint32_t v, uint32_t z)
 Find the center of the 3 intersections, choose best if multiple. More...
 
geo::WireID NearestWireIDOnChan (const double WorldLoc[3], uint32_t chan, unsigned int const plane, unsigned int const tpc=0, unsigned int const cstat=0)
 
unsigned int ChannelToAPA (uint32_t chan)
 Get number of the APA containing the given channel. More...
 
void ChannelToAPA (uint32_t chan, unsigned int &apa, unsigned int &cryo)
 
APAView_t APAView (uint32_t chan)
 Get which of the 4 APA views the channel is in. More...
 
unsigned int ChannelsInView (geo::View_t geoview)
 
uint32_t FirstChannelInView (geo::View_t geoview, unsigned int apa, unsigned int cryo)
 
uint32_t FirstChannelInView (geo::View_t geoview, uint32_t chan)
 
uint32_t FirstChannelInView (uint32_t chan)
 
unsigned int ChannelsInAPAView (APAView_t apaview)
 
unsigned int ChannelsPerAPA ()
 

Private Attributes

art::ServiceHandle< geo::GeometryfGeom
 
unsigned int fChannelsPerAPA
 All APAs have this same number of channels. More...
 
unsigned int fAPAsPerCryo
 
uint32_t fFirstU
 
uint32_t fLastU
 
uint32_t fFirstV
 
uint32_t fLastV
 
uint32_t fFirstZ0
 
uint32_t fLastZ0
 
uint32_t fFirstZ1
 
uint32_t fLastZ1
 
double fChannelRange [2]
 

Detailed Description

Definition at line 39 of file APAGeometryAlg.h.

Constructor & Destructor Documentation

dune::apa::APAGeometryAlg::APAGeometryAlg ( fhicl::ParameterSet const &  pset)

Definition at line 39 of file APAGeometryAlg.cxx.

40  {
41  this->reconfigure(pset);
42  this->Init();
43  }
void Init()
Initialize some chanel numbers to speed up other methods.
void reconfigure(fhicl::ParameterSet const &p)
dune::apa::APAGeometryAlg::APAGeometryAlg ( )

Definition at line 47 of file APAGeometryAlg.cxx.

48  {
49  this->Init();
50  }
void Init()
Initialize some chanel numbers to speed up other methods.
dune::apa::APAGeometryAlg::~APAGeometryAlg ( )
virtual

Definition at line 53 of file APAGeometryAlg.cxx.

54  {
55  }

Member Function Documentation

bool dune::apa::APAGeometryAlg::APAChannelsIntersect ( uint32_t  chan1,
uint32_t  chan2,
std::vector< geo::WireIDIntersection > &  IntersectVector 
)

If the channels intersect, get all intersections.

Definition at line 490 of file APAGeometryAlg.cxx.

492  {
493 
494 
495  // Get the WireIDs and view for each channel, make sure views are different
496  geo::WireIDIntersection widIntersect;
497  std::vector< geo::WireID > wids1 = fGeom->ChannelToWire( chan1 );
498  std::vector< geo::WireID > wids2 = fGeom->ChannelToWire( chan2 );
499  geo::View_t view1 = fGeom->View( chan1 );
500  geo::View_t view2 = fGeom->View( chan2 );
501  if( view1 == view2 ){
502  mf::LogWarning("APAChannelsIntersect") << "Comparing two channels in the same view, return false";
503  return false; }
504  if( wids1[0].Cryostat != wids2[0].Cryostat ||
505  this->ChannelToAPA(chan1) != this->ChannelToAPA(chan2) ){
506  throw cet::exception("APAChannelsIntersect") << "Comparing two channels in in different APAs: "
507  << "channel " << chan1 << " in Cryo "
508  << wids1[0].Cryostat << ", APA " << this->ChannelToAPA(chan1)
509  << ", and channel " << chan2 << " in Cryo "
510  << wids2[0].Cryostat << ", APA " << this->ChannelToAPA(chan2)
511  << "\n";
512  return false; }
513 
514 
515  // Loop through wids1 and see if wids2 has any intersecting wires,
516  // given that the WireIDs are in the same TPC
517  for( unsigned int i1 = 0; i1 < wids1.size() ; i1++){
518  for( unsigned int i2 = 0; i2 < wids2.size() ; i2++){
519 
520  // make sure it is reasonable to intersect
521  if( wids1[i1].Plane == wids2[i2].Plane ||
522  wids1[i1].TPC != wids2[i2].TPC || // not reasonable for a *WireID*
523  wids1[i1].Cryostat != wids2[i2].Cryostat ) continue;
524 
525 // std::cout << "Checking: \n WireID 1 = ("
526 // << wids1[i1].Cryostat << "," << wids1[i1].TPC << ","
527 // << wids1[i1].Plane << "," << wids1[i1].Wire
528 // << ") \n WireID 2 = ("
529 // << wids2[i2].Cryostat << "," << wids2[i2].TPC << ","
530 // << wids2[i2].Plane << "," << wids2[i2].Wire << ")" << std::endl;
531 
532  // Check if they even intersect; if they do, push back
533  if( fGeom->WireIDsIntersect( wids1[i1], wids2[i2], widIntersect ) ){
534 
535 // std::cout << "we have an intersect" << std::endl;
536 
537  IntersectVector.push_back( widIntersect );
538  }
539  }
540  }
541 
542  // Of all considered configurations, there are never more than
543  // 4 intersections per channel pair
544  if( IntersectVector.size() > 4 ){
545  mf::LogWarning("APAChannelsIntersect") << "Got " << IntersectVector.size()
546  << " intersections for channels "
547  << chan1 << " and " << chan2
548  << " - never expect more than 4, so far"; }
549 
550 
551  // With increasing IntersectVector index, the WireID
552  // vector indices of the intersecting wireIDs increase.
553  // This matches the direction in which ChannelToWire
554  // builds its output WireID vector in the APA/35t Alg
555  std::sort( IntersectVector.begin(), IntersectVector.end() );
556 
557  // return true if any intersection points were found
558  if( IntersectVector.size() == 0 ) return false;
559  else return true;
560 
561  }
enum geo::_plane_proj View_t
Enumerate the possible plane projections.
std::vector< geo::WireID > ChannelToWire(raw::ChannelID_t const channel) const
Returns a list of wires connected to the specified TPC channel.
unsigned int ChannelToAPA(uint32_t chan)
Get number of the APA containing the given channel.
art::ServiceHandle< geo::Geometry > fGeom
View_t View(geo::PlaneID const &pid) const
Returns the view (wire orientation) on the channels of specified TPC plane.
bool WireIDsIntersect(WireID const &wid1, WireID const &wid2, geo::Point_t &intersection) const
Computes the intersection between two wires.
MaybeLogger_< ELseverityLevel::ELsev_warning, false > LogWarning
recob::tracking::Plane Plane
Definition: TrackState.h:17
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33
APAView_t dune::apa::APAGeometryAlg::APAView ( uint32_t  chan)

Get which of the 4 APA views the channel is in.

Definition at line 229 of file APAGeometryAlg.cxx.

229  {
230 
231  // it seems trivial to do this for U and V, but this gives a side to
232  // geo::kZ, unlike Geometry::View(c), as is often needed in disambiguation
233 
234  geo::View_t view = fGeom->View( chan );
235  switch(view){
236  default :
237  break;
238  case geo::kU :
239  return kU;
240  case geo::kV :
241  return kV;
242  case geo::kZ :
243  unsigned int modchan = chan % fChannelsPerAPA;
244  // Channel mapping number in the order of U, V, Z0, then Z1
245  if( modchan > fLastV && modchan < fFirstZ1 ) return kZ0;
246  if( modchan > fLastZ0 ) return kZ1;
247  }
248 
249  return kUnknown;
250 
251  }
enum geo::_plane_proj View_t
Enumerate the possible plane projections.
Planes which measure V.
Definition: geo_types.h:130
Z view on the smaller-x side of the APA.
Planes which measure Z direction.
Definition: geo_types.h:132
U view on both sides of the APA.
Planes which measure U.
Definition: geo_types.h:129
Z view on the larger-x side of the APA.
art::ServiceHandle< geo::Geometry > fGeom
View_t View(geo::PlaneID const &pid) const
Returns the view (wire orientation) on the channels of specified TPC plane.
unsigned int fChannelsPerAPA
All APAs have this same number of channels.
V view on both sides of the APA.
unsigned int dune::apa::APAGeometryAlg::ChannelsInAPAView ( APAView_t  apaview)

Definition at line 168 of file APAGeometryAlg.cxx.

168  {
169 
170  switch(apaview){
171  default :
172  return 0;
173  case kU :
174  return fFirstV-fFirstU;
175  case kV :
176  return fFirstZ0-fFirstV;
177  case kZ0 :
178  return fFirstZ1-fFirstZ0;
179  case kZ1 :
180  return fLastZ1-fFirstZ1 + 1;
181  }
182 
183  }
Z view on the smaller-x side of the APA.
U view on both sides of the APA.
Z view on the larger-x side of the APA.
V view on both sides of the APA.
unsigned int dune::apa::APAGeometryAlg::ChannelsInView ( geo::View_t  geoview)

Definition at line 148 of file APAGeometryAlg.cxx.

148  {
149 
150  switch(geoview){
151  default :
152  return 0;
153  case geo::kU :
154  return ChannelsInAPAView( kU );
155  case geo::kV :
156  return ChannelsInAPAView( kV );
157  case geo::kZ :
159  return ChannelsInAPAView( kZ0 );
160  else throw cet::exception("ChannelsInView")
161  << "Both Z sides should have the same amount of channels\n";
162  }
163 
164  }
Planes which measure V.
Definition: geo_types.h:130
Z view on the smaller-x side of the APA.
Planes which measure Z direction.
Definition: geo_types.h:132
U view on both sides of the APA.
Planes which measure U.
Definition: geo_types.h:129
Z view on the larger-x side of the APA.
unsigned int ChannelsInAPAView(APAView_t apaview)
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33
V view on both sides of the APA.
unsigned int dune::apa::APAGeometryAlg::ChannelsPerAPA ( )
inline

Definition at line 88 of file APAGeometryAlg.h.

88 { return fChannelsPerAPA; };
unsigned int fChannelsPerAPA
All APAs have this same number of channels.
unsigned int dune::apa::APAGeometryAlg::ChannelToAPA ( uint32_t  chan)

Get number of the APA containing the given channel.

Definition at line 141 of file APAGeometryAlg.cxx.

141  {
142 
143  return chan / fChannelsPerAPA;
144  }
unsigned int fChannelsPerAPA
All APAs have this same number of channels.
void dune::apa::APAGeometryAlg::ChannelToAPA ( uint32_t  chan,
unsigned int &  apa,
unsigned int &  cryo 
)

Definition at line 125 of file APAGeometryAlg.cxx.

127  {
128 
129  cryo = chan / (fAPAsPerCryo*fChannelsPerAPA);
130 
131  // Number apa uniquely across cryostats so that
132  // apa to recob::Object maps are easy to work with.
133  // If we decide to reset APA number per cryo, uncomment:
134  //chan -= cryo*fAPAsPerCryo*fChannelsPerAPA;
135  apa = chan / fChannelsPerAPA;
136 
137  return;
138  }
unsigned int fChannelsPerAPA
All APAs have this same number of channels.
std::vector< geo::WireID > dune::apa::APAGeometryAlg::ChanSegsPerSide ( uint32_t  chan,
unsigned int  side 
)

Definition at line 255 of file APAGeometryAlg.cxx.

255  {
256 
257  std::vector<geo::WireID> wids = fGeom->ChannelToWire(chan);
258  return this->ChanSegsPerSide(wids, side);
259 
260  }
std::vector< geo::WireID > ChanSegsPerSide(uint32_t chan, unsigned int side)
std::vector< geo::WireID > ChannelToWire(raw::ChannelID_t const channel) const
Returns a list of wires connected to the specified TPC channel.
art::ServiceHandle< geo::Geometry > fGeom
std::vector< geo::WireID > dune::apa::APAGeometryAlg::ChanSegsPerSide ( std::vector< geo::WireID wids,
unsigned int  side 
)

Definition at line 265 of file APAGeometryAlg.cxx.

266  {
267  // Given a vector of wireIDs and an APA side, return
268  // the wireIDs the the tpc side where tpc%2 = side
269 
270  std::vector<geo::WireID> thisSide;
271 
272  for(size_t i = 0; i < wids.size(); i++)
273  if( wids[i].TPC % 2 == side ) thisSide.push_back(wids[i]);
274 
275  return thisSide;
276  }
uint32_t dune::apa::APAGeometryAlg::FirstChannelInView ( geo::View_t  geoview,
unsigned int  apa,
unsigned int  cryo 
)

Definition at line 187 of file APAGeometryAlg.cxx.

189  {
190 
191 
192  switch(geoview){
193  default :
194  return 0 + (uint32_t)(apa + cryo*fAPAsPerCryo)*fChannelsPerAPA;
195  case geo::kU :
196  return fFirstU + (uint32_t)(apa + cryo*fAPAsPerCryo)*fChannelsPerAPA;
197  case geo::kV :
198  return fFirstV + (uint32_t)(apa + cryo*fAPAsPerCryo)*fChannelsPerAPA;
199  case geo::kZ :
200  //TODO: would need tpc number for the rest of this
201  return fFirstZ0 + (uint32_t)(apa + cryo*fAPAsPerCryo)*fChannelsPerAPA;
202  }
203 
204  }
Planes which measure V.
Definition: geo_types.h:130
Planes which measure Z direction.
Definition: geo_types.h:132
Planes which measure U.
Definition: geo_types.h:129
unsigned int fChannelsPerAPA
All APAs have this same number of channels.
uint32_t dune::apa::APAGeometryAlg::FirstChannelInView ( geo::View_t  geoview,
uint32_t  chan 
)

Definition at line 217 of file APAGeometryAlg.cxx.

218  {
219 
220  unsigned int apa, cryo;
221  this->ChannelToAPA( chan, apa, cryo );
222  return this->FirstChannelInView( geoview, apa, cryo );
223 
224  }
unsigned int ChannelToAPA(uint32_t chan)
Get number of the APA containing the given channel.
uint32_t FirstChannelInView(geo::View_t geoview, unsigned int apa, unsigned int cryo)
uint32_t dune::apa::APAGeometryAlg::FirstChannelInView ( uint32_t  chan)

Definition at line 207 of file APAGeometryAlg.cxx.

207  {
208 
209  geo::View_t geoview = fGeom->View(chan);
210  unsigned int apa, cryo;
211  this->ChannelToAPA( chan, apa, cryo );
212  return this->FirstChannelInView( geoview, chan );
213 
214  }
enum geo::_plane_proj View_t
Enumerate the possible plane projections.
unsigned int ChannelToAPA(uint32_t chan)
Get number of the APA containing the given channel.
art::ServiceHandle< geo::Geometry > fGeom
View_t View(geo::PlaneID const &pid) const
Returns the view (wire orientation) on the channels of specified TPC plane.
uint32_t FirstChannelInView(geo::View_t geoview, unsigned int apa, unsigned int cryo)
void dune::apa::APAGeometryAlg::Init ( )

Initialize some chanel numbers to speed up other methods.

Definition at line 63 of file APAGeometryAlg.cxx.

64  {
65 
66  // find the number of channels per APA
67  uint32_t channel = 0;
68  // old logic -- segfaults if there is just one APA
69  //while( fGeom->ChannelToWire(channel+1)[0].TPC < 2 ) channel++;
70  //fChannelsPerAPA = channel + 1;
71 
73  for (channel=0; channel < fGeom->Nchannels(); ++channel)
74  {
75  if ( fGeom->ChannelToWire(channel)[0].TPC > 1 )
76  {
78  break;
79  }
80  }
81 
82  // Step through channel c and find the view boundaries, until
83  // outside of first APA - these help optimize ChannelToAPAView
84  // (very dependent on the conventions implimented in the channel map)
85  fFirstU = 0;
86  uint32_t c = 1;
87  geo::WireID wid = (fGeom->ChannelToWire(c))[0];
88  geo::WireID lastwid;
89  while( wid.TPC < 2 ){
90 
91  if( fGeom->View(c) == geo::kV && fGeom->View(c-1) == geo::kU ){
92  fLastU = c-1;
93  fFirstV = c; }
94 
95  if( fGeom->View(c) == geo::kZ && fGeom->View(c-1) == geo::kV ){
96  fLastV = c-1;
97  fFirstZ0 = c; }
98 
99  if( wid.TPC == lastwid.TPC + 1 ){
100  fLastZ0 = c-1;
101  fFirstZ1 = c; }
102 
103  lastwid = wid;
104  c++;
105  if (c >= fGeom->Nchannels()) break;
106  wid = (fGeom->ChannelToWire(c))[0]; // for the while condition
107 
108  }
109 
110  fLastZ1 = c - 1;
111 
112 
113  if( fLastZ1 + 1 != fChannelsPerAPA ) throw cet::exception("APAGeometryAlg")
114  << "Channel boundaries are inconsistent.\n";
115 
116  // some other things that will be needed
117  fAPAsPerCryo = fGeom->NTPC(0)/2;
120 
121  }
Planes which measure V.
Definition: geo_types.h:130
std::vector< geo::WireID > ChannelToWire(raw::ChannelID_t const channel) const
Returns a list of wires connected to the specified TPC channel.
Planes which measure Z direction.
Definition: geo_types.h:132
uint8_t channel
Definition: CRTFragment.hh:201
geo::Length_t WirePitch(geo::PlaneID const &planeid) const
Returns the distance between two consecutive wires.
unsigned int Nchannels() const
Returns the number of TPC readout channels in the detector.
Planes which measure U.
Definition: geo_types.h:129
art::ServiceHandle< geo::Geometry > fGeom
View_t View(geo::PlaneID const &pid) const
Returns the view (wire orientation) on the channels of specified TPC plane.
unsigned int NTPC(unsigned int cstat=0) const
Returns the total number of TPCs in the specified cryostat.
TPCID_t TPC
Index of the TPC within its cryostat.
Definition: geo_types.h:406
unsigned int fChannelsPerAPA
All APAs have this same number of channels.
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33
bool dune::apa::APAGeometryAlg::LineSegChanIntersect ( TVector3  xyzStart,
TVector3  xyzEnd,
uint32_t  chan,
std::vector< geo::WireID > &  widsCrossed,
bool  ExtendLine = true 
)

If a line given by start/end points intersects a channel.

Definition at line 317 of file APAGeometryAlg.cxx.

320  {
321 
322  // This assumes a smooth wire numbering, and that the line seg is contained in a tpc.
323  // Meant for use with the approximate line calculated
324  // by matching cluster endpoints in disambiguation.
325 
326  // Find tpc, use midpoint in case start/end is on a boundary
327  unsigned int tpc, cryo;
328  double xyzMid[3];
329  xyzMid[0] = (xyzStart[0]+xyzEnd[0])/2;
330  xyzMid[1] = (xyzStart[1]+xyzEnd[1])/2;
331  xyzMid[2] = (xyzStart[2]+xyzEnd[2])/2;
332  fGeom->PositionToTPC(xyzMid, tpc, cryo);
333 
334  // Find the nearest wire number to the line segment endpoints
335  std::vector<geo::WireID> wids = fGeom->ChannelToWire(chan);
336  unsigned int startW = fGeom->NearestWire( xyzStart, wids[0].Plane, tpc, cryo );
337  unsigned int endW = fGeom->NearestWire( xyzEnd, wids[0].Plane, tpc, cryo );
338 
339  if( startW > endW ) std::swap(startW, endW);
340 
341 
342  // Loop through wireIDs and check for intersection, if in the right TPC
343  for( size_t w = 0; w < wids.size(); w++ ){
344  if( wids[w].TPC != tpc ) continue;
345  if( wids[w].Cryostat != cryo ) throw cet::exception("LineSegChanIntersect")
346  << "Channel and line not in the same crostat.\n";
347 
348  // If the current wire id wire number is inbetween the start/end
349  // point wires, the line segment intersects the wireID at some point.
350 
351  // TODO: for now, extend range, but that is application specific. fix asap
352  // The longer we make the range, the more conservative it is, so it is safe
353  // to extend the range a bit to get hits at the ends of the line
354  unsigned int ext = 0;
355  if ( ExtendLine) ext = 10;
356 
357  if( fGeom->ValueInRange( wids[w].Wire*1., (startW-ext)*1., (endW+ext)*1. ) ) widsCrossed.push_back(wids[w]);
358 
359  }
360 
361  if( widsCrossed.size() > 0 ) return true;
362  else return false;
363 
364  }
std::vector< geo::WireID > ChannelToWire(raw::ChannelID_t const channel) const
Returns a list of wires connected to the specified TPC channel.
geo::TPCGeo const & PositionToTPC(geo::Point_t const &point) const
Returns the TPC at specified location.
void swap(Handle< T > &a, Handle< T > &b)
geo::WireID::WireID_t NearestWire(geo::Point_t const &point, geo::PlaneID const &planeid) const
Returns the index of wire closest to position in the specified TPC.
art::ServiceHandle< geo::Geometry > fGeom
bool ValueInRange(double value, double min, double max) const
Returns whether a value is within the specified range.
recob::tracking::Plane Plane
Definition: TrackState.h:17
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33
geo::WireID dune::apa::APAGeometryAlg::NearestWireIDOnChan ( const double  WorldLoc[3],
uint32_t  chan,
unsigned int const  plane,
unsigned int const  tpc = 0,
unsigned int const  cstat = 0 
)

Definition at line 282 of file APAGeometryAlg.cxx.

287  {
288 
289  std::vector<geo::WireID> cWids = fGeom->ChannelToWire( chan );
290 
291  if( cWids[0].Cryostat != cstat )
292  throw cet::exception("APAGeometryAlg") << "Channel " << chan
293  << "not in cryostat " << cstat << "\n";
294  if( std::floor( cWids[0].TPC / 2 ) != std::floor( tpc / 2 ) )
295  throw cet::exception("APAGeometryAlg") << "Channel " << chan
296  << "not in APA " << std::floor(tpc/2) << "\n";
297 
298  // special case for vertical wires
299  if(fGeom->View(chan)==geo::kZ) return fGeom->ChannelToWire(chan)[0];
300 
301  unsigned int xyzWire = fGeom->NearestWireID( WorldLoc, plane, tpc, cstat ).Wire;
302 
303  // The desired wire ID will be the only channel
304  // segment within half the channel range.
305  geo::WireID wid;
306  for(size_t i=0; i<cWids.size(); i++){
307  if( cWids[i].TPC != tpc ) continue;
308  if( std::abs((int)cWids[i].Wire - (int)xyzWire) < fChannelRange[plane]/2 ) wid=cWids[i];
309  }
310 
311  return wid;
312 
313  }
std::vector< geo::WireID > ChannelToWire(raw::ChannelID_t const channel) const
Returns a list of wires connected to the specified TPC channel.
Planes which measure Z direction.
Definition: geo_types.h:132
WireID_t Wire
Index of the wire within its plane.
Definition: geo_types.h:580
T abs(T value)
art::ServiceHandle< geo::Geometry > fGeom
View_t View(geo::PlaneID const &pid) const
Returns the view (wire orientation) on the channels of specified TPC plane.
geo::WireID NearestWireID(geo::Point_t const &point, geo::PlaneID const &planeid) const
Returns the ID of wire closest to position in the specified TPC.
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33
void dune::apa::APAGeometryAlg::reconfigure ( fhicl::ParameterSet const &  p)

Definition at line 58 of file APAGeometryAlg.cxx.

59  {
60  }
std::vector< double > dune::apa::APAGeometryAlg::ThreeChanPos ( uint32_t  u,
uint32_t  v,
uint32_t  z 
)

Find the center of the 3 intersections, choose best if multiple.

Definition at line 369 of file APAGeometryAlg.cxx.

370  {
371 
372  // Say we've associated a U, V, and Z channel -- perhaps by associating hits
373  // or cluster endpoints -- these don't necessarily all intersect, but they
374  // are hopefully pretty close. Find the center of the 3 intersections.
375 
376  // get data needed along the way
377  std::vector< geo::WireIDIntersection > UVIntersects;
378  this->APAChannelsIntersect( u, v, UVIntersects );
379  std::vector< double > UVzToZ(UVIntersects.size());
380  geo::WireID Zwid = fGeom->ChannelToWire(z)[0];
381  unsigned int cryo = Zwid.Cryostat;
382  unsigned int tpc = Zwid.TPC;
383  std::vector<geo::WireID> Uwids = fGeom->ChannelToWire(u);
384  std::vector<geo::WireID> Vwids = fGeom->ChannelToWire(v);
385  std::vector<geo::WireID> UwidsInTPC, VwidsInTPC;
386  for(size_t i=0; i<Uwids.size(); i++) if( Uwids[i].TPC==tpc ) UwidsInTPC.push_back(Uwids[i]);
387  for(size_t i=0; i<Vwids.size(); i++) if( Vwids[i].TPC==tpc ) VwidsInTPC.push_back(Vwids[i]);
388  double Zcent[3] = {0.};
389  fGeom->WireIDToWireGeo( Zwid ).GetCenter(Zcent);
390 
391  std::cout << "Zcent = " << Zcent[2] << ", UVintersects zpos = ";
392  for(size_t uv=0; uv<UVIntersects.size(); uv++){
393  std::cout << UVIntersects[uv].z << ", ";
394  }
395  std::cout << "\n";
396 
397  /////////////////////////////////////
398  /////////////////////////////////////
399 
400  //std::cout << "U = " << u << " V = " << v << ", " << UVIntersects.size() << std::endl;
401 
402  if( UVIntersects.size() == 0 ){
403  if( UwidsInTPC.size() > 1 || VwidsInTPC.size() > 1 )
404  throw cet::exception("ThreeChanPos") << "U/V channels don't intersect, bad return.\n";
405 
406  // Now assume there are only one of each u and v wireIDs on in this TPC
407  mf::LogWarning("ThreeChanPos") << "No U/V intersect, exceptional channels. See if U or V intersects Z\n";
408  std::vector<double> yzCenter(2,0.);
409  geo::WireID Uwid = UwidsInTPC[0];
410  geo::WireID Vwid = VwidsInTPC[0];
411  geo::WireIDIntersection UZInt, VZInt;
412  bool checkUZ = fGeom->WireIDsIntersect( Uwid, Zwid, UZInt );
413  bool checkVZ = fGeom->WireIDsIntersect( Vwid, Zwid, VZInt );
414  if( !checkUZ && !checkVZ )
415  throw cet::exception("NoChanIntersect") << "No channels intersect, bad return.\n";
416  if( checkUZ && !checkVZ ){ yzCenter[0] = UZInt.y; yzCenter[1] = UZInt.z; }
417  if( checkVZ && !checkUZ ){ yzCenter[0] = VZInt.y; yzCenter[1] = VZInt.z; }
418  if( checkUZ && checkVZ ){ yzCenter[0] = (VZInt.y+UZInt.y)/2; yzCenter[1] = (VZInt.z+UZInt.z)/2; }
419  return yzCenter;
420  }
421 
422  /////////////////////////////////////
423  /////////////////////////////////////
424 
425 
426  // In case the uv channels intersect twice on the same side, choose the best case.
427  // Note: this will not happen for APAs with UV angle at about 36, but will for 45
428  std::cout << "UVzToZ = ";
429  for( size_t widI = 0; widI < UVIntersects.size(); widI++ ){
430  UVzToZ[widI] = std::abs( UVIntersects[widI].z - Zcent[2] );
431  std::cout << UVzToZ[widI] << ", ";
432  }
433  std::cout<<"\n";
434 
435  unsigned int bestWidI = 0;
436  double minZdiff = fGeom->Cryostat(cryo).TPC(tpc).Length(); // start it out at maximum z
437  for( unsigned int widI = 0; widI < UVIntersects.size(); widI++ ){
438 
439  //std::cout << "widI = " << widI << std::endl;
440 
441  if( UVIntersects[widI].TPC == tpc && UVzToZ[widI] < minZdiff ){
442  minZdiff = UVzToZ[widI];
443  bestWidI = widI;
444  //std::cout << "bestWidI = " << bestWidI << std::endl;
445  }
446  }
447  geo::WireIDIntersection ChosenUVInt = UVIntersects[bestWidI];
448 
449  // Now having the UV intersection, get the UZ and VZ
450  double UVInt[3] = {0.};
451  UVInt[1] = ChosenUVInt.y; UVInt[2] = ChosenUVInt.z;
452  geo::WireID Uwid = this->NearestWireIDOnChan( UVInt, u, 0, tpc, cryo );
453  geo::WireID Vwid = this->NearestWireIDOnChan( UVInt, v, 1, tpc, cryo );
454  geo::WireIDIntersection UZInt, VZInt;
455  bool checkUZ = fGeom->WireIDsIntersect( Uwid, Zwid, UZInt );
456  bool checkVZ = fGeom->WireIDsIntersect( Vwid, Zwid, VZInt );
457 
458  std::cout << "UZint.z = " << UZInt.z << " (" << checkUZ << "), VZint.z = " << VZInt.z << " (" << checkVZ << ")\n";
459 
460  // find the center
461  std::vector<double> yzCenter(2,0.);
462 
463 
464  if( !checkUZ || !checkVZ ){
465  //throw cet::exception("ThreeChanPos") << "WireIDs were expected to intersect.\n";
466 
467  std::cout << "ChosenUVint.y = " << ChosenUVInt.y << "ChosenUVint.z = " << ChosenUVInt.z << std::endl;
468 
469  //temporary case
470  yzCenter[0] = ChosenUVInt.y;
471  yzCenter[1] = ChosenUVInt.z;
472 
473 
474  } else {
475 
476  yzCenter[0] = (ChosenUVInt.y + UZInt.y + VZInt.y)/3;
477  yzCenter[1] = (ChosenUVInt.z + UZInt.z + VZInt.z)/3;
478 
479  }
480 
481 
482  return yzCenter;
483 
484  }
double z
z position of intersection
Definition: geo_types.h:805
geo::WireID NearestWireIDOnChan(const double WorldLoc[3], uint32_t chan, unsigned int const plane, unsigned int const tpc=0, unsigned int const cstat=0)
std::vector< geo::WireID > ChannelToWire(raw::ChannelID_t const channel) const
Returns a list of wires connected to the specified TPC channel.
double Length() const
Length is associated with z coordinate [cm].
Definition: TPCGeo.h:115
T abs(T value)
CryostatGeo const & Cryostat(geo::CryostatID const &cryoid) const
Returns the specified cryostat.
art::ServiceHandle< geo::Geometry > fGeom
bool WireIDsIntersect(WireID const &wid1, WireID const &wid2, geo::Point_t &intersection) const
Computes the intersection between two wires.
const TPCGeo & TPC(unsigned int itpc) const
Return the itpc&#39;th TPC in the cryostat.
Definition: CryostatGeo.cxx:93
double y
y position of intersection
Definition: geo_types.h:804
MaybeLogger_< ELseverityLevel::ELsev_warning, false > LogWarning
bool APAChannelsIntersect(uint32_t chan1, uint32_t chan2, std::vector< geo::WireIDIntersection > &IntersectVector)
If the channels intersect, get all intersections.
void GetCenter(double *xyz, double localz=0.0) const
Fills the world coordinate of a point on the wire.
Definition: WireGeo.cxx:73
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33
QTextStream & endl(QTextStream &s)
WireGeo const & WireIDToWireGeo(geo::WireID const &wireid) const

Member Data Documentation

unsigned int dune::apa::APAGeometryAlg::fAPAsPerCryo
private

Definition at line 96 of file APAGeometryAlg.h.

double dune::apa::APAGeometryAlg::fChannelRange[2]
private

Definition at line 108 of file APAGeometryAlg.h.

unsigned int dune::apa::APAGeometryAlg::fChannelsPerAPA
private

All APAs have this same number of channels.

Definition at line 95 of file APAGeometryAlg.h.

uint32_t dune::apa::APAGeometryAlg::fFirstU
private

Definition at line 99 of file APAGeometryAlg.h.

uint32_t dune::apa::APAGeometryAlg::fFirstV
private

Definition at line 101 of file APAGeometryAlg.h.

uint32_t dune::apa::APAGeometryAlg::fFirstZ0
private

Definition at line 103 of file APAGeometryAlg.h.

uint32_t dune::apa::APAGeometryAlg::fFirstZ1
private

Definition at line 105 of file APAGeometryAlg.h.

art::ServiceHandle<geo::Geometry> dune::apa::APAGeometryAlg::fGeom
private

Definition at line 88 of file APAGeometryAlg.h.

uint32_t dune::apa::APAGeometryAlg::fLastU
private

Definition at line 100 of file APAGeometryAlg.h.

uint32_t dune::apa::APAGeometryAlg::fLastV
private

Definition at line 102 of file APAGeometryAlg.h.

uint32_t dune::apa::APAGeometryAlg::fLastZ0
private

Definition at line 104 of file APAGeometryAlg.h.

uint32_t dune::apa::APAGeometryAlg::fLastZ1
private

Definition at line 106 of file APAGeometryAlg.h.


The documentation for this class was generated from the following files: