-
Notifications
You must be signed in to change notification settings - Fork 54
Feature/adding blip to caf #871
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
base: develop
Are you sure you want to change the base?
Conversation
|
Adding a comment to each of these to track all 4 related PR. The changes to sbnobj, and sbnanaobj are fully independent of any other changes, so they can be approved first. sbndcode changes rely on sbnobj, so it will have to wait for the first approval. A later simple PR will delete the (now duplicated) class files in the BlipUtils folder here sbncode changes rely on both sbnobj and sbnanaobj, so that will have to wait for both of the first two approvals. |
|
Hi @Jjm321814 can you remove the CRT changes from this PR and instead merge in |
|
Merged in develop. Testing Compile now |
|
Compiled fine |
|
Hi Jacob, I will review properly after tomorrow's reconstruction meeting. In the meanwhile there are still changes to CRT and LightPropagation files that shouldn't be necessary. Can you remove these? Thanks! |
PetrilloAtWork
left a comment
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Left three comments.
| sbnobj::Common_CRT | ||
| #sbndcode_CosmicIdUtils | ||
| sbndcode_BlipUtils | ||
| sbndobj_BlipDataTypes |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I requested the change of this library name.
sbndcode/BlipRecoSBND/CMakeLists.txt
Outdated
| sbndcode_RecoUtils | ||
| sbndcode_OpDetSim | ||
| sbndcode_BlipUtils | ||
| sbndobj_BlipDataTypes |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Please fix the indentation: copy the line above verbatim, then adjust the library name.
| #sbndcode_CosmicIdUtils | ||
| ) | ||
|
|
||
| art_dictionary(DICTIONARY_LIBRARIES sbndcode_BlipUtils) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
No dictionary is produced any more: just remove this art_dictionary line.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I still want to import BlipUtils.h to the Alg folder though, so don't I still need this dictionary to access functions defined there?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Yes, the dictionary does not help with that.
You will need to add sbndcode_BlipUtils in the LIBRARIES list of CMakeLists.txt in the algorithm folder, but the dictionary won't be involved.
…ixed a mislabel in the fcl
…s. This is useful if we were to rerurn older sbndcode reco1 output over a new reco2 version. Thats a really limited usecase but it could come up
henrylay97
left a comment
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
A few small comments
| //#include "larcorealg/Geometry/GeometryCore.h" | ||
| //#include "larcorealg/Geometry/WireReadoutGeom.h" |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
In general can we remove commented headers :)
| sbndcode_CRTUtils | ||
| sbnobj::Common_CRT | ||
| #sbndcode_CosmicIdUtils | ||
| sbnobj::SBND_Blip |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
indentation
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I'm not sure how to fix this -- It looks normal in all my emacs, vim, and vscode. Can you please check out the branch to check if its just a github display issue?
| const int kMaxEDeps = 10000; | ||
| const int kMaxTrkPts = 2000; | ||
| const int kMaxTrkPts = 2000; | ||
| const int kNplanes = blip::kNplanes; |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Can we align equals and values inkeeping with previous lines.
| physics.producers.pandoraShowerRazzle.SimChannelLabel: "simtpc2d:simpleSC" | ||
| physics.producers.pandoraTrackDazzle.SimChannelLabel: "simtpc2d:simpleSC" | ||
| physics.producers.cnnid.WireLabel: "simtpc2d:gauss" | ||
| #physics.producers.cnnid.WireLabel: "simtpc2d:gauss" |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Reiterating this comment - I see there is now a line below that points to "simtpc2d:dnnsp" so L129 can just be removed.
| sbnd_reco1_producers.specialblipgaushit.HitFinderToolVec.CandidateHitsPlane0.MinHitHeight: 8 | ||
| sbnd_reco1_producers.specialblipgaushit.HitFinderToolVec.CandidateHitsPlane1.MinHitHeight: 15 | ||
| sbnd_reco1_producers.specialblipgaushit.HitFinderToolVec.CandidateHitsPlane2.MinHitHeight: 5 |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
There will be conflicts with #887 just be careful during merging that things don't get reverted.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I'll merge that branch in here to make sure they match
PetrilloAtWork
left a comment
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
First, thank you for taking the effort of reducing the TVector3 footprint in our environment.
One problem you are trying to solve here is that you were strongly recommended to drop TVector3 for GenVector but then you have to interface with existing code that still uses the former, and you don't want to undertake the effort of refactoring the whole library.
Conversion routines from and to geo::Point_t and geo::Vector_t are available in the geo::vect namespace. Many of them are straightforward, with some exceptions (notably, the "averaging" of points). I will leave comments around to facilitate the adoption of these tools. If you find that the specification of geo::vect is too frequent in a function, you may also consider to add in the function using namespace geo::vect;, but use it sparsely since it has the potential to change the meaning of existing code.
Another remark is that the "sigma" of the blips seems to me to have an odd definition. This is preexisting the pull request.
| pinfo.position.SetXYZ(0.5*(pinfo.startPoint.X()+pinfo.endPoint.X()), | ||
| 0.5*(pinfo.startPoint.Y()+pinfo.endPoint.Y()), | ||
| 0.5*(pinfo.startPoint.Z()+pinfo.endPoint.Z()) ); | ||
|
|
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
This is an unweighted middle point, so:
| pinfo.position.SetXYZ(0.5*(pinfo.startPoint.X()+pinfo.endPoint.X()), | |
| 0.5*(pinfo.startPoint.Y()+pinfo.endPoint.Y()), | |
| 0.5*(pinfo.startPoint.Z()+pinfo.endPoint.Z()) ); | |
| pinfo.position = geo::vect::middlePoint({ pinfo.startPoint, pinfo.endPoint }); |
(#include "larcorealg/Geometry/geo_vector_utils.h")
Note the braces around the list of points (you can have any number of them).
This function is a simpler interface to MiddlePointAccumulator, which I'll mention below.
| tblip.Position.SetXYZ( w1*tblip.Position.X() + w2*pinfo.position.X(), | ||
| w1*tblip.Position.Y() + w2*pinfo.position.Y(), | ||
| w1*tblip.Position.Z() + w2*pinfo.position.Z()); |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
This is the hardest one, being weighted. I would use the accumulator directly:
| tblip.Position.SetXYZ( w1*tblip.Position.X() + w2*pinfo.position.X(), | |
| w1*tblip.Position.Y() + w2*pinfo.position.Y(), | |
| w1*tblip.Position.Z() + w2*pinfo.position.Z()); | |
| geo::vect::MiddlePointAccumulator mpalg; | |
| mpalg.add(tblip.Position, w1); | |
| mpalg.add(pinfo.Position, w2); | |
| tblip.Position = mpalg.middlePoint(); |
(#include "larcorealg/Geometry/geo_vector_utils.h")
| // together a blip that happened much later but in the same spot) | ||
| if( fabs(blip_i.Time - blip_j.Time) > 5 ) continue; | ||
| float d = (blip_i.Position-blip_j.Position).Mag(); | ||
| float d = TMath::Sqrt((blip_i.Position-blip_j.Position).Mag2()); |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
You can directly use .R():
| float d = TMath::Sqrt((blip_i.Position-blip_j.Position).Mag2()); | |
| float d = (blip_i.Position-blip_j.Position).R(); |
Also note that for determining the closest point, you don't really need to take the square root of the square distance, because the minimum distance is also the minimum distance square, so you can look for the latter and save a tiny bit of computing time.
An awkward feature of the new vectors is that they have Mag2() and R(), but not Mag() and R2(). 🤷
| blip_i.Position.SetXYZ( w1*blip_i.Position.X() + w2*blip_j.Position.X(), | ||
| w1*blip_i.Position.Y() + w2*blip_j.Position.Y(), | ||
| w1*blip_i.Position.Z() + w2*blip_j.Position.Z()); |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
As above:
| blip_i.Position.SetXYZ( w1*blip_i.Position.X() + w2*blip_j.Position.X(), | |
| w1*blip_i.Position.Y() + w2*blip_j.Position.Y(), | |
| w1*blip_i.Position.Z() + w2*blip_j.Position.Z()); | |
| geo::vect::MiddlePointAccumulator mpalg; | |
| mpalg.add(blip_i.Position, w1); | |
| mpalg.add(blip_j.Position, w2); | |
| blip_i.Position = mpalg.middlePoint(); |
However, here we are in a loop context, so I recommend that the accumulation is done first, and then the results saved. You can declare your accumulators (MiddlePointAccumulator, StatCollector) in the outer loop, dealing with blip_i (and feeding in there the blip i), and then fill them in the inner loop:
- initialization in outer loop
for(size_t i=0; i<vtb.size(); i++){ if( isGrouped.at(i) ) continue; else isGrouped.at(i) = true; auto& blip_i = vtb.at(i); geo::vect::MiddlePointAccumulator position; lar::util::StatCollector<double> time; float const weight_i = blip_i.DepElectrons; position.add(blip_i.Position, weight_i); time.add(blip_i.time, weight_i); // ... for(size_t j=i+1; j<vtb.size(); j++){
- filling in here:
Suggested change
blip_i.Position.SetXYZ( w1*blip_i.Position.X() + w2*blip_j.Position.X(), w1*blip_i.Position.Y() + w2*blip_j.Position.Y(), w1*blip_i.Position.Z() + w2*blip_j.Position.Z()); float const weight_j = blip_j.DepElectrons; position.add(blip_j.Position.X(), weight_j); time.add(blip_j.time, weight_j); // ... - putting back into
blip_iat the end of the loop:}//loop over blip_j blip_i.Position = position.middlePoint(); blip_i.time = time.Average(); // ... blip_i.ID = vtb_merged.size();
| // TO-DO: if this track starts or ends at a TPC boundary, | ||
| // we should extend p1 or p2 to outside the AV to avoid blind spots | ||
| TVector3 bp = newBlip.Position; | ||
| TVector3 bp(newBlip.Position.X(), newBlip.Position.Y(), newBlip.Position.Z()); |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Conversion geo::Point_t → TVector3:
| TVector3 bp(newBlip.Position.X(), newBlip.Position.Y(), newBlip.Position.Z()); | |
| TVector3 bp = geo::vect::toTVector3(newBlip.Position); |
(#include "larcorealg/Geometry/geo_vectors_utils_TVector.h"); or using the more generic convertTo():
| TVector3 bp(newBlip.Position.X(), newBlip.Position.Y(), newBlip.Position.Z()); | |
| TVector3 bp = geo::vect::convertTo<TVector3>(newBlip.Position); |
(for this one: #include "larcorealg/Geometry/geo_vectors_utils.h").
| for(auto& v : wirex ) newblip.Position.SetXYZ( newblip.Position.X() + v.X() * fact, | ||
| newblip.Position.Y() + v.Y() * fact, | ||
| newblip.Position.Z() + v.Z() * fact); |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
This should work:
| for(auto& v : wirex ) newblip.Position.SetXYZ( newblip.Position.X() + v.X() * fact, | |
| newblip.Position.Y() + v.Y() * fact, | |
| newblip.Position.Z() + v.Z() * fact); | |
| geo::vect::MiddlePointAccumulator position; | |
| position.add(begin(wirex), end(wirex)); | |
| newblip.Position = position.middlePoint(); |
I think this will work, still, given the following, it would be better having wirex redefined as std::vector<geo::Point_t>, especially considering that it is used only for this computation, and much better to have two accumulators initialised before the loop (a StatCollector2D sc for y and z — x is 0 —, and optionally a MiddlePointAccumulator for the position — but that one can be extracted by the other two), collect the statistics adding the coordinates in the loop, and finally computing the σ (sqrt(sc.VarianceX() + sc.VarianceY())) and average position ({ 0, sc.AverageX(), sc.AverageY() }). Also, the definition of σ feels off to me: it's the average of distances, while the usual RMS is the square root of the average of the squared distances (as in my example with VarianceX()/Y()).
| #sbndcode_CosmicIdUtils | ||
| ) | ||
|
|
||
| art_dictionary(DICTIONARY_LIBRARIES sbndcode_BlipUtils) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Yes, the dictionary does not help with that.
You will need to add sbndcode_BlipUtils in the LIBRARIES list of CMakeLists.txt in the algorithm folder, but the dictionary won't be involved.


Moved the blip-related structs and classes to sbnobj.
This PR must not be approved until SBNSoftware/sbnobj#155 is approved/released! Otherwise it will break blip production.
I also had to make a few CRT changes to successfully compile off the current sbnobj file.
https://sbn-docdb.fnal.gov/cgi-bin/sso/ShowDocument?docid=44445