WireCell / wire-cell-toolkit

Toolkit for Liquid Argon TPC Simulation and Reconstruction
https://wirecell.github.io/
Other
7 stars 22 forks source link

Simulating cathode non-planarity #329

Closed wenqiang-gu closed 2 months ago

wenqiang-gu commented 3 months ago

I was contacted by some SBN folks (M. Mooney et al.) regarding the needs to simulate cathode non-planarity. @brettviren @HaiwangYu Could you comment on the following?

Basically, they aim to model the X position of the cathode as a function of the Y and Z positions: X_cathode = f(Y,Z). In ICARUS, this non-planarity could contribute to about 10% of the systematic uncertainty in the $\nu_{\mu}$ selection(as mentioned by M. Mooney).

Since changes in X_cathode could cause some fraction of energy deposits to be assigned to a different TPC drift volume, the SBN team would like to simulate this effect. They can provide a 2D map characterizing the function X_cathode = f(Y,Z). It's also noted that the effect is more significant at the boundary of the cryostat, while at the center (small Y,Z) the effect is negligible. This could reduce some computing time.

Currently, in Drifter, the drift volume determination is based solely on the X position: https://github.com/WireCell/wire-cell-toolkit/blob/1a00d641b90d500238c4696bc204f45905c2a454/gen/inc/WireCellGen/Drifter.h#L172

It appears that we can implement a function to determine a drift volume with depo's X,Y,Z position? For example,

xr.inside_bulk(depo->pos().x(), depo()->pos().y(), depo->pos().z() );

And inside the function Xregion::inside_bulk(double x, double y, double z), we can correct the position of cathode? For example,

bool Gen::Drifter::Xregion::inside_bulk(double x, double y, double z) const
{
    double cathode = f(y,z); // SBN's 2D map 
    return (response < x and x < cathode) or (cathode < x and x < response);
}
brettviren commented 3 months ago

Hi @wenqiang-gu.

I would like Drifter not to hold this new feature.

Really, depo filtering of any kind does not belong inside Drifter. In hindsight I see now that the existing Xregion selection is better placed in some upstream IDepoSetFilter. Maybe one day we clean this up. But, for now it can stay.

So, for this new filtering, I suggest a new IDepoSetFilter be developed to implement this "bent cathode" partitioning of the depo set. The implementation (BentCathodeFilter or some better name) can be generic and live inside WCT and follow your f(x,y) design.

A change in configuration for the Drifter's will be needed. Each Drifter must have its "cathode" plane placed beyond the furthest part of the real, bent cathode. This means the Xregion of a pair of Drifter's operating on either side of one physical cathode will overlap. But, that is okay as they rely on the two BentCathodeFilter's operating on ether side to have consistent configurations to properly partition the depo set. That is, we rely on the BentCathodeFilter to assure there are no depos between the bent cathode and the Xregion cathode for each Drifter. The Xregion cathode can even be defined at "infinity".

wenqiang-gu commented 3 months ago

@brettviren Please let me know if I understand correctly. Below is an example of a sim+nf+sp configuration. A DepoSetDrifter is applied before the "fanpipe" of all anodes.

local simnfsp_pipes = [
  g.pipeline([
               splusn_pipes[n], // from sim.jsonnet
               nf_pipes[n],
               sp_pipes[n],
             ],
             ‘simnfsp_pipe_%d' % n)
  for n in std.range(0, std.length(tools.anodes) - 1);

local fanpipe = f.fanpipe('DepoSetFanout', simnfsp_pipes, 'FrameFanin’);

local graph = g.pipeline([deposet, deposet_drifter, fanpipe, sink]);

Now if I understand your proposal, another IDepoSetFilter (e.g.,"bentCahodeFilter") can be applied to update DepoSet for each anode. So we can move both "bentCahodeFilter" and DepoSetDrifter into the fanpipe, for example:

local simnfsp_pipes = [
  g.pipeline([
           bent_cathode_deposet_filiters[n], // IDepoSetFIlter
           deposet_drifters[n], // IDepoSetFilter
               splusn_pipes[n],
               nf_pipes[n],
               sp_pipes[n],
             ],
             ‘simnfsp_pipe_%d' % n)
  for n in std.range(0, std.length(tools.anodes) - 1);

local fanpipe = f.fanpipe('DepoSetFanout', simnfsp_pipes, 'FrameFanin’);

local graph = g.pipeline([deposet, fanpipe, sink]);
brettviren commented 3 months ago

Hi @wenqiang-gu

Really good question as it changes my mind drastically.

I think we can greatly simplify things from what you describe.

The main thing is we should keep our data flow graphs simple with a single DepoSetDrifter for the whole detector. We should avoid having to make, eg, one per each "pipeline".

I see now the key to doing this is to remove Drifter::Xregion from Drifter (making Drifter simpler) and extend the new Xregion to handle both the current planar representation and the new "bent" representation.

This actually brings me back closer to your original idea. Effectively we can make Drifter support "bent" cathodes but through an improved and removed Xregion. Sorry for not seeing this better path until now!

There are a few details to describe. Too much to type into GitHub. It may be easier for me to implement them than describe them in text. Actually, I volunteer to do the work. In any case we could chat on zoom, maybe after today's WC meeting?

For configuration, what I envision is that the new Xregion will check if the value of any of its cathode/response/anode "planes" is a scalar value (as currently we give) or if it is an object. This would allow user to give configuration for the new Xregion like:

{
  cathode: {
    x: [...]
    y: [...]
    z: [...]
  },
  response: 10*wc.cm,
  anode: 0.0,
}

In this case, Xregion sees cathode as an object instead of a scalar value and loads the the three cathode.{x,y,z} arrays into a interpolation (probably bicubic).

The inside_*(double x) methods need to change to inside_*(const Point& pt). This will actually make Drifter even slightly simpler.

With this change, Drifter configuration can be backwards compatible while transparently supporting arbitrary shaped (well, must be single-valued) "planes" defined by interpolated samples of their surface.

wenqiang-gu commented 3 months ago

Hi @brettviren , that sounds good to me -- please proceed with the implementation.

If you need a configuration for testing, here is a standalone WCT simulation in PDSP: https://github.com/WireCell/wire-cell-toolkit/blob/master/cfg/pgrapher/experiment/pdsp/wct-sim-deposet-drifter.jsonnet However, please remove the "deposet_rotate" from the pgraph, it's just for test purpose:

local graph = g.pipeline([depos, plainbagger,/* deposet_rotate,*/ setdrifter, parallel_graph, sink]);

Also, do you think it would be beneficial to ask Mike Mooney to present at next week's WireCell meeting on this bent cathode issue? This could help ensure we're aligned on their request.

brettviren commented 2 months ago

Fixed in #336