Bugfix/blips on empty events #901
Conversation
|
https://docs.google.com/presentation/d/1igXPcS_Qst99layiHuMEZiUxz4EHMxGHutf7Mr-jbT0/edit?slide=id.g3bbdbff95f2_0_1#slide=id.g3bbdbff95f2_0_1 here is a useful overview set of slides for these changes. Sorry for the many commits most are just adding/removing debug. For the hardcoded gaushit parameter the change in this PR is sufficient for any standard running someone does/for Fall production. Longterm I think these should not be hardcoded and just use the fcl-provided hit collection, and we assume the user knows what they are doing. The existing MC truth back tracking works. I just added a few additional breaks in MC truth blips. I found a way to get labeling I am totally happy with but the performance hit is too large for anything other than a small production aimed specifically at blips. I am happy with the current state of the MCTruth matching in this PR now. I promise some future slides looking at how the MC matching is doing in Fall production pre/post this PR in the next ~week |
henrylay97
left a comment
henrylay97
left a comment
There was a problem hiding this comment.
Hi Jacob, thanks for the PR. I have a few comments.
My main concern is around the fix for the empty event problem. Rather than finding the root cause of the problem, the solution seems to be to provide a second vector that you hope won't be empty as well.
I would rather the actual cause of the segfault be found and prevented. This would also remove any need for hard coded module labels.
The rest looks good to me, happy to discuss more offline!
| const double dir((tpcgeom.DriftSign() == geo::DriftSign::Negative) ? +1.0 : -1.0); | ||
| float x_ticks_coefficient = kDriftVelocity*kTickPeriod; | ||
| float goofy_offset = -xyz.X() / (dir * x_ticks_coefficient); |
There was a problem hiding this comment.
There is lots of bad mis-alignment in these files, should be fixable with a quick de-tabify offline :)
There was a problem hiding this comment.
Went through and generally fixed tab alignment. There are a few lines of white space with awkward tabbing still, but all the content looks aligned properly on my end.
| //kXTicksOffsets[cstat][tpc][pl] += fTimeOffset[pl]; | ||
| std::cout << " offsetting plane " << pl << " by " << fTimeOffset[pl] << " ticks " << std::endl; | ||
|
|
||
| //kXTicksOffsets[cstat][tpc][pl] += fTimeOffset[pl]; |
There was a problem hiding this comment.
Why is this commented and not just removed?
There was a problem hiding this comment.
Good point. The timeOffset is applied independently later to the hits directly. This can be removed.
| std::vector<geo::WireID> wids = wireReadoutGeom->Get().ChannelToWire( (*(simchanlist[chIndex])).Channel() ); //Not sure why this is a vector, but it should have len 1 | ||
| const geo::PlaneID& planeID = wids[0].planeID(); | ||
| if(int(planeID.Plane) != fCaloPlane) continue; //only take calorimetry plane IDE values | ||
| unsigned int MaxTDCTick = 4294967294; // 1 under 32 bit unsigned int max |
There was a problem hiding this comment.
Why 1 less? And you could use UINT_MAX here for better readability and to keep it consistent for any system.
There was a problem hiding this comment.
Didn't know that existed. Updated it now
| if(evt.getByLabel("specialblipgaushit",hitHandleGH)) art::fill_ptr_vector(hitlistGH, hitHandleGH); | ||
| else if(evt.getByLabel("gaushit",hitHandleGH)) art::fill_ptr_vector(hitlistGH, hitHandleGH); |
There was a problem hiding this comment.
Module labels should never be hard coded. I appreciate this was true before this PR but given you're touching the code at this point it needs fixing.
I also don't understand the use of two different hit collections - only one of them will have been an input to the BlipReco and so only one of them is a valid input for the truth matching. If the code crashes on an empty vector somewhere then the fix should be a catch to stop that happening rather than the providing of a different vector.
There was a problem hiding this comment.
I have removed the hitHandleGH, and will remove the other spots that gaushit gets a special default privileging. It will just be the user's responsibility to make sure they have ran all the necessary inputs
| art::FindManyP<recob::Track> fmtrk(hitHandle,evt,fTrkProducer); | ||
| art::FindManyP<recob::Track> fmtrkGH(hitHandleGH,evt,fTrkProducer); | ||
| art::FindMany<simb::MCParticle,anab::BackTrackerHitMatchingData> fmhh(hitHandleGH,evt,"gaushitTruthMatch"); | ||
| art::FindMany<simb::MCParticle, anab::BackTrackerHitMatchingData> fmhh(hitHandleGH,evt,"blipgaushitTruthMatch"); |
There was a problem hiding this comment.
Added a new fcl param that defaults to blipgaushitTruthMatch
| std::map< int, int > map_gh; | ||
| // if input collection is already gaushit, this is trivial | ||
| if( fHitProducer == "gaushit" ) { | ||
| if( fHitProducer == "gaushit" || fHitProducer == "specialblipgaushit") { |
There was a problem hiding this comment.
This was some functionality to allow the user to provide a hitcollection that never had truth matching or PandoraTracking ran on it. That was the whole purpose of the hitHandleGH, and I think that was more useful for uB reprocessing and just got included when the blips were ported over to SBND. I've left in map_gh, but this branching is eliminated and it just assumes the trivial case where our hitcollection was ran through truth matching and tracking.
|
|
||
| // find associated track | ||
| if( fHitProducer == "gaushit" && fmtrk.isValid() ) { | ||
| if( fHitProducer == "specialblipgaushit" && fmtrk.isValid() ) { |
There was a problem hiding this comment.
Again same comment, this feels like a solution that creates a lot of future problems. Why should this block only be run if the HitProducer is a certain value?
There was a problem hiding this comment.
Changed to just check if fmtrk is a valid handle.
| printf("G4 IDs contib \n"); | ||
| for(int g4ID : hc.G4IDs) printf(" %i ", g4ID); | ||
| printf("\n"); |
There was a problem hiding this comment.
Is BlipAna run as standard in production?
For example, in the calib ntuples? If so, the print statements need to be hidden behind a debug flag.
There was a problem hiding this comment.
BlipAna is run in standard production for calib ntuples.
The function these print statements live in is only called with a debug flag
if( fDebugMode ) PrintClusterInfo(clust);
… some branches to assume the hitHandle the user provided was run through truthmatching and pandoraTrack
|
Thanks for the suggestions, most are in. Let me check compilation and a simple output in a typical Fall production run |
|
Merged in develop and the compiler works and entries in the reco2 files look good. After getting the normal successful processing table I got this error |
A few bugs existed in old versions of blip reco
The first bug was related to a hardcoded hitHandleGH label, which looked at gaushit rather than specialBlipGausHit. On ~emtpy events the lower threshold (specialBlipGausHit) will have a few entries, but the code uses hitHandleGH for truth matching, and that higher threshold hit collection can be empty. That leads to crashes.
For now the hardcoded label is just updated to accommodate specialBlipGausHit.
The second bug is not critical, as access hitclusters this way is done through a blip, so the ID is the same as the blip you are looking at and the clusters included must be matched (into this blip). The hitcluster collection values are updated after the struct assignment (copy constructor), so I updated the blip alg to also update the blip.Cluster[i].blipID, and isMatched.
The third bug I am still digging into. I was worried that the switch from sim::energy_deposition to sim::channel->IDE was causing issues, but in the new arrangement (with an optional interface to use either) gives consistent results between them.
I have some notes I can build into slides for the first two but would like to spend a little more time on the third bug.