boundaryhits sometimes counts incident electrons more than once

cameronc137 commented 5 years ago

Environment: (where does this bug occur, have you tried other environments)

branch (master for latest released): develop
revision (HEAD for most recent): HEAD
OS or system: Centos 6.9
Special ROOT or Geant4 versions? Newest of both

Steps to reproduce: (give a step by step account of how to trigger the bug)

Execute the macros/runexample_optical.mac (preferably with a thousand or so events)
Run the pe analysis code (from analysis/bin/pe after running make install in the analysis/build directory) and move the output file to the side.
Edit the pe.cc code and remove the && Qcounted==0 and related checks in the quartz, reflector air, and light guide air logical checks (around line 155, 160, and 173 and 182).
Rerun the pe analysis code and now compare the electron hits in the quartz/PE spectra.
(You can also see this issue by scanning the raw data and looking for electron hits in det==50001, but this is more tedious or requires putting cout statements in the PE code).

Expected Result: (what do you expect when you execute the steps above)

Whether or not I check if an electron has been counted yet in the detector should not matter, since the detector type for quartz and reflector or light guide air are set to boundaryhits. Boundaryhits should not register a hit more than once for a given electron passing through a medium.

Actual Result: (what do you get when you execute the steps above)

Instead, maybe 10% of the time, an electron will in fact step multiple (~5 or 6) times inside the quartz volume, erroneously increasing the number of counts in my PE analysis script by a factor of 2 or so above what it should be.

Probably these are just the deltas, and because the deltas a born inside of the quartz bar they are not subject to a boundary crossing before they start hitting in the material, hence this behavior.

Is there a way to amend the boundaryhits detector type to accomodate for particles born inside of the volume as well? Something like:

if mtrid==0 then behave like a boundary hit, and if mtrid>0 then only register the first hit.

cameronc137 commented 5 years ago

Is this still the case after the updates in #182? I'll check soon

wdconinc commented 5 years ago

I don't think this should have changed. I looked at this issue yesterday but wish there was a better demonstration than running the full pe code... Something like: if you run T->Scan("hit.trid","hit.pid == 0 && hit.det == soandso") then I see more than one entry per hit.trid for event soandso.

cameronc137 commented 5 years ago

I had to update the beam generator to allow macro commands with 3 vectors to work: 259fb69e6e9edd7ec73945695554f52d5d8f97e5

