All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Public Member Functions | Private Member Functions | Private Attributes | List of all members
wcls::SimChannelSink Class Reference

#include <SimChannelSink.h>

Inheritance diagram for wcls::SimChannelSink:
wcls::IArtEventVisitor

Public Member Functions

 SimChannelSink ()
 
virtual void produces (art::ProducesCollector &collector)
 IArtEventVisitor. More...
 
virtual void visit (art::Event &event)
 Implement to visit an Art event. More...
 
virtual bool operator() (const WireCell::IDepo::pointer &indepo, WireCell::IDepo::pointer &outdepo)
 IDepoFilter. More...
 
virtual WireCell::Configuration default_configuration () const
 IConfigurable. More...
 
virtual void configure (const WireCell::Configuration &config)
 
- Public Member Functions inherited from wcls::IArtEventVisitor
virtual ~IArtEventVisitor ()
 

Private Member Functions

void save_as_simchannel (const WireCell::IDepo::pointer &depo)
 

Private Attributes

WireCell::IDepo::pointer m_depo
 
std::vector
< WireCell::IAnodePlane::pointer > 
m_anodes
 
WireCell::IRandom::pointer m_rng
 
std::map< unsigned int,
sim::SimChannel
m_mapSC
 
std::string m_artlabel
 
double m_readout_time
 
double m_tick
 
double m_start_time
 
double m_nsigma
 
double m_drift_speed
 
double m_u_to_rp
 
double m_v_to_rp
 
double m_y_to_rp
 
double m_u_time_offset
 
double m_v_time_offset
 
double m_y_time_offset
 
double m_g4_ref_time
 
bool m_use_energy
 
bool m_use_extra_sigma
 

Detailed Description

Definition at line 23 of file SimChannelSink.h.

Constructor & Destructor Documentation

SimChannelSink::SimChannelSink ( )

Definition at line 17 of file SimChannelSink.cxx.

18  : m_depo(nullptr)
19 {
20  m_mapSC.clear();
21 }
WireCell::IDepo::pointer m_depo
std::map< unsigned int, sim::SimChannel > m_mapSC

Member Function Documentation

void SimChannelSink::configure ( const WireCell::Configuration &  config)
virtual

Definition at line 47 of file SimChannelSink.cxx.

48 {
49  m_anodes.clear();
50  auto anodes_tn = cfg["anodes_tn"];
51  for (auto anode_tn: anodes_tn) {
52  auto anode = Factory::find_tn<IAnodePlane>(anode_tn.asString());
53  m_anodes.push_back(anode);
54  }
55 
56  // Initialize a SimChannel map. Two assumptions are followed in
57  // larsim::BackTracker::FindCimChannel()
58  // 1) fill all channels even with empty SimChannel
59  // 2) SimChannels are sorted in channel number
60  // In SimChannelSink::visit(), the map is reinitialized insted of
61  // calling map::clear()
62  for (auto& anode: m_anodes){
63  for (auto& channel: anode->channels()){
64  // std::cout << "-> channel: " << channel << std::endl;
65  m_mapSC.emplace(channel, sim::SimChannel(channel));
66  }
67  }
68  // std::cout << "m_mapSC.size(): " << m_mapSC.size() << std::endl;
69 
70  if (m_anodes.empty()) {
71  const std::string anode_tn = cfg["anode"].asString();
72  if (anode_tn.empty()) {
73  THROW(ValueError() << errmsg{"SimChannelSink requires an anode plane"});
74  }
75  auto anode = Factory::find_tn<IAnodePlane>(anode_tn);
76  m_anodes.push_back(anode);
77  }
78 
79  const std::string rng_tn = cfg["rng"].asString();
80  if (rng_tn.empty()) {
81  THROW(ValueError() << errmsg{"SimChannelSink requires a noise source"});
82  }
83  m_rng = Factory::find_tn<IRandom>(rng_tn);
84 
85  m_artlabel = cfg["artlabel"].asString();
86  m_artlabel = get(cfg,"artlabel",m_artlabel);
87  m_start_time = get(cfg,"start_time",-1.6*units::ms);
88  m_readout_time = get(cfg,"readout_time",4.8*units::ms);
89  m_tick = get(cfg,"tick",0.5*units::us);
90  m_nsigma = get(cfg,"nsigma",3.0);
91  m_drift_speed = get(cfg,"drift_speed",1.098*units::mm/units::us);
92  m_u_to_rp = get(cfg,"uboone_u_to_rp",94*units::mm); // compatible interface
93  m_v_to_rp = get(cfg,"uboone_v_to_rp",97*units::mm);
94  m_y_to_rp = get(cfg,"uboone_y_to_rp",100*units::mm);
95  // compatible interface for protoDUNE-SP
96  if (cfg.isMember("u_to_rp")) m_u_to_rp = get(cfg,"u_to_rp",90.58*units::mm);
97  if (cfg.isMember("v_to_rp")) m_v_to_rp = get(cfg,"v_to_rp",95.29*units::mm);
98  if (cfg.isMember("y_to_rp")) m_y_to_rp = get(cfg,"y_to_rp",100.0*units::mm);
99 
100  m_u_time_offset = get(cfg,"u_time_offset",0.0*units::us);
101  m_v_time_offset = get(cfg,"v_time_offset",0.0*units::us);
102  m_y_time_offset = get(cfg,"y_time_offset",0.0*units::us);
103  m_g4_ref_time = get(cfg,"g4_ref_time",-4050*units::us);
104  m_use_energy = get(cfg,"use_energy",false);
105  m_use_extra_sigma = get(cfg,"use_extra_sigma",false);
106 
107 }
Energy deposited on a readout channel by simulated tracks.
Definition: SimChannel.h:145
std::string m_artlabel
WireCell::IRandom::pointer m_rng
std::vector< WireCell::IAnodePlane::pointer > m_anodes
std::map< unsigned int, sim::SimChannel > m_mapSC
WireCell::Configuration SimChannelSink::default_configuration ( ) const
virtual

