test_GeometryDune35t.cxx
Go to the documentation of this file.
1 // test_GeometryDune35t.cxx
2 
3 // David Adams
4 // OCtober 2016
5 //
6 // Test the DUNE 35t geometry including channel mapping.
7 
8 #include "test_GeometryDune.cxx"
9 
10 //**********************************************************************
11 
12 void setExpectedValues(ExpectedValues& ev) {
13  ev.gname ="dune35t_geo";
14  ev.fullname = "dune35t4apa_v6";
15  // Geometry counts.
16  ev.ncry = 1;
17  ev.ntpc = 8;
18  ev.npla = 3;
19  ev.napa = 4;
20  ev.nrop = 4;
21  // Signal type and view for each TPC plane.
22  resize(ev.view, ev.ncry, ev.ntpc, ev.npla, geo::kUnknown);
23  resize(ev.sigType, ev.ncry, ev.ntpc, ev.npla, geo::kMysteryType);
24  for ( Index icry=0; icry<ev.ncry; ++icry ) {
25  for ( Index itpc=0; itpc<ev.ntpc; ++itpc ) {
26  for ( Index ipla=0; ipla<ev.npla; ++ipla ) {
27  ev.sigType[icry][itpc][ipla] = (ipla < 2) ? geo::kInduction : geo::kCollection;
28  }
29  ev.view[icry][itpc][0] = geo::kU;
30  ev.view[icry][itpc][1] = geo::kV;
31  ev.view[icry][itpc][2] = geo::kZ;
32  }
33  }
34  // Wire count for each TPC plane.
35  resize(ev.nwirPerPlane, ev.ntpc, ev.npla, 0);
36  ev.nwirPerPlane[0][0] = 359;
37  ev.nwirPerPlane[0][1] = 345;
38  ev.nwirPerPlane[0][2] = 112;
39  ev.nwirPerPlane[2][0] = 194;
40  ev.nwirPerPlane[2][1] = 188;
41  ev.nwirPerPlane[2][2] = 112;
42  ev.nwirPerPlane[4][0] = 236;
43  ev.nwirPerPlane[4][1] = 228;
44  ev.nwirPerPlane[4][2] = 112;
45  for ( Index ipla=0; ipla<ev.npla; ++ipla ) ev.nwirPerPlane[1][ipla] = ev.nwirPerPlane[0][ipla];
46  for ( Index ipla=0; ipla<ev.npla; ++ipla ) ev.nwirPerPlane[3][ipla] = ev.nwirPerPlane[2][ipla];
47  for ( Index ipla=0; ipla<ev.npla; ++ipla ) ev.nwirPerPlane[5][ipla] = ev.nwirPerPlane[4][ipla];
48  for ( Index ipla=0; ipla<ev.npla; ++ipla ) ev.nwirPerPlane[6][ipla] = ev.nwirPerPlane[0][ipla];
49  for ( Index ipla=0; ipla<ev.npla; ++ipla ) ev.nwirPerPlane[7][ipla] = ev.nwirPerPlane[0][ipla];
50  resize(ev.nchaPerRop, ev.nrop, 0);
51  ev.nchaPerRop[0] = 144;
52  ev.nchaPerRop[1] = 144;
53  ev.nchaPerRop[2] = 112;
54  ev.nchaPerRop[3] = 112;
55  for ( Index irop=0; irop<ev.nrop; ++irop ) ev.nchaPerApa += ev.nchaPerRop[irop];
56  ev.nchatot = ev.napa*ev.nchaPerApa;
57  resize(ev.chacry, ev.nchatot, InvalidIndex);
58  resize(ev.chaapa, ev.nchatot, InvalidIndex);
59  resize(ev.charop, ev.nchatot, InvalidIndex);
60  Index icha = 0;
61  for ( Index icry=0; icry<ev.ncry; ++icry ) {
62  for ( Index iapa=0; iapa<ev.napa; ++iapa ) {
63  for ( Index irop=0; irop<ev.nrop; ++irop ) {
64  for ( Index kcha=0; kcha<ev.nchaPerRop[irop]; ++kcha ) {
65  ev.chacry[icha] = icry;
66  ev.chaapa[icha] = iapa;
67  ev.charop[icha] = irop;
68  ++icha;
69  }
70  }
71  }
72  }
73  resize(ev.firstchan, ev.ncry, ev.napa, ev.nrop, raw::InvalidChannelID);
74  Index chan = 0;
75  for ( Index icry=0; icry<ev.ncry; ++icry ) {
76  for ( Index iapa=0; iapa<ev.napa; ++iapa ) {
77  for ( Index irop=0; irop<ev.nrop; ++irop ) {
78  ev.firstchan[icry][iapa][irop] = chan;
79  chan += ev.nchaPerRop[irop];
80  }
81  }
82  }
83  #include "set35tSpacePoints.dat"
84  // Optical detectors.
85  ev.nopdet = 8;
86  ev.nopdetcha.push_back( 8);
87  ev.nopdetcha.push_back(12);
88  ev.nopdetcha.push_back( 2);
89  ev.nopdetcha.push_back(12);
90  ev.nopdetcha.push_back( 8);
91  ev.nopdetcha.push_back(12);
92  ev.nopdetcha.push_back( 8);
93  ev.nopdetcha.push_back(12);
94  ev.nopcha = ev.nopdet*12;
95  resize(ev.opdetcha, ev.nopdet, 12, 0);
96  for ( Index iopt=0; iopt<ev.nopdet; ++iopt ) {
97  Index icha = 12*iopt;
98  for ( Index ioch=0; ioch<ev.nopdetcha[iopt]; ++ ioch ) {
99  ev.opdetcha[iopt][ioch] = icha++;
100  }
101  }
102 }
103 
104 //**********************************************************************
105 
107  const string myname = "setExpectedValuespacePoints: ";
108  string ofname = "set35tSpacePoints.dat";
109  cout << myname << "Writing " << ofname << endl;
110  const CryostatGeo& crygeo = pgeo->Cryostat(0);
111  double origin[3] = {0.0};
112  double crypos[3] = {0.0};
113  crygeo.LocalToWorld(origin, crypos);
114  double cxlo = crypos[0] - crygeo.HalfWidth();
115  double cxhi = crypos[0] + crygeo.HalfWidth();
116  double cylo = crypos[1] - crygeo.HalfHeight();
117  double cyhi = crypos[1] + crygeo.HalfHeight();
118  double czlo = crypos[2] - 0.5*crygeo.Length();
119  double czhi = crypos[2] + 0.5*crygeo.Length();
120  cout << "Cryostat limits: "
121  << "(" << cxlo << ", " << cylo << ", " << czlo << "), "
122  << "(" << cxhi << ", " << cyhi << ", " << czhi << ")" << endl;
123  vector<double> zfs = {0.2, 0.3, 0.4, 0.5, 0.6, 0.68 };
124  vector<double> yfs = {0.2, 0.5, 0.8};
125  vector<double> xfs = {0.1, 0.3, 0.6 };
126  ofstream fout(ofname.c_str());
127  for ( double zf : zfs ) {
128  double z = czlo + zf*(czhi-czlo);
129  for ( double yf : yfs ) {
130  double y = cylo + yf*(cyhi-cylo);
131  for ( double xf : xfs ) {
132  double x = cxlo + xf*(cxhi-cxlo);
133  double xyz[3] = {x, y, z};
134  TPCID tpcid = pgeo->FindTPCAtPosition(xyz);
135  unsigned int itpc = tpcid.TPC;
136  const TPCGeo& tpcgeo = pgeo->TPC(tpcid);
137  unsigned int npla = tpcgeo.Nplanes();
138  fout << " ev.posXyz.push_back(SpacePoint(" << x << ", " << y << ", " << z << "));" << endl;
139  for ( unsigned int ipla=0; ipla<npla; ++ipla ) {
140  PlaneID plaid(tpcid, ipla);
141  double xwire = pgeo->WireCoordinate(x, y, plaid);
142  fout << " ev.posTpc.push_back(" << itpc << ");" << endl;
143  fout << " ev.posPla.push_back(" << ipla << ");" << endl;
144  fout << " ev.posWco.push_back(" << xwire << ");" << endl;
145  } // end loop over planes
146  } // end loop over x
147  } // end loop over y
148  } // end loop over z
149 }
150 
151 //**********************************************************************
geo::Length_t WireCoordinate(double YPos, double ZPos, geo::PlaneID const &planeid) const
Returns the index of the nearest wire to the specified position.
Who knows?
Definition: geo_types.h:147
Planes which measure V.
Definition: geo_types.h:130
Unknown view.
Definition: geo_types.h:136
unsigned int Nplanes() const
Number of planes in this tpc.
Definition: TPCGeo.h:165
void setExpectedValues(ExpectedValues &ev)
The data type to uniquely identify a Plane.
Definition: geo_types.h:472
Geometry information for a single TPC.
Definition: TPCGeo.h:38
Planes which measure Z direction.
Definition: geo_types.h:132
unsigned int Index
Geometry information for a single cryostat.
Definition: CryostatGeo.h:43
void resize(Vector< T > &vec1, Index n1, const V &val)
const Index InvalidIndex
geo::TPCID FindTPCAtPosition(double const worldLoc[3]) const
Returns the ID of the TPC at specified location.
Planes which measure U.
Definition: geo_types.h:129
constexpr ChannelID_t InvalidChannelID
ID of an invalid channel.
Definition: RawTypes.h:32
Signal from induction planes.
Definition: geo_types.h:145
double HalfWidth() const
Half width of the cryostat [cm].
CryostatGeo const & Cryostat(geo::CryostatID const &cryoid) const
Returns the specified cryostat.
The geometry of one entire detector, as served by art.
Definition: Geometry.h:196
The data type to uniquely identify a TPC.
Definition: geo_types.h:386
double HalfHeight() const
Half height of the cryostat [cm].
void setExpectedValuesSpacePoints(Geometry *pgeo)
void LocalToWorld(const double *cryo, double *world) const
Transform point from local cryostat frame to world frame.
Definition: CryostatGeo.h:387
for(std::string line;std::getline(inFile, line);)
Definition: regex_t.cc:37
list x
Definition: train.py:276
TPCGeo const & TPC(unsigned int const tpc=0, unsigned int const cstat=0) const
Returns the specified TPC.
TPCID_t TPC
Index of the TPC within its cryostat.
Definition: geo_types.h:406
constexpr Point origin()
Returns a origin position with a point of the specified type.
Definition: geo_vectors.h:227
double Length() const
Length of the cryostat [cm].
Definition: CryostatGeo.h:107
QTextStream & endl(QTextStream &s)
Signal from collection planes.
Definition: geo_types.h:146