← Back to branch summary

~durga/maus/rel709

« back to all changes in this revision

Viewing changes to src/common_cpp/Recon/SciFi/SciFiClusterRec.cc

  • Committer: Durga Rajaram
  • Date: 2013-08-27 04:36:50 UTC
  • mfrom: (659.1.73 rc)
  • Revision ID: durga@fnal.gov-20130827043650-me0hgsbzlzikdoik
Tags: MAUS-v0.7.0
MAUS-v0.7.0

Show diffs side-by-side

added added

removed removed

Lines of Context:
src/common_cpp/Recon/SciFi/SciFiClusterRec.cc
15
15
 *
16
16
 */
17
17
#include "src/common_cpp/Recon/SciFi/SciFiClusterRec.hh"
18
 
#include <algorithm>
19
 
 
20
 
#include "src/common_cpp/Utils/ThreeVectorUtils.hh"
21
 
 
22
 
#include "Geant4/G4ThreeVector.hh"
23
 
#include "Geant4/G4RotationMatrix.hh"
24
 
 
25
18
 
26
19
namespace MAUS {
27
20
 
28
 
SciFiClusterRec::SciFiClusterRec() {}
 
21
SciFiClusterRec::SciFiClusterRec():_size_exception(0),
 
22
                                   _min_npe(0.) {}
29
23
 
30
 
SciFiClusterRec::SciFiClusterRec(int cluster_exception, double min_npe,
31
 
                                 std::vector<const MiceModule*> modules)
32
 
                                 :_size_exception(cluster_exception),
33
 
                                  _min_npe(min_npe),
34
 