IConfigurable.

Definition at line 23 of file SimChannelSink.cxx.

24 {
25  Configuration cfg;
26  cfg["anodes_tn"] = Json::arrayValue;
27  cfg["anode"] = "AnodePlane"; // either `anodes_tn` or `anode`
28  cfg["rng"] = "Random";
29  cfg["artlabel"] = "simpleSC";
30  cfg["start_time"] = -1.6*units::ms;
31  cfg["readout_time"] = 4.8*units::ms;
32  cfg["tick"] = 0.5*units::us;
33  cfg["nsigma"] = 3.0;
34  cfg["drift_speed"] = 1.098*units::mm/units::us;
35  cfg["uboone_u_to_rp"] = 94*units::mm;
36  cfg["uboone_v_to_rp"] = 97*units::mm;
37  cfg["uboone_y_to_rp"] = 100*units::mm;
38  cfg["u_time_offset"] = 0.0*units::us;
39  cfg["v_time_offset"] = 0.0*units::us;
40  cfg["y_time_offset"] = 0.0*units::us;
41  cfg["g4_ref_time"] = -4050*units::us; // uboone: -4050us, pdsp: -250us
42  cfg["use_energy"] = false;
43  cfg["use_extra_sigma"] = false;
44  return cfg;
45 }
bool SimChannelSink::operator() ( const WireCell::IDepo::pointer &  indepo,
WireCell::IDepo::pointer &  outdepo 
)
virtual

IDepoFilter.

Definition at line 284 of file SimChannelSink.cxx.

286 {
287  outdepo = indepo;
288  if (indepo) {
289  m_depo = indepo;
291  }
292 
293  return true;
294 }
void save_as_simchannel(const WireCell::IDepo::pointer &depo)
WireCell::IDepo::pointer m_depo
void SimChannelSink::produces ( art::ProducesCollector &  collector)
virtual

IArtEventVisitor.

Reimplemented from wcls::IArtEventVisitor.

Definition at line 109 of file SimChannelSink.cxx.

110 {
111  collector.produces< std::vector<sim::SimChannel> >(m_artlabel);
112 }
std::string m_artlabel
void SimChannelSink::save_as_simchannel ( const WireCell::IDepo::pointer &  depo)
private

Definition at line 114 of file SimChannelSink.cxx.

