From 401f71220a3b0e3956bd8c9bbaa6e206ceb3ef7f Mon Sep 17 00:00:00 2001 From: albertoesmp Date: Mon, 25 Sep 2023 13:45:39 +0200 Subject: [PATCH] Solved ray-box intersection issue for early aborts --- scripts/debug/hda_pulse_records_plotter.py | 217 +++++++++++++++++- src/alg/raycast/GroveKDTreeRaycaster.h | 8 + src/alg/raycast/KDGroveRaycaster.cpp | 10 +- src/alg/raycast/KDGroveRaycaster.h | 4 + src/alg/raycast/KDTreeRaycaster.cpp | 157 ++++++++++++- src/alg/raycast/KDTreeRaycaster.h | 14 ++ src/alg/raycast/Raycaster.h | 4 + src/dataanalytics/HDA_GlobalVars.cpp | 117 ++++++++++ src/dataanalytics/HDA_GlobalVars.h | 38 ++- src/dataanalytics/HDA_GlobalVarsReporter.h | 40 ++++ src/dataanalytics/HDA_PulseRecorder.cpp | 2 +- src/dataanalytics/HDA_PulseRecorder.h | 45 +++- src/scanner/MultiScanner.cpp | 3 +- src/scanner/MultiScanner.h | 3 +- src/scanner/Scanner.h | 3 +- src/scanner/ScanningDevice.cpp | 27 ++- src/scanner/ScanningDevice.h | 3 +- src/scanner/SingleScanner.cpp | 3 +- src/scanner/SingleScanner.h | 3 +- .../detector/DynFullWaveformPulseRunnable.cpp | 7 + .../detector/DynFullWaveformPulseRunnable.h | 3 + .../detector/FullWaveformPulseRunnable.cpp | 59 ++++- .../detector/FullWaveformPulseRunnable.h | 6 +- src/scene/Scene.cpp | 43 ++++ src/scene/Scene.h | 4 + src/scene/primitives/AABB.cpp | 5 +- src/scene/primitives/AABB.h | 21 ++ 27 files changed, 811 insertions(+), 38 deletions(-) diff --git a/scripts/debug/hda_pulse_records_plotter.py b/scripts/debug/hda_pulse_records_plotter.py index 524871fd4..ecff9d91c 100755 --- a/scripts/debug/hda_pulse_records_plotter.py +++ b/scripts/debug/hda_pulse_records_plotter.py @@ -94,7 +94,20 @@ def read_records(path, sep=','): 'radius_step': subray_sim[:, 1], 'circle_steps': subray_sim[:, 2], 'circle_step': subray_sim[:, 3], - 'divergence_angle_rad': subray_sim[:, 4] + 'divergence_angle_rad': subray_sim[:, 4], + 'ray_dir_norm': subray_sim[:, 5], + 'subray_dir_norm': subray_sim[:, 6], + 'ray_subray_angle_rad': subray_sim[:, 7], + 'ray_subray_sign_check': subray_sim[:, 8], + 'subray_tmin': subray_sim[:, 9], + 'subray_tmax': subray_sim[:, 10], + 'subray_rt_pos_dir': subray_sim[:, 11], # Ray-tracing positive dir. + 'subray_rt_neg_dir': subray_sim[:, 12], # Ray-tracing negative dir. + 'subray_rt_par_dir': subray_sim[:, 13], # Ray-tracing parallel dir. + 'subray_rt_no_second': subray_sim[:, 14], # Ray-tracing no second half + 'subray_rt_no_first': subray_sim[:, 15], # Ray-tracing no first half + 'subray_rt_both': subray_sim[:, 16], # Ray-tracing both sides + 'subray_rt_both2': subray_sim[:, 17] # Ray-tracing both, second try } @@ -120,6 +133,8 @@ def plot_records(arec, brec, outdir): do_incidence_angle_plots(arec, brec, outdir) do_by_incidence_angle_plots(arec, brec, outdir) do_subray_hit_plots(arec, brec, outdir) + do_ray_subray_plots(arec, brec, outdir) + do_raytracing_plots(arec, brec, outdir) def validate_record(key, rec, recid): @@ -614,7 +629,205 @@ def do_subray_hit_plots(arec, brec, outdir): # TODO Rethink : Implement fig.tight_layout() # Save figure to file and remove it from memory - fig.savefig(os.path.join(outdir, 'subray_hit_plots.png')) + fig.savefig(os.path.join(outdir, 'subray_hit.png')) + fig.clear() + plt.close(fig) + + +def do_ray_subray_plots(arec, brec, outdir): + # Validate ray subray data + if( + not validate_record('ray_dir_norm', arec, 'a') or + not validate_record('subray_dir_norm', arec, 'a') or + not validate_record('ray_subray_angle_rad', arec, 'a') or + not validate_record('ray_subray_sign_check', arec, 'a') or + not validate_record('subray_tmin', arec, 'a') or + not validate_record('subray_tmax', arec, 'a') or + not validate_record('ray_dir_norm', brec, 'b') or + not validate_record('subray_dir_norm', brec, 'b') or + not validate_record('ray_subray_angle_rad', brec, 'b') or + not validate_record('ray_subray_sign_check', brec, 'b') or + not validate_record('subray_tmin', brec, 'b') or + not validate_record('subray_tmax', brec, 'b') + ): + print('Cannot do ray-subray plots') + return + # Do the ray subray plots + fig = init_figure() # Initialize figure + # CASE A + ax = fig.add_subplot(3, 4, 1) # Initialize ray norm subplot + do_incidence_angle_subplot( + fig, ax, arec['ray_dir_norm'], + xlabel='Ray direction norm' + ) + ax = fig.add_subplot(3, 4, 2) # Initialize subray norm subplot + do_incidence_angle_subplot( + fig, ax, arec['subray_dir_norm'], + xlabel='Subay direction norm' + ) + ax = fig.add_subplot(3, 4, 5) # Initialize ray-subray angle subplot + do_incidence_angle_subplot( + fig, ax, arec['ray_subray_angle_rad']*180/np.pi, + xlabel='Ray-subray angle (deg)' + ) + ax = fig.add_subplot(3, 4, 6) # Initialize ray-subray sign check subplot + do_incidence_angle_subplot( + fig, ax, arec['ray_subray_sign_check'], + xlabel='Sign equality check' + ) + ax = fig.add_subplot(3, 4, 9) # Initialize subray tmin subplot + do_incidence_angle_subplot( + fig, ax, arec['subray_tmin'], + xlabel='Subray $t_{\\mathrm{min}}$' + ) + ax = fig.add_subplot(3, 4, 10) # Initialize subray tmin subplot + do_incidence_angle_subplot( + fig, ax, arec['subray_tmax'], + xlabel='Subray $t_{\\mathrm{max}}$' + ) + # CASE B + ax = fig.add_subplot(3, 4, 3) # Initialize ray norm subplot + do_incidence_angle_subplot( + fig, ax, brec['ray_dir_norm'], + xlabel='Ray direction norm' + ) + ax = fig.add_subplot(3, 4, 4) # Initialize subray norm subplot + do_incidence_angle_subplot( + fig, ax, brec['subray_dir_norm'], + xlabel='Subay direction norm' + ) + ax = fig.add_subplot(3, 4, 7) # Initialize ray-subray angle subplot + do_incidence_angle_subplot( + fig, ax, brec['ray_subray_angle_rad']*180/np.pi, + xlabel='Ray-subray angle (deg)' + ) + ax = fig.add_subplot(3, 4, 8) # Initialize ray-subray sign check subplot + do_incidence_angle_subplot( + fig, ax, brec['ray_subray_sign_check'], + xlabel='Sign equality check' + ) + ax = fig.add_subplot(3, 4, 11) # Initialize subray tmin subplot + do_incidence_angle_subplot( + fig, ax, brec['subray_tmin'], + xlabel='Subray $t_{\\mathrm{min}}$' + ) + ax = fig.add_subplot(3, 4, 12) # Initialize subray tmin subplot + do_incidence_angle_subplot( + fig, ax, brec['subray_tmax'], + xlabel='Subray $t_{\\mathrm{max}}$' + ) + fig.tight_layout() + # Save figure to file and remove it from memory + fig.savefig( + os.path.join(outdir, 'ray_subray.png') + ) + fig.clear() + plt.close(fig) + + +def do_raytracing_plots(arec, brec, outdir): + # TODO Rethink : Implement + if( + not validate_record('subray_hit', arec, 'a') or + not validate_record('subray_rt_pos_dir', arec, 'a') or + not validate_record('subray_rt_neg_dir', arec, 'a') or + not validate_record('subray_rt_par_dir', arec, 'a') or + not validate_record('subray_rt_no_second', arec, 'a') or + not validate_record('subray_rt_no_first', arec, 'a') or + not validate_record('subray_rt_both', arec, 'a') or + not validate_record('subray_rt_both2', arec, 'a') or + not validate_record('subray_hit', brec, 'b') or + not validate_record('subray_rt_pos_dir', brec, 'b') or + not validate_record('subray_rt_neg_dir', brec, 'b') or + not validate_record('subray_rt_par_dir', brec, 'b') or + not validate_record('subray_rt_no_second', brec, 'b') or + not validate_record('subray_rt_no_first', brec, 'b') or + not validate_record('subray_rt_both', brec, 'b') or + not validate_record('subray_rt_both2', brec, 'b') + ): + print('Cannot do ray-tracing plots') + return + # Do the ray-tracing plots + fig = init_figure() # Initialize figure + # CASE A + ax = fig.add_subplot(4, 7, 1) # Initialize positive direction hist + do_subray_hit_subplot_hist( + fig, ax, arec['subray_hit'], arec['subray_rt_pos_dir'], + ylabel='Absolute' + ) + ax = fig.add_subplot(4, 7, 2) # Initialize negative direction hist + do_subray_hit_subplot_hist( + fig, ax, arec['subray_hit'], arec['subray_rt_neg_dir'], + ) + ax = fig.add_subplot(4, 7, 3) # Initialize parallel direction hist + do_subray_hit_subplot_hist( + fig, ax, arec['subray_hit'], arec['subray_rt_par_dir'], + ) + ax = fig.add_subplot(4, 7, 4) # Initialize no second half hist + do_subray_hit_subplot_hist( + fig, ax, arec['subray_hit'], arec['subray_rt_no_second'], + ) + ax = fig.add_subplot(4, 7, 5) # Initialize no second half hist + do_subray_hit_subplot_hist( + fig, ax, arec['subray_hit'], arec['subray_rt_no_first'], + ) + ax = fig.add_subplot(4, 7, 6) # Initialize both sides hist + do_subray_hit_subplot_hist( + fig, ax, arec['subray_hit'], arec['subray_rt_both'], + ) + ax = fig.add_subplot(4, 7, 7) # Initialize both sides hist + do_subray_hit_subplot_hist( + fig, ax, arec['subray_hit'], arec['subray_rt_both2'], + ) + ax = fig.add_subplot(4, 7, 8) # Initialize positive direction hist + do_subray_hit_subplot_hist( + fig, ax, arec['subray_hit'], arec['subray_rt_pos_dir'], + relative=True, + xlabel='Positive direction', + ylabel='Relative (100%)' + ) + ax = fig.add_subplot(4, 7, 9) # Initialize negative direction hist + do_subray_hit_subplot_hist( + fig, ax, arec['subray_hit'], arec['subray_rt_neg_dir'], + relative=True, + xlabel='Negative direction', + ) + ax = fig.add_subplot(4, 7, 10) # Initialize parallel direction hist + do_subray_hit_subplot_hist( + fig, ax, arec['subray_hit'], arec['subray_rt_par_dir'], + relative=True, + xlabel='Parallel direction', + ) + ax = fig.add_subplot(4, 7, 11) # Initialize no second half hist + do_subray_hit_subplot_hist( + fig, ax, arec['subray_hit'], arec['subray_rt_no_second'], + relative=True, + xlabel='No second half', + ) + ax = fig.add_subplot(4, 7, 12) # Initialize no second half hist + do_subray_hit_subplot_hist( + fig, ax, arec['subray_hit'], arec['subray_rt_no_first'], + relative=True, + xlabel='No first half', + ) + ax = fig.add_subplot(4, 7, 13) # Initialize both sides hist + do_subray_hit_subplot_hist( + fig, ax, arec['subray_hit'], arec['subray_rt_both'], + relative=True, + xlabel='Both sides', + ) + ax = fig.add_subplot(4, 7, 14) # Initialize both sides hist + do_subray_hit_subplot_hist( + fig, ax, arec['subray_hit'], arec['subray_rt_both2'], + relative=True, + xlabel='Both sides ($2^{\mathrm{nd}}$ try)', + ) + # CASE B + fig.tight_layout() + # Save figure to file and remove it from memory + fig.savefig( + os.path.join(outdir, 'subray_ray_tracing.png') + ) fig.clear() plt.close(fig) diff --git a/src/alg/raycast/GroveKDTreeRaycaster.h b/src/alg/raycast/GroveKDTreeRaycaster.h index 85ceac12e..031fb7bb1 100644 --- a/src/alg/raycast/GroveKDTreeRaycaster.h +++ b/src/alg/raycast/GroveKDTreeRaycaster.h @@ -80,8 +80,16 @@ class GroveKDTreeRaycaster : double tmin, double tmax, bool groundOnly +#ifdef DATA_ANALYTICS + ,std::vector &subraySimRecord, + bool const isSubray=false +#endif ) override {return KDTreeRaycaster::search( rayOrigin, rayDir, tmin, tmax, groundOnly +#ifdef DATA_ANALYTICS + ,subraySimRecord, + isSubray +#endif );} // *** GROVE DYNAMIC TREE METHODS *** // diff --git a/src/alg/raycast/KDGroveRaycaster.cpp b/src/alg/raycast/KDGroveRaycaster.cpp index 002fcd6fb..007e57159 100644 --- a/src/alg/raycast/KDGroveRaycaster.cpp +++ b/src/alg/raycast/KDGroveRaycaster.cpp @@ -26,13 +26,21 @@ RaySceneIntersection * KDGroveRaycaster::search( double tmin, double tmax, bool groundOnly +#ifdef DATA_ANALYTICS + ,std::vector &subraySimRecord, + bool const isSubray +#endif ){ std::map out; size_t const m = grove->getNumTrees(); - RaySceneIntersection *bestRSI = nullptr, *rsi = nullptr; + RaySceneIntersection *bestRSI = nullptr, *rsi; for(size_t i = 0 ; i < m ; ++i){ rsi = grove->getTreeShared(i)->search( rayOrigin, rayDir, tmin, tmax, groundOnly +#ifdef DATA_ANALYTICS + ,subraySimRecord, + isSubray +#endif ); if( bestRSI==nullptr || diff --git a/src/alg/raycast/KDGroveRaycaster.h b/src/alg/raycast/KDGroveRaycaster.h index 943dcfbef..8dc4a18b0 100644 --- a/src/alg/raycast/KDGroveRaycaster.h +++ b/src/alg/raycast/KDGroveRaycaster.h @@ -56,6 +56,10 @@ class KDGroveRaycaster : public Raycaster{ double tmin, double tmax, bool groundOnly +#ifdef DATA_ANALYTICS + ,std::vector &subraySimRecord, + bool const isSubray=false +#endif ) override; // *** KDGROVE RELATED METHODS *** // diff --git a/src/alg/raycast/KDTreeRaycaster.cpp b/src/alg/raycast/KDTreeRaycaster.cpp index 70bdc048a..ccd860ce7 100644 --- a/src/alg/raycast/KDTreeRaycaster.cpp +++ b/src/alg/raycast/KDTreeRaycaster.cpp @@ -1,6 +1,11 @@ #include "KDTreeRaycaster.h" #include "logging.hpp" +#ifdef DATA_ANALYTICS +#include +using helios::analytics::HDA_GV; +#endif + using namespace std; map KDTreeRaycaster::searchAll( @@ -32,6 +37,10 @@ RaySceneIntersection* KDTreeRaycaster::search( double const tmin, double const tmax, bool const groundOnly +#ifdef DATA_ANALYTICS + ,std::vector &subraySimRecord, + bool const isSubray +#endif ){ // Prepare search KDTreeRaycasterSearch search(rayDir, rayOrigin, groundOnly); @@ -43,12 +52,38 @@ RaySceneIntersection* KDTreeRaycaster::search( search.rayOriginArray.push_back(rayOrigin.z); // Do recursive search +#ifdef DATA_ANALYTICS + std::size_t positiveDirectionCount = 0, negativeDirectionCount = 0; + std::size_t noSecondHalfCount = 0, noFirstHalfCount = 0; + std::size_t parallelDirectionCount = 0, bothSidesCount = 0; + std::size_t bothSidesSecondTryCount = 0; +#endif Primitive* prim = this->search_recursive( this->root.get(), tmin-epsilon, tmax+epsilon, search +#ifdef DATA_ANALYTICS + ,positiveDirectionCount, + negativeDirectionCount, + parallelDirectionCount, + noSecondHalfCount, + noFirstHalfCount, + bothSidesCount, + bothSidesSecondTryCount, + isSubray +#endif ); +#ifdef DATA_ANALYTICS + // Write counts to subray simulation record + subraySimRecord[11] = (double) positiveDirectionCount; + subraySimRecord[12] = (double) negativeDirectionCount; + subraySimRecord[13] = (double) parallelDirectionCount; + subraySimRecord[14] = (double) noSecondHalfCount; + subraySimRecord[15] = (double) noFirstHalfCount; + subraySimRecord[16] = (double) bothSidesCount; + subraySimRecord[17] = (double) bothSidesSecondTryCount; +#endif // Handle search output if (prim == nullptr) return nullptr; @@ -170,13 +205,22 @@ Primitive* KDTreeRaycaster::search_recursive( double const tmin, double const tmax, KDTreeRaycasterSearch &search +#ifdef DATA_ANALYTICS + ,size_t &positiveDirectionCount, + size_t &negativeDirectionCount, + size_t ¶llelDirectionCount, + size_t &noSecondHalfCount, + size_t &noFirstHalfCount, + size_t &bothSidesCount, + size_t &bothSidesSecondTryCount, + bool const isSubray +#endif ) const { if(node==nullptr) return nullptr; // Null nodes cannot contain primitives Primitive* hitPrim = nullptr; // ######### BEGIN If node is a leaf, perform ray-primitive intersection on all primitives in the leaf's bucket ########### if (node->splitAxis == -1) { - //logging::DEBUG("leaf node:"); for (auto prim : *node->primitives) { double const newDistance = prim->getRayIntersectionDistance( search.rayOrigin, search.rayDir @@ -198,10 +242,43 @@ Primitive* KDTreeRaycaster::search_recursive( // (in other leaves, with other primitives) that are *closer* to // the ray originWaypoint. If this was the case, the returned intersection // would be wrong. - - if( - (newDistance > 0 && newDistance < search.closestHitDistance) && - (newDistance >= tmin && newDistance <= tmax) +#ifdef DATA_ANALYTICS + bool const nonNegativeDistance = newDistance > 0; + bool const closerThanClosest = + newDistance < search.closestHitDistance; + bool const tminCheck = newDistance >= tmin; + bool const tmaxCheck = newDistance <= tmax; + if(!nonNegativeDistance){ // For any traced ray + HDA_GV.incrementRaycasterLeafNegativeDistancesCount(); + } else if(!closerThanClosest){ + HDA_GV.incrementRaycasterLeafFurtherThanClosestCount(); + } else if(!tminCheck){ + HDA_GV.incrementRaycasterLeafFailedTminCheckCount(); + } + else if(!tmaxCheck){ + HDA_GV.incrementRaycasterLeafFailedTmaxCheckCount(); + } + if(isSubray){ // For any traced ray that is a subray + if(!nonNegativeDistance){ + HDA_GV.incrementSubrayLeafNegativeDistancesCount(); + } else if(!closerThanClosest){ + HDA_GV.incrementSubrayLeafFurtherThanClosestCount(); + } else if(!tminCheck){ + HDA_GV.incrementSubrayLeafFailedTminCheckCount(); + } + else if(!tmaxCheck){ + HDA_GV.incrementSubrayLeafFailedTmaxCheckCount(); + } + + } + if( + nonNegativeDistance && closerThanClosest && + tminCheck && tmaxCheck +#else + if( + newDistance > 0 && newDistance < search.closestHitDistance && + newDistance >= tmin && newDistance <= tmax +#endif ){ if( !search.groundOnly || @@ -229,6 +306,9 @@ Primitive* KDTreeRaycaster::search_recursive( // Case 1: Ray goes in positive direction - it passes through the left side first, then through the right: if (search.rayDirArray[a] > 0) { +#ifdef DATA_ANALYTICS + ++positiveDirectionCount; +#endif first = node->left; second = node->right; @@ -237,6 +317,9 @@ Primitive* KDTreeRaycaster::search_recursive( } // Case 2: Ray goes in negative direction - it passes through the right side first, then through the left: else if (search.rayDirArray[a] < 0) { +#ifdef DATA_ANALYTICS + ++negativeDirectionCount; +#endif first = node->right; second = node->left; @@ -245,6 +328,9 @@ Primitive* KDTreeRaycaster::search_recursive( } // Case 3: Ray goes parallel to the split plane - it passes through only one side, depending on it's originWaypoint: else { +#ifdef DATA_ANALYTICS + ++parallelDirectionCount; +#endif first = (search.rayOriginArray[a] < node->splitPos) ? node->left : node->right; second = (search.rayOriginArray[a] < node->splitPos) ? @@ -257,22 +343,75 @@ Primitive* KDTreeRaycaster::search_recursive( // thit >= tmax means that the ray crosses the split plane *after it has already left the volume*. // In this case, it never enters the second half. if (thit >= tmax) { - hitPrim = search_recursive(first, tmin, tmax, search); +#ifdef DATA_ANALYTICS + ++noSecondHalfCount; +#endif + hitPrim = search_recursive( + first, tmin, tmax, search +#ifdef DATA_ANALYTICS + ,positiveDirectionCount, + negativeDirectionCount, + parallelDirectionCount, + noSecondHalfCount, + noFirstHalfCount, + bothSidesCount, + bothSidesSecondTryCount +#endif + ); } // thit <= tmin means that the ray crosses the split plane *before it enters the volume*. // In this case, it never enters the first half: else if (thit <= tmin) { - hitPrim = search_recursive(second, tmin, tmax, search); +#ifdef DATA_ANALYTICS + ++noFirstHalfCount; +#endif + hitPrim = search_recursive( + second, tmin, tmax, search +#ifdef DATA_ANALYTICS + ,positiveDirectionCount, + negativeDirectionCount, + parallelDirectionCount, + noSecondHalfCount, + noFirstHalfCount, + bothSidesCount, + bothSidesSecondTryCount +#endif + ); } // Otherwise, the ray crosses the split plane within the volume. // This means that it passes through both sides: else { - hitPrim = search_recursive(first, tmin, thit+epsilon, search); +#ifdef DATA_ANALYTICS + ++bothSidesCount; +#endif + hitPrim = search_recursive( + first, tmin, thit+epsilon, search +#ifdef DATA_ANALYTICS + ,positiveDirectionCount, + negativeDirectionCount, + parallelDirectionCount, + noSecondHalfCount, + noFirstHalfCount, + bothSidesCount, + bothSidesSecondTryCount +#endif + ); if (hitPrim == nullptr) { - hitPrim = search_recursive(second, thit-epsilon, tmax, search); + hitPrim = search_recursive( + second, thit-epsilon, tmax, search +#ifdef DATA_ANALYTICS + ,positiveDirectionCount, + negativeDirectionCount, + parallelDirectionCount, + noSecondHalfCount, + noFirstHalfCount, + bothSidesCount, + bothSidesSecondTryCount +#endif + ); } } diff --git a/src/alg/raycast/KDTreeRaycaster.h b/src/alg/raycast/KDTreeRaycaster.h index 108f33dbd..c42b27b21 100644 --- a/src/alg/raycast/KDTreeRaycaster.h +++ b/src/alg/raycast/KDTreeRaycaster.h @@ -140,6 +140,10 @@ class KDTreeRaycasterSearch{ double const tmin, double const tmax, bool const groundOnly +#ifdef DATA_ANALYTICS + ,std::vector &subraySimRecord, + bool const isSubray=false +#endif ) override; protected: @@ -183,5 +187,15 @@ class KDTreeRaycasterSearch{ double const tmin, double const tmax, KDTreeRaycasterSearch &search +#ifdef DATA_ANALYTICS + ,size_t &positiveDirectionCount, + size_t &negativeDirectionCount, + size_t ¶llelDirectionCount, + size_t &noSecondHalfCount, + size_t &noFirstHalfCount, + size_t &bothSidesCount, + size_t &bothSidesSecondTryCount, + bool const isSubray=false +#endif ) const; }; \ No newline at end of file diff --git a/src/alg/raycast/Raycaster.h b/src/alg/raycast/Raycaster.h index 9d6d03462..875363b9f 100644 --- a/src/alg/raycast/Raycaster.h +++ b/src/alg/raycast/Raycaster.h @@ -62,5 +62,9 @@ class Raycaster{ double tmin, double tmax, bool groundOnly +#ifdef DATA_ANALYTICS + ,std::vector &subraySimRecord, + bool const isSubray=false +#endif ) = 0; }; \ No newline at end of file diff --git a/src/dataanalytics/HDA_GlobalVars.cpp b/src/dataanalytics/HDA_GlobalVars.cpp index f3964f8e6..d618cc407 100644 --- a/src/dataanalytics/HDA_GlobalVars.cpp +++ b/src/dataanalytics/HDA_GlobalVars.cpp @@ -53,6 +53,16 @@ HDA_GlobalVars & HDA_GlobalVars::incrementNonIntersectiveSubraysCount(){ return *this; } +HDA_GlobalVars & +HDA_GlobalVars::incrementNonIntersectiveSubraysDueToNullTimeCount(){ + std::unique_lock lock( + nonIntersectiveSubraysDueToNullTimeCount_mutex + ); + ++nonIntersectiveSubraysDueToNullTimeCount; + return *this; + +} + HDA_GlobalVars & HDA_GlobalVars::incrementSubrayIntersectionCount() { std::unique_lock lock(subrayIntersectionCount_mutex); ++subrayIntersectionCount; @@ -65,6 +75,62 @@ HDA_GlobalVars & HDA_GlobalVars::incrementSubrayNonIntersectionCount() { return *this; } +HDA_GlobalVars & +HDA_GlobalVars::incrementRaycasterLeafNegativeDistancesCount(){ + std::unique_lock lock( + raycasterLeafNegativeDistancesCount_mutex + ); + ++raycasterLeafNegativeDistancesCount; + return *this; +} + +HDA_GlobalVars & +HDA_GlobalVars::incrementRaycasterLeafFurtherThanClosestCount(){ + std::unique_lock lock( + raycasterLeafFurtherThanClosestCount_mutex + ); + ++raycasterLeafFurtherThanClosestCount; + return *this; +} + +HDA_GlobalVars & HDA_GlobalVars::incrementRaycasterLeafFailedTminCheckCount(){ + std::unique_lock lock(raycasterLeafFailedTminCheckCount_mutex); + ++raycasterLeafFailedTminCheckCount; + return *this; +} + +HDA_GlobalVars & HDA_GlobalVars::incrementRaycasterLeafFailedTmaxCheckCount(){ + std::unique_lock lock(raycasterLeafFailedTmaxCheckCount_mutex); + ++raycasterLeafFailedTmaxCheckCount; + return *this; +} + +HDA_GlobalVars & +HDA_GlobalVars::incrementSubrayLeafNegativeDistancesCount(){ + std::unique_lock lock(subrayLeafNegativeDistancesCount_mutex); + ++subrayLeafNegativeDistancesCount; + return *this; +} + +HDA_GlobalVars & +HDA_GlobalVars::incrementSubrayLeafFurtherThanClosestCount(){ + std::unique_lock lock(subrayLeafFurtherThanClosestCount_mutex); + ++subrayLeafFurtherThanClosestCount; + return *this; +} + +HDA_GlobalVars & HDA_GlobalVars::incrementSubrayLeafFailedTminCheckCount(){ + std::unique_lock lock(subrayLeafFailedTminCheckCount_mutex); + ++subrayLeafFailedTminCheckCount; + return *this; +} + +HDA_GlobalVars & HDA_GlobalVars::incrementSubrayLeafFailedTmaxCheckCount(){ + std::unique_lock lock(subrayLeafFailedTmaxCheckCount_mutex); + ++subrayLeafFailedTmaxCheckCount; + return *this; +} + // *** READ METHODS *** // // ********************** // std::size_t HDA_GlobalVars::getGeneratedSubraysCount(){ @@ -96,6 +162,13 @@ std::size_t HDA_GlobalVars::getNonIntersectiveSubraysCount() { return nonIntersectiveSubraysCount; } +std::size_t HDA_GlobalVars::getNonIntersectiveSubraysDueToNullTimeCount(){ + std::unique_lock lock( + nonIntersectiveSubraysDueToNullTimeCount_mutex + ); + return nonIntersectiveSubraysDueToNullTimeCount; +} + std::size_t HDA_GlobalVars::getSubrayIntersectionCount() { std::unique_lock lock(subrayIntersectionCount_mutex); return subrayIntersectionCount; @@ -111,6 +184,50 @@ std::size_t HDA_GlobalVars::getIntensityComputationsCount(){ return intensityComputationsCount; } +std::size_t HDA_GlobalVars::getRaycasterLeafNegativeDistancesCount(){ + std::unique_lock lock( + raycasterLeafNegativeDistancesCount_mutex + ); + return raycasterLeafNegativeDistancesCount; +} + +std::size_t HDA_GlobalVars::getRaycasterLeafFurtherThanClosestCount(){ + std::unique_lock lock( + raycasterLeafFurtherThanClosestCount_mutex + ); + return raycasterLeafFurtherThanClosestCount; +} + +std::size_t HDA_GlobalVars::getRaycasterLeafFailedTminCheckCount(){ + std::unique_lock lock(raycasterLeafFailedTminCheckCount_mutex); + return raycasterLeafFailedTminCheckCount; +} + +std::size_t HDA_GlobalVars::getRaycasterLeafFailedTmaxCheckCount(){ + std::unique_lock lock(raycasterLeafFailedTmaxCheckCount_mutex); + return raycasterLeafFailedTmaxCheckCount; +} + +unsigned long HDA_GlobalVars::getSubrayLeafNegativeDistancesCount(){ + std::unique_lock lock(subrayLeafNegativeDistancesCount_mutex); + return subrayLeafNegativeDistancesCount; +} + +unsigned long HDA_GlobalVars::getSubrayLeafFurtherThanClosestCount(){ + std::unique_lock lock(subrayLeafFurtherThanClosestCount_mutex); + return subrayLeafFurtherThanClosestCount; +} + +unsigned long HDA_GlobalVars::getSubrayLeafFailedTminCheckCount(){ + std::unique_lock lock(subrayLeafFailedTminCheckCount_mutex); + return subrayLeafFailedTminCheckCount; +} + +unsigned long HDA_GlobalVars::getSubrayLeafFailedTmaxCheckCount(){ + std::unique_lock lock(subrayLeafFailedTmaxCheckCount_mutex); + return subrayLeafFailedTmaxCheckCount; +} + }} diff --git a/src/dataanalytics/HDA_GlobalVars.h b/src/dataanalytics/HDA_GlobalVars.h index 024f97e91..b8d7d20ff 100644 --- a/src/dataanalytics/HDA_GlobalVars.h +++ b/src/dataanalytics/HDA_GlobalVars.h @@ -19,9 +19,18 @@ class HDA_GlobalVars{ std::size_t generatedSubraysCount; std::size_t intersectiveSubraysCount; std::size_t nonIntersectiveSubraysCount; + std::size_t nonIntersectiveSubraysDueToNullTimeCount; std::size_t subrayIntersectionCount; std::size_t subrayNonIntersectionCount; std::size_t intensityComputationsCount; + std::size_t raycasterLeafNegativeDistancesCount; + std::size_t raycasterLeafFurtherThanClosestCount; + std::size_t raycasterLeafFailedTminCheckCount; + std::size_t raycasterLeafFailedTmaxCheckCount; + unsigned long subrayLeafNegativeDistancesCount; + unsigned long subrayLeafFurtherThanClosestCount; + unsigned long subrayLeafFailedTminCheckCount; + unsigned long subrayLeafFailedTmaxCheckCount; protected: // *** CONCURRENCY HANDLING ATTRIBUTES *** // @@ -31,9 +40,18 @@ class HDA_GlobalVars{ std::mutex generatedSubraysCount_mutex; std::mutex intersectiveSubraysCount_mutex; std::mutex nonIntersectiveSubraysCount_mutex; + std::mutex nonIntersectiveSubraysDueToNullTimeCount_mutex; std::mutex subrayIntersectionCount_mutex; std::mutex subrayNonIntersectionCount_mutex; std::mutex intensityComputationsCount_mutex; + std::mutex raycasterLeafNegativeDistancesCount_mutex; + std::mutex raycasterLeafFurtherThanClosestCount_mutex; + std::mutex raycasterLeafFailedTminCheckCount_mutex; + std::mutex raycasterLeafFailedTmaxCheckCount_mutex; + std::mutex subrayLeafNegativeDistancesCount_mutex; + std::mutex subrayLeafFurtherThanClosestCount_mutex; + std::mutex subrayLeafFailedTminCheckCount_mutex; + std::mutex subrayLeafFailedTmaxCheckCount_mutex; public: // *** CONSTRUCTION / DESTRUCTION *** // @@ -44,6 +62,7 @@ class HDA_GlobalVars{ generatedSubraysCount(0), intersectiveSubraysCount(0), nonIntersectiveSubraysCount(0), + nonIntersectiveSubraysDueToNullTimeCount(0), subrayIntersectionCount(0), subrayNonIntersectionCount(0), intensityComputationsCount(0) @@ -57,9 +76,18 @@ class HDA_GlobalVars{ HDA_GlobalVars & incrementGeneratedSubraysCount(); HDA_GlobalVars & incrementIntersectiveSubraysCount(); HDA_GlobalVars & incrementNonIntersectiveSubraysCount(); + HDA_GlobalVars & incrementNonIntersectiveSubraysDueToNullTimeCount(); HDA_GlobalVars & incrementSubrayIntersectionCount(); HDA_GlobalVars & incrementSubrayNonIntersectionCount(); HDA_GlobalVars & incrementIntensityComputationsCount(); + HDA_GlobalVars & incrementRaycasterLeafNegativeDistancesCount(); + HDA_GlobalVars & incrementRaycasterLeafFurtherThanClosestCount(); + HDA_GlobalVars & incrementRaycasterLeafFailedTminCheckCount(); + HDA_GlobalVars & incrementRaycasterLeafFailedTmaxCheckCount(); + HDA_GlobalVars & incrementSubrayLeafNegativeDistancesCount(); + HDA_GlobalVars & incrementSubrayLeafFurtherThanClosestCount(); + HDA_GlobalVars & incrementSubrayLeafFailedTminCheckCount(); + HDA_GlobalVars & incrementSubrayLeafFailedTmaxCheckCount(); // *** READ METHODS *** // // ********************** // @@ -68,10 +96,18 @@ class HDA_GlobalVars{ std::size_t getGeneratedRaysAfterEarlyAbortCount(); std::size_t getIntersectiveSubraysCount(); std::size_t getNonIntersectiveSubraysCount(); + std::size_t getNonIntersectiveSubraysDueToNullTimeCount(); std::size_t getSubrayIntersectionCount(); std::size_t getSubrayNonIntersectionCount(); std::size_t getIntensityComputationsCount(); - + std::size_t getRaycasterLeafNegativeDistancesCount(); + std::size_t getRaycasterLeafFurtherThanClosestCount(); + std::size_t getRaycasterLeafFailedTminCheckCount(); + std::size_t getRaycasterLeafFailedTmaxCheckCount(); + unsigned long getSubrayLeafNegativeDistancesCount(); + unsigned long getSubrayLeafFurtherThanClosestCount(); + unsigned long getSubrayLeafFailedTminCheckCount(); + unsigned long getSubrayLeafFailedTmaxCheckCount(); }; }} diff --git a/src/dataanalytics/HDA_GlobalVarsReporter.h b/src/dataanalytics/HDA_GlobalVarsReporter.h index 4f03dc431..19ace808f 100644 --- a/src/dataanalytics/HDA_GlobalVarsReporter.h +++ b/src/dataanalytics/HDA_GlobalVarsReporter.h @@ -27,6 +27,16 @@ class HDA_GlobalVarsReporter{ // Initialize string stream to build the print std::stringstream ss; + // Extract variables for the sake of convenience + size_t const raycasterLeafFailedTminCheckCount = + HDA_GV.getRaycasterLeafFailedTminCheckCount(); + size_t const raycasterLeafFailedTmaxCheckCount = + HDA_GV.getRaycasterLeafFailedTmaxCheckCount(); + size_t const subrayLeafFailedTminCheckCount = + HDA_GV.getSubrayLeafFailedTminCheckCount(); + size_t const subrayLeafFailedTmaxCheckCount = + HDA_GV.getSubrayLeafFailedTmaxCheckCount(); + // Print global variables ss << "HDA GLOBAL VARS REPORT:\n\n"; ss << "Generated rays before early abort: " @@ -38,12 +48,42 @@ class HDA_GlobalVarsReporter{ << gv.getIntersectiveSubraysCount() << "\n"; ss << "Non-intersective subrays: " << gv.getNonIntersectiveSubraysCount() << "\n"; + ss << "\tNon-intersective subrays due to null time: " + << gv.getNonIntersectiveSubraysDueToNullTimeCount() << "\n"; ss << "Subray intersections: " << gv.getSubrayIntersectionCount() << "\n"; ss << "Subray non-intersections: " << gv.getSubrayNonIntersectionCount() << "\n"; ss << "Number of computed intensities: " << gv.getIntensityComputationsCount() << "\n"; + ss << "Raycaster fails-on-leaves distribution:\n" + << "\tNegative distances: " + << HDA_GV.getRaycasterLeafNegativeDistancesCount() << "\n" + << "\tFurther than current closest: " + << HDA_GV.getRaycasterLeafFurtherThanClosestCount() << "\n" + << "\tFailed tmin checks: " + << raycasterLeafFailedTminCheckCount << "\n" + << "\tFailed tmax checks: " + << raycasterLeafFailedTmaxCheckCount << "\n" + << "\tFailed t checks: " + << ( + raycasterLeafFailedTminCheckCount + + raycasterLeafFailedTmaxCheckCount + ) << "\n" + << "\tSubrays only:\n" + << "\t\tNegative distances: " + << HDA_GV.getSubrayLeafNegativeDistancesCount() << "\n" + << "\t\tFurther than current closest: " + << HDA_GV.getSubrayLeafFurtherThanClosestCount() << "\n" + << "\t\tFailed tmin checks: " + << subrayLeafFailedTminCheckCount << "\n" + << "\t\tFailed tmax checks: " + << subrayLeafFailedTmaxCheckCount << "\n" + << "\t\tFailed t checks: " + << ( + raycasterLeafFailedTminCheckCount + + raycasterLeafFailedTmaxCheckCount + ) << "\n"; // Print through info logging level system std::string text = ss.str(); logging::INFO(ss.str()); diff --git a/src/dataanalytics/HDA_PulseRecorder.cpp b/src/dataanalytics/HDA_PulseRecorder.cpp index e8518898d..0b27f94e5 100644 --- a/src/dataanalytics/HDA_PulseRecorder.cpp +++ b/src/dataanalytics/HDA_PulseRecorder.cpp @@ -56,7 +56,7 @@ void HDA_PulseRecorder::recordIntensityCalculation( } } -void HDA_PulseRecorder::recordSubraySimuilation( +void HDA_PulseRecorder::recordSubraySimulation( std::vector const &record ){ std::unique_lock lock(subraySimMutex); diff --git a/src/dataanalytics/HDA_PulseRecorder.h b/src/dataanalytics/HDA_PulseRecorder.h index fa862732d..c8c67e5f1 100644 --- a/src/dataanalytics/HDA_PulseRecorder.h +++ b/src/dataanalytics/HDA_PulseRecorder.h @@ -60,6 +60,49 @@ class HDA_PulseRecorder : public HDA_Recorder{ * [3] -> Circle step * * [4] -> Divergence angle (in rad) + * + * [5] -> Ray direction norm + * + * [6] -> Subray direction norm + * + * [7] -> Angle between ray and subray (in rad) + * + * [8] -> Ray-subray sign check (1 if sign match, 0 otherwise) + * + * [9] -> Min time for subray intersection + * + * [10] -> Max time for subray intersection + * + * [11] -> Count of subray as positive direction during ray tracing + * + * [12] -> Count of subray as negative direction during ray tracing + * + * [13] -> Count of subray as parallel direction during ray tracing + * + * [14] -> Count of subray as not touching the second half of a node during + * ray tracing + * + * [15] -> Count of subray as not touching the first half of a node during + * ray tracing + * + * [16] -> Count of subray as touching both sides of a node during ray + * tracing + * + * [17] -> Count of subrays as touching both sides of a node during a + * second ray tracing try (which uses thit-epsilon for tmin, while the + * first try uses thit+epsilon) + * + * [18] -> Subray direction (x component) + * + * [19] -> Subray direction (y component) + * + * [20] -> Subray direction (z component) + * + * [21] -> Ray direction (x component) + * + * [22] -> Ray direction (y component) + * + * [23] -> Ray direction (z component) */ std::shared_ptr>> subraySim; @@ -129,7 +172,7 @@ class HDA_PulseRecorder : public HDA_Recorder{ /** * @brief Handle all the records for the current subray simulation. */ - virtual void recordSubraySimuilation(std::vector const &record); + virtual void recordSubraySimulation(std::vector const &record); /** * @brief Like * HDA_PulseRecorder::recordSubraySimulation(std::vector) diff --git a/src/scanner/MultiScanner.cpp b/src/scanner/MultiScanner.cpp index e60eaa2a5..db051cbd4 100644 --- a/src/scanner/MultiScanner.cpp +++ b/src/scanner/MultiScanner.cpp @@ -167,7 +167,8 @@ void MultiScanner::computeSubrays( std::map &reflections, vector &intersects #ifdef DATA_ANALYTICS - ,bool &subrayHit + ,bool &subrayHit, + std::vector &subraySimRecord #endif )> handleSubray, vector const &tMinMax, diff --git a/src/scanner/MultiScanner.h b/src/scanner/MultiScanner.h index f04c47341..e4a4edac4 100644 --- a/src/scanner/MultiScanner.h +++ b/src/scanner/MultiScanner.h @@ -166,7 +166,8 @@ class MultiScanner : public Scanner{ std::map &reflections, vector &intersects #ifdef DATA_ANALYTICS - ,bool &subrayHit + ,bool &subrayHit, + std::vector &subraySimRecord #endif )> handleSubray, vector const &tMinMax, diff --git a/src/scanner/Scanner.h b/src/scanner/Scanner.h index 78bb769eb..e8bbc5db3 100644 --- a/src/scanner/Scanner.h +++ b/src/scanner/Scanner.h @@ -451,7 +451,8 @@ class Scanner : public Asset { std::map &reflections, vector &intersects #ifdef DATA_ANALYTICS - ,bool &subrayHit + ,bool &subrayHit, + std::vector &subraySimRecord #endif )> handleSubray, vector const &tMinMax, diff --git a/src/scanner/ScanningDevice.cpp b/src/scanner/ScanningDevice.cpp index 5e61ce83e..214c73624 100644 --- a/src/scanner/ScanningDevice.cpp +++ b/src/scanner/ScanningDevice.cpp @@ -227,7 +227,8 @@ void ScanningDevice::computeSubrays( std::map &reflections, vector &intersects #ifdef DATA_ANALYTICS - ,bool &subrayHit + ,bool &subrayHit, + std::vector &subraySimRecord #endif )> handleSubray, std::vector const &tMinMax, @@ -250,7 +251,8 @@ void ScanningDevice::computeSubrays( double const subrayDivergenceAngle_rad = radiusStep * radiusStep_rad; // Rotate subbeam into divergence step (towards outer rim of the beam cone): - Rotation r1 = Rotation(Directions::right, subrayDivergenceAngle_rad); + Rotation r1 = Rotation(Directions::right, subrayDivergenceAngle_rad); // TODO Restore + //Rotation r1 = Rotation(Directions::right, 0.0); // TODO Remove // Calculate circle step width: int circleSteps = (int)(PI_2 * radiusStep); @@ -264,6 +266,11 @@ void ScanningDevice::computeSubrays( // # Loop over sub-rays along the circle for (int circleStep = 0; circleStep < circleSteps; circleStep++){ +#ifdef DATA_ANALYTICS + std::vector subraySimRecord( + 24, std::numeric_limits::quiet_NaN() + ); +#endif handleSubray( tMinMax, circleStep, @@ -274,18 +281,18 @@ void ScanningDevice::computeSubrays( reflections, intersects #ifdef DATA_ANALYTICS - ,subrayHit + ,subrayHit, + subraySimRecord #endif ); #ifdef DATA_ANALYTICS HDA_GV.incrementGeneratedSubraysCount(); - pulseRecorder->recordSubraySimuilation(std::vector({ - (double)subrayHit, - (double) radiusStep, - (double) circleSteps, - (double) circleStep, - subrayDivergenceAngle_rad - })); + subraySimRecord[0] = (double) subrayHit; + subraySimRecord[1] = (double) radiusStep; + subraySimRecord[2] = (double) circleSteps; + subraySimRecord[3] = (double) circleStep; + subraySimRecord[4] = subrayDivergenceAngle_rad; + pulseRecorder->recordSubraySimulation(subraySimRecord); #endif } } diff --git a/src/scanner/ScanningDevice.h b/src/scanner/ScanningDevice.h index 4e919c1d5..11111ce37 100644 --- a/src/scanner/ScanningDevice.h +++ b/src/scanner/ScanningDevice.h @@ -309,7 +309,8 @@ class ScanningDevice : public Asset { std::map &reflections, vector &intersects #ifdef DATA_ANALYTICS - ,bool &subrayHit + ,bool &subrayHit, + std::vector &subraySimRecord #endif )> handleSubray, std::vector const &tMinMax, diff --git a/src/scanner/SingleScanner.cpp b/src/scanner/SingleScanner.cpp index 1a7c11b80..72d9f3b8e 100644 --- a/src/scanner/SingleScanner.cpp +++ b/src/scanner/SingleScanner.cpp @@ -203,7 +203,8 @@ void SingleScanner::computeSubrays( std::map &reflections, vector &intersects #ifdef DATA_ANALYTICS - ,bool &subrayHit + ,bool &subrayHit, + std::vector &subraySimRecord #endif )> handleSubray, vector const &tMinMax, diff --git a/src/scanner/SingleScanner.h b/src/scanner/SingleScanner.h index 34bb88bcb..1f7134053 100644 --- a/src/scanner/SingleScanner.h +++ b/src/scanner/SingleScanner.h @@ -139,7 +139,8 @@ class SingleScanner : public Scanner{ std::map &reflections, vector &intersects #ifdef DATA_ANALYTICS - ,bool &subrayHit + ,bool &subrayHit, + std::vector &subraySimRecord #endif )> handleSubray, vector const &tMinMax, diff --git a/src/scanner/detector/DynFullWaveformPulseRunnable.cpp b/src/scanner/detector/DynFullWaveformPulseRunnable.cpp index f71b04695..34a0542a5 100644 --- a/src/scanner/detector/DynFullWaveformPulseRunnable.cpp +++ b/src/scanner/detector/DynFullWaveformPulseRunnable.cpp @@ -7,10 +7,17 @@ DynFullWaveformPulseRunnable::findIntersection( vector const &tMinMax, glm::dvec3 const &o, glm::dvec3 const &v +#ifdef DATA_ANALYTICS + ,std::vector &subraySimRecord +#endif ) const{ // Consider scene AABB intersection check was done before at operator() // Thus, proceed to raycasting directly return shared_ptr(raycaster->search( o, v, tMinMax[0], tMinMax[1], false +#ifdef DATA_ANALYTICS + ,subraySimRecord, + true // Subray flag +#endif )); } diff --git a/src/scanner/detector/DynFullWaveformPulseRunnable.h b/src/scanner/detector/DynFullWaveformPulseRunnable.h index 9bea8c46e..772d44a07 100644 --- a/src/scanner/detector/DynFullWaveformPulseRunnable.h +++ b/src/scanner/detector/DynFullWaveformPulseRunnable.h @@ -50,5 +50,8 @@ class DynFullWaveformPulseRunnable : public FullWaveformPulseRunnable { vector const &tMinMax, glm::dvec3 const &o, glm::dvec3 const &v +#ifdef DATA_ANALYTICS + ,std::vector &subraySimRecord +#endif ) const override; }; \ No newline at end of file diff --git a/src/scanner/detector/FullWaveformPulseRunnable.cpp b/src/scanner/detector/FullWaveformPulseRunnable.cpp index aac64abba..c2675480e 100644 --- a/src/scanner/detector/FullWaveformPulseRunnable.cpp +++ b/src/scanner/detector/FullWaveformPulseRunnable.cpp @@ -18,6 +18,9 @@ #ifdef DATA_ANALYTICS #include using helios::analytics::HDA_GV; + +#include +#include #endif using namespace std; @@ -52,11 +55,26 @@ void FullWaveformPulseRunnable::operator()( vector tMinMax = scene.getAABB()->getRayIntersection( pulse.getOriginRef(), beamDir - ); + ); // TODO Restore + // TODO Remove --- + /*std::shared_ptr aabb = std::make_shared(*scene.getAABB()); + aabb->vertices[0].pos.x -= 1.0; + aabb->vertices[0].pos.y -= 1.0; + aabb->vertices[0].pos.z -= 1.0; + aabb->vertices[1].pos.x += 1.0; + aabb->vertices[1].pos.y += 1.0; + aabb->vertices[1].pos.z += 1.0; + aabb->bounds[0] = aabb->vertices[0].pos; + aabb->bounds[1] = aabb->vertices[1].pos; + vector tMinMax = aabb->getRayIntersection( + pulse.getOriginRef(), + beamDir + );*/ + // --- TODO Remove #ifdef DATA_ANALYTICS HDA_GV.incrementGeneratedRaysBeforeEarlyAbortCount(); #endif - if (tMinMax.empty()) { + if (tMinMax.empty() || (tMinMax[0] < 0 && tMinMax[1] < 0)) { logging::DEBUG("Early abort - beam does not intersect with the scene"); scanner->setLastPulseWasHit(false, pulse.getDeviceIndex()); return; @@ -129,7 +147,8 @@ void FullWaveformPulseRunnable::computeSubrays( std::map &reflections, vector &intersects #ifdef DATA_ANALYTICS - ,bool &subrayHit + ,bool &subrayHit, + std::vector &subraySimRecord #endif ) -> void { handleSubray( @@ -143,6 +162,7 @@ void FullWaveformPulseRunnable::computeSubrays( intersects #ifdef DATA_ANALYTICS ,subrayHit, + subraySimRecord, calcIntensityRecords #endif ); @@ -169,6 +189,7 @@ void FullWaveformPulseRunnable::handleSubray( vector &intersects #ifdef DATA_ANALYTICS ,bool &subrayHit, + std::vector &subraySimRecord, std::vector> &calcIntensityRecords #endif ){ @@ -182,6 +203,24 @@ void FullWaveformPulseRunnable::handleSubray( glm::dvec3 subrayDirection = pulse.getAttitude().applyTo(r2) .applyTo(Directions::forward); +#ifdef DATA_ANALYTICS + glm::dvec3 rayDirection = pulse.computeDirection(); + subraySimRecord[5] = glm::l2Norm(rayDirection); // Ray norm + subraySimRecord[6] = glm::l2Norm(subrayDirection); // Subray norm + subraySimRecord[7] = glm::angle( // Angle between ray and subray + rayDirection, + subrayDirection + ); + subraySimRecord[8] = (rayDirection[0] < 0) == (subrayDirection[0] < 0); + subraySimRecord[9] = tMinMax[0]; + subraySimRecord[10] = tMinMax[1]; + subraySimRecord[18] = subrayDirection.x; + subraySimRecord[19] = subrayDirection.y; + subraySimRecord[20] = subrayDirection.z; + subraySimRecord[21] = rayDirection.x; + subraySimRecord[22] = rayDirection.y; + subraySimRecord[23] = rayDirection.z; +#endif glm::dvec3 subrayOrigin(pulse.getOrigin()); bool rayContinues = true; @@ -192,6 +231,9 @@ void FullWaveformPulseRunnable::handleSubray( tMinMax, subrayOrigin, subrayDirection +#ifdef DATA_ANALYTICS + ,subraySimRecord +#endif ); if (intersect != nullptr && intersect->prim != nullptr) { @@ -718,8 +760,17 @@ shared_ptr FullWaveformPulseRunnable::findIntersection( vector const &tMinMax, glm::dvec3 const &o, glm::dvec3 const &v +#ifdef DATA_ANALYTICS + ,std::vector &subraySimRecord +#endif ) const { - return scene.getIntersection(tMinMax, o, v, false); + return scene.getIntersection( + tMinMax, o, v, false +#ifdef DATA_ANALYTICS + ,subraySimRecord, + true // Subray flag +#endif + ); } void FullWaveformPulseRunnable::captureFullWave( diff --git a/src/scanner/detector/FullWaveformPulseRunnable.h b/src/scanner/detector/FullWaveformPulseRunnable.h index 0ca586486..dfa6be5d4 100644 --- a/src/scanner/detector/FullWaveformPulseRunnable.h +++ b/src/scanner/detector/FullWaveformPulseRunnable.h @@ -114,7 +114,8 @@ class FullWaveformPulseRunnable : public AbstractPulseRunnable { vector &intersects #ifdef DATA_ANALYTICS ,bool &subrayHit, - std::vector> &calcIntensityRecords + std::vector &subraySimRecord, + std::vector> &calcIntensityRecords #endif ); /** @@ -273,6 +274,9 @@ class FullWaveformPulseRunnable : public AbstractPulseRunnable { vector const &tMinMax, glm::dvec3 const &o, glm::dvec3 const &v +#ifdef DATA_ANALYTICS + ,std::vector &subraySimRecord +#endif ) const; /** * @brief Detect full waveform peaks diff --git a/src/scene/Scene.cpp b/src/scene/Scene.cpp index b7d943d39..f404ef274 100644 --- a/src/scene/Scene.cpp +++ b/src/scene/Scene.cpp @@ -23,6 +23,11 @@ using namespace glm; #include #include +#ifdef DATA_ANALYTICS +#include +using helios::analytics::HDA_GV; +#endif + using SurfaceInspector::maths::Plane; using SurfaceInspector::maths::PlaneFitter; @@ -91,6 +96,16 @@ bool Scene::finalizeLoading(bool const safe) { // Store original bounding box (CRS coordinates): this->bbox_crs = AABB::getForPrimitives(primitives); + // TODO Rethink : Does this solve the issue ? --- + this->bbox_crs->vertices[0].pos.x -= 1.0; + this->bbox_crs->vertices[0].pos.y -= 1.0; + this->bbox_crs->vertices[0].pos.z -= 1.0; + this->bbox_crs->vertices[1].pos.x += 1.0; + this->bbox_crs->vertices[1].pos.y += 1.0; + this->bbox_crs->vertices[1].pos.z += 1.0; + this->bbox_crs->bounds[0] = this->bbox_crs->vertices[0].pos; + this->bbox_crs->bounds[1] = this->bbox_crs->vertices[1].pos; + // --- TODO Rethink : Does this solve the issue ? glm::dvec3 diff = this->bbox_crs->getMin(); stringstream ss; ss << "CRS bounding box (by vertices): " << this->bbox_crs->toString() @@ -110,6 +125,16 @@ bool Scene::finalizeLoading(bool const safe) { // Get new bounding box of translated scene: this->bbox = AABB::getForPrimitives(primitives); + // TODO Rethink : Does this solve the issue ? --- + this->bbox->vertices[0].pos.x -= 1.0; + this->bbox->vertices[0].pos.y -= 1.0; + this->bbox->vertices[0].pos.z -= 1.0; + this->bbox->vertices[1].pos.x += 1.0; + this->bbox->vertices[1].pos.y += 1.0; + this->bbox->vertices[1].pos.z += 1.0; + this->bbox->bounds[0] = this->bbox->vertices[0].pos; + this->bbox->bounds[1] = this->bbox->vertices[1].pos; + // --- TODO Rethink : Does this solve the issue ? ss << "Actual bounding box (by vertices): " << this->bbox->toString(); logging::INFO(ss.str()); @@ -171,7 +196,14 @@ Scene::getIntersection( bool const groundOnly ) const { vector tMinMax = bbox->getRayIntersection(rayOrigin, rayDir); +#ifdef DATA_ANALYTICS + std::vector subraySimRecord_fake(24, 0); // Not really recorded + return getIntersection( + tMinMax, rayOrigin, rayDir, groundOnly, subraySimRecord_fake + ); +#else return getIntersection(tMinMax, rayOrigin, rayDir, groundOnly); +#endif } std::shared_ptr Scene::getIntersection( @@ -179,13 +211,24 @@ std::shared_ptr Scene::getIntersection( glm::dvec3 const &rayOrigin, glm::dvec3 const &rayDir, bool const groundOnly +#ifdef DATA_ANALYTICS + ,std::vector &subraySimRecord, + bool const isSubray +#endif ) const{ if (tMinMax.empty()) { logging::DEBUG("tMinMax is empty"); +#ifdef DATA_ANALYTICS + HDA_GV.incrementNonIntersectiveSubraysDueToNullTimeCount(); +#endif return nullptr; } return shared_ptr(raycaster->search( rayOrigin, rayDir, tMinMax[0], tMinMax[1], groundOnly +#ifdef DATA_ANALYTICS + ,subraySimRecord, + isSubray +#endif )); } diff --git a/src/scene/Scene.h b/src/scene/Scene.h index 99f7afb82..2481cdb34 100644 --- a/src/scene/Scene.h +++ b/src/scene/Scene.h @@ -217,6 +217,10 @@ class Scene : public Asset { glm::dvec3 const &rayOrigin, glm::dvec3 const &rayDir, bool const groundOnly +#ifdef DATA_ANALYTICS + ,std::vector &subraySimRecord, + bool const isSubray=false +#endif ) const ; /** * @brief Obtain all intersections between the ray and the scene, if any diff --git a/src/scene/primitives/AABB.cpp b/src/scene/primitives/AABB.cpp index 56d2bb250..8154d7a14 100644 --- a/src/scene/primitives/AABB.cpp +++ b/src/scene/primitives/AABB.cpp @@ -163,14 +163,15 @@ std::vector AABB::getRayIntersection( const int ysign = dir.y < 0; const int zsign = dir.z < 0; + // As long as divide by zero yields +-infinity, this logic should work double tmin, tmax; tmin = (bounds[xsign].x - orig.x) / dir.x; tmax = (bounds[1 - xsign].x - orig.x) / dir.x; const double tymin = (bounds[ysign].y - orig.y) / dir.y; const double tymax = (bounds[1 - ysign].y - orig.y) / dir.y; - if (tmin > tymax || tymin > tmax) return std::vector(); - if (tymin > tmin)tmin = tymin; + if (tmin > tymax || tymin > tmax) return {}; + if (tymin > tmin) tmin = tymin; if (tymax < tmax) tmax = tymax; const double tzmin = (bounds[zsign].z - orig.z) / dir.z; diff --git a/src/scene/primitives/AABB.h b/src/scene/primitives/AABB.h index 7fc6478ca..c1e06ddd8 100644 --- a/src/scene/primitives/AABB.h +++ b/src/scene/primitives/AABB.h @@ -128,6 +128,27 @@ class AABB : public Primitive { const glm::dvec3& intersectionPoint ) override; /** + * WARNING! The returned + * intersection times should be interpreted as follows, where \f$t_*\f$ + * is the minimum time and \f$t^*\f$ the maximum: + * + * If \f$t_* \geq 0, t^* \geq 0\f$ then the ray intersects the box twice, + * one time to enter the box, another time to leave it. + * + * If \f$t_* < 0, t^* < 0\f$ then the ray does not intersect the + * box, the line will do, but not the ray which is not allowed to + * go backward (see Convex Optimization by Stephen Boyd and Lieven + * Vandenberghe, where a ray is defined as + * \f$\left\{o + tv : t \geq 0\right\}\f$). + * + * If either \f$t_* < 0\f$ or \f$t^* < 0\f$ then the ray intersects the + * box once, it starts inside the box and the intersection occurs when + * the ray leaves the box. + * + * If the returned vector is empty, then there is no line-box intersection + * at all. + * + * * @see Primitive::getRayIntersection */ std::vector getRayIntersection(