                                  _modules(modules) {}
 
24
SciFiClusterRec::SciFiClusterRec(int cluster_exception,
 
25
                                 double min_npe,
 
26
                                 const std::map<int, SciFiPlaneGeometry> &geometry_map):
 
27
                                   _size_exception(cluster_exception),
 
28
                                   _min_npe(min_npe),
 
29
                                   _geometry_map(geometry_map) {}
35
30
 
36
31
SciFiClusterRec::~SciFiClusterRec() {}
37
32
 
57
52
 
58
53
std::vector<SciFiDigit*> SciFiClusterRec::get_seeds(SciFiEvent &evt) {
59
54
  std::vector<SciFiDigit*> seeds_in_event;
60
 
  for ( unsigned int dig = 0; dig < evt.digits().size(); dig++ ) {
 
55
  for ( size_t dig = 0; dig < evt.digits().size(); ++dig ) {
61
56
    if ( evt.digits()[dig]->get_npe() > _min_npe/2.0 )
62
57
      seeds_in_event.push_back(evt.digits()[dig]);
63
58
  }
64
59
  return seeds_in_event;
65
60
}
66
61
 
67
 
void SciFiClusterRec::make_clusters(SciFiEvent &evt, std::vector<SciFiDigit*>   &seeds) {
68
 
  int seeds_size = seeds.size();
69
 
  for ( int i = 0; i < seeds_size; i++ ) {
70
 
    if ( !seeds[i]->is_used() ) {
 
62
void SciFiClusterRec::make_clusters(SciFiEvent &evt, std::vector<SciFiDigit*> &seeds) {
 
63
  size_t seeds_size = seeds.size();
 
64
  for ( size_t i = 0; i < seeds_size; i++ ) {
 
65
    if ( !(seeds[i]->is_used()) ) {
71
66
      SciFiDigit* neigh = NULL;
72
67
      SciFiDigit* seed = seeds[i];
73
68
 
74
69
      double pe = seed->get_npe();
75
70
      // Look for a neighbour.
76
 
      for ( int j = i+1; j < seeds_size; j++ ) {
 
71
      for ( size_t j = i+1; j < seeds_size; ++j ) {
77
72
        if ( are_neighbours(seeds[i], seeds[j]) ) {
78
73
          neigh = seeds[j];
79
74
        }
96
91
}
97
92
 
98
93
void SciFiClusterRec::process_cluster(SciFiCluster *clust) {
99
 
  // Get the MiceModule of the plane...
100
94
  int tracker = clust->get_tracker();
101
95
  int station = clust->get_station();
102
96
  int plane   = clust->get_plane();
103
 
  const MiceModule* this_plane = NULL;
104
 
  this_plane = find_plane(tracker, station, plane);
105
 
  assert(this_plane != NULL);
106
 
  // compute it's direction & position in TRF...
107
 
  ThreeVector trf_dir(0., 1., 0.);
108
 
  ThreeVector trf_pos(0., 0., 0.);
109
 
  double alpha;
110
 
  construct(clust, this_plane, trf_dir, trf_pos, alpha);
111
 
 
112
 
  clust->set_direction(trf_dir);
113
 
  clust->set_position(trf_pos);
114
 
 
115
 
  clust->set_alpha(alpha);
116
 
  int id = 15*tracker + 3*(station-1) + (plane);
 
97
 
 
98
  int id =  3*(station-1) + (plane+1);
 
99
  id = ( tracker == 0 ? -id : id );
117
100
  clust->set_id(id);
118
 
  /*
119
 
   std::cerr << "----------Clustering--------- \n"
120
 
            << "Site ID: " << id << "\n"
121
 
            << "Tracker " << tracker << ", station " << station << ", plane " << plane << "\n"
122
 
            << "Fibre direction: " << dir << "\n"
123
 
            << "Position: " << tracker_ref_frame_pos << "\n";
124
 
  */
125
 
}
126
 
 
127
 
void SciFiClusterRec::construct(SciFiCluster *clust,
128
 
                                const MiceModule* this_plane,
129
 
                                ThreeVector &dir,
130
 
                                ThreeVector &tracker_ref_frame_pos,
131
 
                                double &alpha) {
132
 
  ThreeVector perp(-1., 0., 0.);
133
 
 
134
 
  CLHEP::HepRotation zflip;
135
 
  const Hep3Vector rowx(-1., 0., 0.);
136
 
  const Hep3Vector rowy(0., 1., 0.);
137
 
  const Hep3Vector rowz(0., 0., -1.);
138
 
  zflip.setRows(rowx, rowy, rowz);
139
 
  G4RotationMatrix trot(this_plane->globalRotation());
140
 
  // Rotations of the planes in the Tracker Reference Frame.
141
 
  if ( clust->get_tracker() == 0 ) {
142
 
    trot= trot*zflip;
143
 
    dir  *= trot;
144
 
    perp *= trot;
145
 
  } else if ( clust->get_tracker() == 1 ) {
146
 
    dir  *= trot;
147
 
    perp *= trot;
 
101
 
 
102
  std::map<int, SciFiPlaneGeometry>::iterator iterator;
 
103
  iterator = _geometry_map.find(id);
 
104
  // Throw if the plane isn't found.
 
105
  if ( iterator == _geometry_map.end() ) {
 
106
    throw(Squeal(Squeal::nonRecoverable,
 
107
    "Failed to find SciFi plane in _geometry_map.",
 
108
    "SciFiClusterRec::process_cluster"));
148
109
  }
 
110
  SciFiPlaneGeometry this_plane = (*iterator).second;
 
111
  ThreeVector plane_direction   = this_plane.Direction;
 
112
  ThreeVector plane_position    = this_plane.Position;
 
113
  double Pitch                  = this_plane.Pitch;
 
114
  double CentralFibre           = this_plane.CentralFibre;
 
115
  // alpha is the distance to the central fibre.
 
116
  double alpha   = clust->get_channel()-CentralFibre;
 
117
  double dist_mm = Pitch * 7.0 / 2.0 * alpha;
149
118
 
150
 
  double Pitch = this_plane->propertyDouble("Pitch");
151
 
  double CentralFibre = this_plane->propertyDouble("CentralFibre");
152
 
  double dist_mm = Pitch * 7.0 / 2.0 * (clust->get_channel() - CentralFibre);
153
 
  ThreeVector plane_position = clhep_to_root(this_plane->globalPosition());
 
119
  ThreeVector perp = plane_direction.Orthogonal();
154
120
  ThreeVector position = dist_mm * perp + plane_position;
155
 
  ThreeVector reference = get_reference_frame_pos(clust->get_tracker());
156
 
 
157
 
  tracker_ref_frame_pos = position - reference;
158
 
 
159
 
  // ThreeVector tracker_ref_frame_pos;
160
 
  if ( clust->get_tracker() == 0 ) {
161
 
    tracker_ref_frame_pos = - (position - reference);
162
 
  } else if ( clust->get_tracker() == 1 ) {
163
 
    tracker_ref_frame_pos = position - reference;
164
 
  }
165
 
 
166
 
  alpha = clust->get_channel() - CentralFibre;
167
 
  if ( clust->get_tracker() == 1 ) {
168
 
    alpha = -alpha;
169
 
  }
170
 
}
171
 
 
172
 
const MiceModule* SciFiClusterRec::find_plane(int tracker, int station, int plane) {
173
 
  const MiceModule* scifi_plane = NULL;
174
 
  for ( unsigned int j = 0; !scifi_plane && j < _modules.size(); j++ ) {
175
 
    // Find the right module
176
 
    if ( _modules[j]->propertyExists("Tracker", "int") &&
177
 
         _modules[j]->propertyExists("Station", "int") &&
178
 
         _modules[j]->propertyExists("Plane", "int")  &&
179
 
         _modules[j]->propertyInt("Tracker") ==
180
 
         tracker &&
181
 
         _modules[j]->propertyInt("Station") ==
182
 
         station &&
183
 
         _modules[j]->propertyInt("Plane") ==
184
 
         plane ) {
185
 
         // Save the module
186
 
      scifi_plane = _modules[j];
187
 
    }
188
 
  }
189
 
  return scifi_plane;
190
 
}
191
 
 
192
 
ThreeVector SciFiClusterRec::get_reference_frame_pos(int tracker) {
193
 
  // Reference plane is plane 0, station 1 of current tracker.
194
 
  int station = 1;
195
 
  int plane   = 0;
196
 
  const MiceModule* reference_plane = NULL;
197
 
  reference_plane = find_plane(tracker, station, plane);
198
 
 
199
 
  assert(reference_plane != NULL);
200
 
  ThreeVector reference_pos = clhep_to_root(reference_plane->globalPosition());
201
 
 
202
 
  return reference_pos;
 
121
 
 
122
  clust->set_direction(plane_direction);
 
123
  clust->set_position(position);
 
124
  clust->set_alpha(alpha);
203
125
}
204
126
 
205
127
bool SciFiClusterRec::are_neighbours(SciFiDigit *seed_i, SciFiDigit *seed_j) {
206
128
  bool neigh = false;
207
129
 
208
 
  if ( !seed_j->is_used() && // seed is unused
 
130
  if ( !(seed_j->is_used()) && // seed is unused
209
131
       seed_j->get_spill() == seed_i->get_spill() && // same spill
210
132
       seed_j->get_event() == seed_i->get_event() && // same event
211
133
       seed_j->get_tracker() == seed_i->get_tracker() && // same tracker