114  {
115  // Binning tbins(m_readout_time/m_tick, m_start_time, m_start_time+m_readout_time);
116 
117  /* Start the gate ealier for the depos between the response
118  * plane and the anode plane. Those depos are anti-drifted
119  * to the reponse plane, so the start time is earlier.
120  * c.f. jsonnet config in wirecell toolkit: params.sim.ductor
121  */
122  double response_plane = 10.0 * units::cm;
123  double response_time_offset = response_plane / m_drift_speed;
124  int response_nticks = (int)(response_time_offset / m_tick);
125  Binning tbins(m_readout_time/m_tick + response_nticks, m_start_time - response_time_offset, m_start_time+m_readout_time);
126 
127  if(!depo) return;
128 
129  int ctr = 0;
130  while(ctr<1){
131  ctr++;
132  // if(ctr % 10000==0){ std::cout<<"COUNTER "<<ctr<<std::endl;}
133  for (auto anode: m_anodes) {
134  for(auto face : anode->faces()){
135  auto boundbox = face->sensitive();
136  if(!boundbox.inside(depo->pos())) continue;
137 
138  // int plane = -1;
139  // for(Pimpos* pimpos : {uboone_u, uboone_v, uboone_y}){
140  for (auto plane : face->planes()) {
141  // plane++;
142  int iplane = plane->planeid().index();
143  if (iplane<0) continue;
144  const Pimpos* pimpos = plane->pimpos();
145  auto& wires = plane->wires();
146 
147  const double center_time = depo->time();
148  const double center_pitch = pimpos->distance(depo->pos());
149 
150  double sigma_L = depo->extent_long();
151  if (m_use_extra_sigma) {
152  int nrebin=1;
153  double time_slice_width = nrebin * m_drift_speed * m_tick; // units::mm
154  double add_sigma_L = 1.428249 * time_slice_width / nrebin / (m_tick/units::us); // units::mm
155  sigma_L = sqrt( pow(depo->extent_long(),2) + pow(add_sigma_L,2) );// / time_slice_width;
156  }
157  Gen::GausDesc time_desc(center_time, sigma_L / m_drift_speed);
158  // Gen::GausDesc time_desc(center_time, depo->extent_long() / m_drift_speed);
159  {
160  double nmin_sigma = time_desc.distance(tbins.min());
161  double nmax_sigma = time_desc.distance(tbins.max());
162 
163  double eff_nsigma = depo->extent_long() / m_drift_speed>0?m_nsigma:0;
164  if (nmin_sigma > eff_nsigma || nmax_sigma < -eff_nsigma) {
165  break;
166  }
167  }
168 
169  // auto ibins = pimpos->impact_binning();
170  auto wbins = pimpos->region_binning(); // wire binning
171 
172  double sigma_T = depo->extent_tran();
173  if (m_use_extra_sigma) {
174  double add_sigma_T = wbins.binsize();
175  if (iplane==0) add_sigma_T *= (0.402993*0.3);
176  else if (iplane==1) add_sigma_T *= (0.402993*0.5);
177  else if (iplane==2) add_sigma_T *= (0.188060*0.2);
178  sigma_T = sqrt( pow(depo->extent_tran(),2) + pow(add_sigma_T,2) ); // / wbins.binsize();
179  }
180  Gen::GausDesc pitch_desc(center_pitch, sigma_T);
181  // Gen::GausDesc pitch_desc(center_pitch, depo->extent_tran());
182  {
183  double nmin_sigma = pitch_desc.distance(wbins.min());
184  double nmax_sigma = pitch_desc.distance(wbins.max());
185 
186  double eff_nsigma = depo->extent_tran()>0?m_nsigma:0;
187  if (nmin_sigma > eff_nsigma || nmax_sigma < -eff_nsigma) {
188  break;
189  }
190  }
191 
192  auto gd = std::make_shared<Gen::GaussianDiffusion>(depo, time_desc, pitch_desc);
193  gd->set_sampling(tbins, wbins, m_nsigma, 0, 1);
194 
195  double xyz[3];
196  int id = -10000;
197  double energy = 100.0;
198  if(depo->prior()){
199  id = depo->prior()->id();
200  if(m_use_energy){ energy = depo->prior()->energy(); }
201  }
202  else{
203  id = depo->id();
204  if(m_use_energy){ energy = depo->energy(); }
205  }
206 
207  const auto patch = gd->patch();
208  const int poffset_bin = gd->poffset_bin();
209  const int toffset_bin = gd->toffset_bin();
210  const int np = patch.rows();
211  const int nt = patch.cols();
212 
213  int min_imp = 0;
214  int max_imp = wbins.nbins();
215 
216  for (int pbin = 0; pbin != np; pbin++){
217  int abs_pbin = pbin + poffset_bin;
218  if (abs_pbin < min_imp || abs_pbin >= max_imp) continue;
219 
220  auto iwire = wires[abs_pbin];
221  int channel = iwire->channel();
222  // int channel = abs_pbin;
223  // if(plane == 1){ channel = abs_pbin+2400; }
224  // if(plane == 2){ channel = abs_pbin+4800; }
225 
226  auto channelData = m_mapSC.find(channel);
227  sim::SimChannel& sc = (channelData == m_mapSC.end())
228  ? (m_mapSC[channel]=sim::SimChannel(channel))
229  : channelData->second;
230 
231  for (int tbin = 0; tbin!= nt; tbin++){
232  int abs_tbin = tbin + toffset_bin;
233  double charge = patch(pbin, tbin);
234  double tdc = tbins.center(abs_tbin);
235 
236  // double wire_response_offset = iwire->center().x() - pimpos->origin().x();
237  if(iplane == 0){
238  tdc = tdc + (m_u_to_rp/m_drift_speed) + m_u_time_offset;
239  // xyz[0] = depo->pos().x()/units::cm - 94*units::mm/units::cm; // m_u_to_rp/units::cm;
240  }
241  if(iplane == 1){
242  tdc = tdc + (m_v_to_rp/m_drift_speed) + m_v_time_offset;
243  // xyz[0] = depo->pos().x()/units::cm - 97*units::mm/units::cm; // m_v_to_rp/units::cm;
244  }
245  if(iplane == 2){
246  tdc = tdc + (m_y_to_rp/m_drift_speed) + m_y_time_offset;
247  // xyz[0] = depo->pos().x()/units::cm - 100*units::mm/units::cm; // m_y_to_rp/units::cm;
248  }
249  // xyz[0] = depo->pos().x()/units::cm + wire_response_offset/units::cm;
250  WireCell::IDepo::pointer orig = depo_chain(depo).back(); // first depo in the chain
251  xyz[0] = orig->pos().x()/units::cm;
252  xyz[1] = orig->pos().y()/units::cm;
253  xyz[2] = orig->pos().z()/units::cm;
254 
255  unsigned int temp_time = (unsigned int) ( (tdc - m_g4_ref_time) / m_tick );
256  charge = abs(charge);
257  if(charge>1){
258  sc.AddIonizationElectrons(id, temp_time, charge, xyz, energy*abs(charge/depo->charge()));
259  }
260  }
261  }
262  } // plane
263  } //face
264  } // anode
265  }
266 }
Energy deposited on a readout channel by simulated tracks.
Definition: SimChannel.h:145
T abs(T value)
for($it=0;$it< $RaceTrack_number;$it++)
std::vector< WireCell::IAnodePlane::pointer > m_anodes
std::map< unsigned int, sim::SimChannel > m_mapSC
void AddIonizationElectrons(TrackID_t trackID, TDC_t tdc, double numberElectrons, double const *xyz, double energy, TrackID_t origTrackID=util::kBogusI)
Add ionization electrons and energy to this channel.
Definition: SimChannel.cxx:55
void SimChannelSink::visit ( art::Event &  event)
virtual