Running macros/runexample_optical.mac appears to still have this issue, though you have to run ~ 1000 events instead of 100 to see anything meaningful (which could be a shower instead of a multi-stepping delta, I'm not sure)

I run T->Scan("hit.z","hit.pid==11 && hit.mtrid>0") and I see a few entries where there are multiple electron hits inside the material and clearly not at the front face (which is around -126mm in this geometry).

cameronc137 commented 5 years ago

Ahhh, updated again to keep the units with 3 vectors in remollGenBeam (so basically I reverted back to the correct macro command functionality and updated the direction setting format in the macros instead): 59ac2f7ab1329941c706875488a4fd657a3f5dcd

wdconinc commented 5 years ago

I see, this is actually an issue with electrons, not opticalphotons (why runexample_optical.mac?). And it only affects the following detectors:

11:56:10 wdconinc@herakles ~/git/remoll (develop *$%=) $ grep -A1 boundaryhits geometry/detector_5open.gdml
         <auxiliary auxtype="DetType" auxvalue="boundaryhits"/> 
         <auxiliary auxtype="DetNo" auxvalue="50001"/>  
--
         <auxiliary auxtype="DetType" auxvalue="boundaryhits"/> 
         <auxiliary auxtype="DetNo" auxvalue="50101"/>  
--
         <auxiliary auxtype="DetType" auxvalue="boundaryhits"/> 
         <auxiliary auxtype="DetNo" auxvalue="50401"/>

So, this here is the output of T->Draw("hit.det:hit.z","hit.pid==11 && hit.mtrid>0"):

However, scanning through them does not show them to be multiple entries per track.

re-root [5] T->Scan("hit.trid:hit.det:hit.z","hit.pid==11 && hit.mtrid>0") 
***********************************************************
*    Row   * Instance *  hit.trid *   hit.det *     hit.z *
***********************************************************
*      247 *        1 *      1258 *     50001 * -127.5232 *
*      247 *        2 *       467 *     50101 * -119.7268 *
*      247 *        3 *       467 *     50101 * -116.9691 *
*      247 *        4 *       467 *     50101 * -128.8717 *
*      247 *      292 *       467 *     50201 * -116.4340 *
*      247 *      293 *       467 *     50201 * -128.5010 *
*      247 *      295 *       467 *     50201 * -129.2281 *
*      247 *      659 *       467 *     50301 * -117.1379 *
*      247 *      660 *       467 *     50301 * -116.3997 *
*      266 *        5 *       819 *     50101 * -118.9939 *
*      266 *      282 *       819 *     50201 * -116.9327 *
*      266 *      709 *       819 *     50301 * -117.5542 *
*      308 *        1 *         3 *     50001 * -127.3463 *
*      311 *        1 *      1757 *     50001 * -113.6626 *
*      453 *        0 *         3 *     50001 * -125.8716 *
*      576 *        5 *        97 *     50101 * -118.5526 *
*      576 *      254 *        97 *     50201 * -117.0388 *
*      576 *     1024 *        97 *     50301 * -117.6461 *
*      777 *        4 *         6 *     50101 * -127.7598 *
*      777 *       41 *         6 *     50201 * -128.1028 *
*      793 *        0 *         2 *     50001 * -124.6902 *
*      902 *        1 *         2 *     50201 * -129.9814 *
*      911 *        4 *       124 *     50101 * -118.5587 *
*      911 *      124 *       124 *     50201 * -117.0277 *
*      911 *      527 *       124 *     50301 * -117.6391 *
Type <CR> to continue or q to quit ==>     
*     1164 *        2 *        11 *     50001 * -127.4970 *
*     1202 *        0 *         2 *     50001 * -125.8198 *
*     1213 *        4 *       115 *     50101 * -118.7665 *
*     1213 *      480 *       115 *     50201 * -117.0045 *
*     1213 *     1168 *       115 *     50301 * -117.6136 *
*     1217 *        1 *      1124 *     50001 * -116.0762 *
*     1261 *        2 *       897 *     50101 * -119.0713 *
*     1261 *      281 *       897 *     50201 * -117.0351 *
*     1261 *      564 *       897 *     50301 * -117.6227 *
*     1281 *        1 *      1159 *     50001 * -116.0681 *
*     1376 *        1 *         4 *     50101 * -117.9944 *
*     1376 *      243 *         4 *     50201 * -117.0824 *
*     1376 *      741 *         4 *     50301 * -117.7119 *
*     1385 *        1 *         7 *     50001 * -126.6894 *
*     1388 *        0 *         2 *     50001 * -123.1895 *
*     1407 *        1 *      1114 *     50001 * -115.7254 *
*     1472 *        0 *         2 *     50001 * -127.1664 *
*     1475 *        0 *         2 *     50001 * -122.7099 *
*     1576 *        1 *      1032 *     50001 * -116.2187 *
*     1588 *        0 *         2 *     50001 * -125.3732 *
*     1588 *        1 *         2 *     50001 * -123.6794 *
*     1767 *        0 *         2 *     50001 * -124.7878 *
*     1832 *        1 *      1139 *     50101 * -118.8377 *
*     1832 *      332 *      1139 *     50201 * -115.6763 *
*     1832 *      333 *      1139 *     50201 * -115.7491 *
Type <CR> to continue or q to quit ==> 
*     1832 *      711 *      1139 *     50301 * -116.3458 *
*     1856 *        0 *         2 *     50001 * -122.6925 *
*     1897 *      137 *       525 *     50201 * -120.6091 *
*     1962 *        1 *         6 *     50001 * -127.7700 *
*     1962 *        3 *         6 *     50101 * -128.2725 *
*     1962 *        6 *         6 *     50101 * -118.6320 *
*     1962 *      307 *         6 *     50201 * -117.0131 *
*     1962 *     1043 *         6 *     50301 * -117.6261 *
*     2000 *        0 *         2 *     50001 * -124.2395 *
*     2131 *        1 *         2 *     50201 * -142.1546 *
*     2156 *        1 *         4 *     50001 * -127.8697 *
*     2156 *        3 *         4 *     50101 * -129.2273 *
*     2156 *        4 *         4 *     50101 * -127.8697 *
*     2156 *        5 *       929 *     50101 * -119.6607 *
*     2156 *        6 *       929 *     50101 * -116.9233 *
*     2156 *       35 *         4 *     50201 * -129.2786 *
*     2156 *      557 *       929 *     50301 * -117.0489 *
*     2156 *      595 *       929 *     50401 * -128.9825 *
*     2156 *     1198 *       929 *     50701 * -169.9439 *
*     2164 *        0 *         6 *     50001 * -124.8726 *
*     2270 *        1 *      1277 *     50001 * -113.5629 *
*     2483 *        0 *         2 *     50001 * -127.2839 *
*     2522 *      297 *      1030 *     50201 * -111.2253 *
*     2540 *        1 *      1178 *     50001 * -115.9638 *
*     2607 *      296 *      1065 *     50201 * -116.9274 *
Type <CR> to continue or q to quit ==> 
*     2795 *      262 *      1157 *     50201 * -115.5775 *
*     2855 *        0 *         2 *     50001 * -124.9402 *
*     3076 *        0 *         2 *     50001 * -124.2084 *
*     3122 *      241 *       956 *     50201 * -115.7711 *
*     3258 *        0 *         2 *     50001 * -126.3509 *
*     3286 *        0 *         2 *     50001 * -125.9705 *
*     3381 *       30 *         8 *     50201 * -185.2009 *
*     3441 *        1 *         7 *     50001 * -125.3510 *
*     3658 *        1 *       153 *     50001 * -125.9399 *
*     3753 *        1 *       685 *     50101 * -119.2936 *
*     3753 *        2 *       685 *     50101 * -115.9125 *
*     3753 *      271 *       685 *     50201 * -116.4720 *
*     3753 *      611 *       685 *     50301 * -115.9924 *
*     3765 *        0 *         2 *     50001 * -126.2919 *
*     3765 *        1 *         2 *     50001 * -125.1620 *
*     3885 *        0 *         2 *     50001 * -127.0162 *
*     3928 *        0 *         2 *     50001 * -122.8099 *
*     4306 *        1 *      1113 *     50101 * -118.0689 *
*     4306 *      300 *      1113 *     50201 * -116.8948 *
*     4306 *      605 *      1113 *     50301 * -117.5665 *
*     4371 *        0 *         2 *     50001 * -125.1471 *
*     4452 *        0 *         2 *     50001 * -127.0611 *
*     4511 *      219 *       997 *     50201 * -115.4482 *
*     4630 *      294 *      1143 *     50201 * -116.8573 *
*     4804 *        2 *       116 *     50101 * -119.7248 *
Type <CR> to continue or q to quit ==> 
*     4804 *      269 *       116 *     50201 * -116.8338 *
*     4804 *      753 *       116 *     50301 * -117.4509 *
*     4869 *        0 *         2 *     50001 * -127.0377 *
*     4936 *        0 *         2 *     50001 * -127.4512 *
*     4987 *      281 *      1184 *     50201 * -117.0949 *
*     5026 *        1 *         3 *     50001 * -122.4368 *
*     5040 *        1 *      1636 *     50001 * -117.3977 *
*     5040 *      313 *       897 *     50201 * -116.5811 *
*     5040 *      353 *      1636 *     50201 * -115.0293 *
*     5265 *        1 *       460 *     50001 * -112.4467 *
*     5704 *        0 *         2 *     50001 * -124.6688 *
*     5723 *        0 *         2 *     50001 * -124.2151 *
*     5881 *        0 *         2 *     50001 * -126.2764 *
*     5900 *        0 *         2 *     50001 * -115.3196 *
*     5977 *        1 *         5 *     50001 * -127.7638 *
*     5977 *        4 *         5 *     50101 * -128.1762 *
*     6042 *        1 *      1363 *     50001 * -121.0184 *
*     6112 *        1 *         4 *     50001 * -127.7445 *
*     6112 *        3 *         4 *     50101 * -127.8416 *
*     6425 *        1 *      1129 *     50001 * -114.3462 *
*     6453 *        6 *         5 *     50101 * -128.6900 *
*     6453 *       27 *         5 *     50201 * -148.9641 *
*     6476 *        2 *      1204 *     50101 * -117.9020 *
*     6476 *      275 *      1204 *     50201 * -117.1041 *
*     6476 *      573 *      1204 *     50301 * -117.7290 *
Type <CR> to continue or q to quit ==> 
*     6495 *        1 *      1274 *     50001 * -114.0991 *
*     6760 *        9 *         2 *     50201 * -131.8565 *
*     6800 *        0 *         2 *     50001 * -126.4940 *
*     6926 *        0 *         2 *     50001 * -122.3507 *
*     6935 *        0 *        10 *     50001 * -125.5220 *
*     6952 *        0 *         2 *     50001 * -127.3510 *
*     7069 *        0 *         2 *     50001 * -126.8147 *
*     7129 *        0 *         2 *     50001 * -125.0415 *
*     7150 *        0 *         2 *     50001 * -124.0440 *
*     7319 *        0 *         2 *     50001 * -124.0280 *
*     7414 *        1 *         6 *     50001 * -127.8223 *
*     7414 *        3 *         6 *     50101 * -128.5112 *
*     7424 *        0 *         2 *     50001 * -124.4654 *
*     7424 *        1 *         2 *     50001 * -124.2455 *
*     7424 *        2 *        10 *     50001 * -124.2434 *
*     7513 *        0 *         2 *     50001 * -125.5268 *
*     7522 *        1 *         2 *     50001 * -125.0093 *
*     7629 *        1 *      1764 *     50001 * -124.9847 *
*     7733 *        0 *         2 *     50001 * -127.2963 *
*     7769 *        1 *      1258 *     50001 * -115.1034 *
*     7787 *        5 *        98 *     50101 * -121.4561 *
*     7787 *      255 *        98 *     50201 * -116.2336 *
*     7787 *      622 *        98 *     50301 * -116.8933 *
*     7799 *        0 *         2 *     50001 * -122.6590 *
*     7807 *        1 *      1236 *     50001 * -114.8734 *
Type <CR> to continue or q to quit ==> 
*     7894 *        4 *         2 *     50201 * -135.7543 *
*     8096 *        1 *       817 *     50001 * -113.1372 *
*     8228 *        0 *         2 *     50001 * -126.9913 *
*     8228 *        1 *         2 *     50001 * -127.2463 *
*     8228 *        2 *         2 *     50001 * -127.2623 *
*     8228 *        3 *         6 *     50001 * -127.2705 *
*     8607 *        1 *      1264 *     50001 * -114.4635 *
*     8607 *        2 *      1294 *     50001 * -114.4319 *
*     8688 *        0 *         2 *     50001 * -123.8131 *
*     8733 *      313 *      1404 *     50201 * -117.1180 *
*     8794 *        0 *         3 *     50001 * -124.5799 *
*     9008 *        0 *         2 *     50001 * -123.8210 *
*     9041 *        0 *         5 *     50001 * -122.7711 *
*     9205 *        0 *         2 *     50001 * -125.6304 *
*     9223 *        0 *         2 *     50001 * -125.0461 *
*     9383 *        1 *         2 *     50201 * -134.9497 *
*     9413 *        1 *      1170 *     50001 * -116.4832 *
*     9550 *      247 *      1038 *     50201 * -114.4763 *
*     9627 *        0 *         2 *     50001 * -125.5247 *
*     9661 *        0 *         2 *     50001 * -125.5932 *
*     9663 *        1 *       802 *     50001 * -116.7205 *
*     9663 *        2 *       802 *     50001 * -116.7293 *
*     9673 *        0 *         2 *     50001 * -124.9405 *
*     9688 *        1 *      2276 *     50001 * -115.9183 *
*     9703 *        0 *         2 *     50001 * -123.6140 *
Type <CR> to continue or q to quit ==> 
*     9800 *      285 *       760 *     50201 * -117.0775 *
*     9800 *      640 *       760 *     50301 * -117.6309 *
*     9816 *        1 *      1166 *     50001 * -114.5949 *
***********************************************************
==> 178 selected entries
(long long) 178

Am I correct in assuming that in event 1588 you don't want the same trid to have two hits?

wdconinc commented 5 years ago

As an aside. I think all this messing with which hits are stored is just a recipe for confusion. If I were to restart (but that would break all existing analysis), we'd store all hits by default. Period. All these exceptions, in both directions, are causing confusion.

wdconinc commented 5 years ago

That event 1588 might actually legitimately cross the boundary twice since the hit.pz changes sign...

JeffersonLab / remoll