Implement to visit an Art event.

Implements wcls::IArtEventVisitor.

Definition at line 268 of file SimChannelSink.cxx.

269 {
270  std::unique_ptr<std::vector<sim::SimChannel> > out(new std::vector<sim::SimChannel>);
271 
272  for(auto& m : m_mapSC){
273  out->emplace_back(m.second);
274  }
275 
276  event.put(std::move(out), m_artlabel);
277  // m_mapSC.clear();
278  for (auto& elem: m_mapSC){
279  elem.second = sim::SimChannel(elem.first);
280  }
281 
282 }
Energy deposited on a readout channel by simulated tracks.
Definition: SimChannel.h:145
std::string m_artlabel
tuple m
now if test mode generate materials, CRT shell, world, gdml header else just generate CRT shell for u...
std::map< unsigned int, sim::SimChannel > m_mapSC

Member Data Documentation

std::vector<WireCell::IAnodePlane::pointer> wcls::SimChannelSink::m_anodes
private

Definition at line 46 of file SimChannelSink.h.

std::string wcls::SimChannelSink::m_artlabel
private

Definition at line 53 of file SimChannelSink.h.

WireCell::IDepo::pointer wcls::SimChannelSink::m_depo
private

Definition at line 44 of file SimChannelSink.h.

double wcls::SimChannelSink::m_drift_speed
private

Definition at line 58 of file SimChannelSink.h.

double wcls::SimChannelSink::m_g4_ref_time
private

Definition at line 65 of file SimChannelSink.h.

std::map<unsigned int,sim::SimChannel> wcls::SimChannelSink::m_mapSC
private

Definition at line 49 of file SimChannelSink.h.

double wcls::SimChannelSink::m_nsigma
private

Definition at line 57 of file SimChannelSink.h.

double wcls::SimChannelSink::m_readout_time
private

Definition at line 54 of file SimChannelSink.h.

WireCell::IRandom::pointer wcls::SimChannelSink::m_rng
private

Definition at line 47 of file SimChannelSink.h.

double wcls::SimChannelSink::m_start_time
private

Definition at line 56 of file SimChannelSink.h.

double wcls::SimChannelSink::m_tick
private

Definition at line 55 of file SimChannelSink.h.

double wcls::SimChannelSink::m_u_time_offset
private

Definition at line 62 of file SimChannelSink.h.

double wcls::SimChannelSink::m_u_to_rp
private

Definition at line 59 of file SimChannelSink.h.

bool wcls::SimChannelSink::m_use_energy
private

Definition at line 66 of file SimChannelSink.h.

bool wcls::SimChannelSink::m_use_extra_sigma
private

Definition at line 67 of file SimChannelSink.h.

double wcls::SimChannelSink::m_v_time_offset
private

Definition at line 63 of file SimChannelSink.h.

double wcls::SimChannelSink::m_v_to_rp
private

Definition at line 60 of file SimChannelSink.h.

double wcls::SimChannelSink::m_y_time_offset
private

Definition at line 64 of file SimChannelSink.h.

double wcls::SimChannelSink::m_y_to_rp
private

Definition at line 61 of file SimChannelSink.h.


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