eic-opticks/CSG/csg_intersect_leaf.h

0001 #pragma once
0002 /**
0003 csg_intersect_leaf.h : distance_leaf and intersect_leaf functions
0004 ===================================================================
0005
0006 Thus header needs to be included before csg_intersect_node.h which needs to be included before csg_intersect_tree.h
0007
0008 distance_leaf_sphere
0009 intersect_leaf_sphere
0010     CSG_SPHERE, robust_quadratic_roots
0011
0012 distance_leaf_zsphere
0013 intersect_leaf_zsphere
0014     CSG_ZSPHERE, robust_quadratic_roots
0015
0016 distance_leaf_convexpolyhedron
0017 intersect_leaf_convexpolyhedron
0018     CSG_CONVEXPOLYHEDRON, plane intersections
0019
0020 MISSING : distance_leaf_cone
0021 intersect_leaf_cone
0022     CSG_CONE, newcone with robust_quadratic_roots, oldcone without
0023
0024 MISSING : distance_leaf_hyperboloid
0025 intersect_leaf_hyperboloid
0026     CSG_HYPERBOLOID, robust_quadratic_roots
0027
0028 distance_leaf_box3
0029 intersect_leaf_box3
0030     CSG_BOX3, plane intersections
0031
0032 distance_leaf_plane
0033 intersect_leaf_plane
0034     CSG_PLANE
0035
0036 MISSING : distance_leaf_phicut
0037 intersect_leaf_phicut
0038     CSG_PHICUT
0039
0040 distance_leaf_slab
0041 intersect_leaf_slab
0042     CSG_SLAB
0043
0044 distance_leaf_cylinder
0045 intersect_leaf_cylinder
0046     CSG_CYLINDER, robust_quadratic_roots_disqualifying
0047
0048 MISSING : distance_leaf_infcylinder
0049 intersect_leaf_infcylinder
0050     CSG_INFCYLINDER, robust_quadratic_roots
0051
0052 MISSING : distance_leaf_disc
0053 intersect_leaf_disc
0054     CSG_DISC, disc still using the pseudo-general flop-heavy approach similar to oldcylinder
0055
0056     * TODO: adopt less-flops approach like newcylinder
0057     * (NOT URGENT AS disc NOT CURRENTLY VERY RELEVANT IN ACTIVE GEOMETRIES)
0058
0059 distance_leaf
0060 intersect_leaf
0061
0062 Bringing over functions from  ~/opticks/optixrap/cu/csg_intersect_primitive.h
0063
0064 **/
0065
0066 #include "csg_intersect_leaf_head.h"
0067 #include "OpticksCSG.h"
0068 #include "squad.h"
0069
0070 #include "CSGNode.h"
0071 #include "CSGPrim.h"
0072
0073 #include "csg_robust_quadratic_roots.h"
0074 #include "csg_classify.h"
0075
0076 #if !defined(PRODUCTION) && defined(DEBUG_RECORD)
0077 #include <csignal>
0078 #endif
0079
0080 #if !defined(PRODUCTION) && defined(DEBUG_CYLINDER)
0081 #include "CSGDebug_Cylinder.hh"
0082 #endif
0083
0084 #if !defined(PRODUCTION) && defined(DEBUG_PIDXYZ)
0085 #include <stdio.h>
0086 #endif
0087
0088
0089 #include "csg_intersect_leaf_sphere.h"
0090 #include "csg_intersect_leaf_zsphere.h"
0091 #include "csg_intersect_leaf_cylinder.h"
0092 #include "csg_intersect_leaf_box3.h"
0093 #include "csg_intersect_leaf_newcone.h"
0094 #include "csg_intersect_leaf_convexpolyhedron.h"
0095 #include "csg_intersect_leaf_hyperboloid.h"
0096
0097 #include "csg_intersect_leaf_halfspace.h"
0098 #include "csg_intersect_leaf_phicut.h"
0099
0100 #if !defined(PRODUCTION) && defined(CSG_EXTRA)
0101 #include "csg_intersect_leaf_plane.h"
0102 #include "csg_intersect_leaf_slab.h"
0103 #include "csg_intersect_leaf_thetacut.h"
0104 #include "csg_intersect_leaf_oldcone.h"
0105 #include "csg_intersect_leaf_oldcylinder.h"
0106 #include "csg_intersect_leaf_infcylinder.h"
0107 #include "csg_intersect_leaf_disc.h"
0108 #endif
0109
0110
0111 /**
0112 distance_leaf
0113 ---------------
0114
0115 For hints on how to implement distance functions for more primitives:
0116
0117 * https://www.iquilezles.org/www/articles/distfunctions/distfunctions.htm
0118 * env-;sdf-
0119
0120 **/
0121
0122 LEAF_FUNC
0123 float distance_leaf( const float3& global_position, const CSGNode* node, const float4* plan, const qat4* itra )
0124 {
0125     const unsigned typecode = node->typecode() ;
0126     const unsigned gtransformIdx = node->gtransformIdx() ;
0127     const bool complement = node->is_complement();
0128
0129     const qat4* q = gtransformIdx > 0 ? itra + gtransformIdx - 1 : nullptr ;  // gtransformIdx is 1-based, 0 meaning None
0130
0131     float3 local_position  = q ? q->right_multiply(global_position,  1.f) : global_position ;
0132     float distance = 0.f ;
0133
0134     switch(typecode)
0135     {
0136         case CSG_SPHERE:           distance = distance_leaf_sphere(            local_position, node->q0 )           ; break ;
0137         case CSG_ZSPHERE:          distance = distance_leaf_zsphere(           local_position, node->q0, node->q1 ) ; break ;
0138         case CSG_CYLINDER:         distance = distance_leaf_cylinder(          local_position, node->q0, node->q1 ) ; break ;
0139         case CSG_BOX3:             distance = distance_leaf_box3(              local_position, node->q0 )           ; break ;
0140         case CSG_CONE:             distance = 0.f                                                                   ; break ;
0141         case CSG_CONVEXPOLYHEDRON: distance = distance_leaf_convexpolyhedron(  local_position, node, plan )         ; break ;
0142         case CSG_HYPERBOLOID:      distance = 0.f                                                                   ; break ;
0143         case CSG_PHICUT:           distance = distance_leaf_phicut(            local_position, node->q0 )           ; break ;
0144 #if !defined(PRODUCTION) && defined(CSG_EXTRA)
0145         case CSG_PLANE:            distance = distance_leaf_plane(             local_position, node->q0 )           ; break ;
0146         case CSG_SLAB:             distance = distance_leaf_slab(              local_position, node->q0, node->q1 ) ; break ;
0147         case CSG_OLDCYLINDER:      distance = distance_leaf_cylinder(          local_position, node->q0, node->q1 ) ; break ;
0148 #endif
0149     }
0150
0151     const float sd = complement ? -distance : distance  ;
0152 #if !defined(PRODUCTION) && defined(DEBUG)
0153     printf("//distance_leaf typecode %d name %s complement %d sd %10.4f \n", typecode, CSG::Name(typecode), complement, sd  );
0154 #endif
0155     return sd ;
0156 }
0157
0158
0159 /**
0160 intersect_leaf : must be purely single node
0161 ----------------------------------------------
0162
0163 Notice that only the inverse CSG transforms are needed on the GPU as these are used to
0164 transform the ray_origin and ray_direction into the local origin and direction in the
0165 local frame of the node.
0166
0167 This is called from both::
0168
0169     csg_intersect_node.h
0170     csg_intersect_tree.h
0171
0172 **/
0173
0174 LEAF_FUNC
0175 void intersect_leaf(bool& valid_isect, float4& isect, const CSGNode* node, const float4* plan, const qat4* itra, const float t_min , const float3& ray_origin , const float3& ray_direction, bool dumpxyz )
0176 {
0177     valid_isect = false ;
0178     isect.x = 0.f ;
0179     isect.y = 0.f ;
0180     isect.z = 0.f ;
0181     isect.w = 0.f ;
0182
0183     const unsigned typecode = node->typecode() ;
0184     const unsigned gtransformIdx = node->gtransformIdx() ;
0185     const bool complement = node->is_complement();
0186
0187     const qat4* q = gtransformIdx > 0 ? itra + gtransformIdx - 1 : nullptr ;  // gtransformIdx is 1-based, 0 meaning None
0188
0189     float3 origin    = q ? q->right_multiply(ray_origin,    1.f) : ray_origin ;
0190     float3 direction = q ? q->right_multiply(ray_direction, 0.f) : ray_direction ;
0191
0192 #if !defined(PRODUCTION) && defined(DEBUG_RECORD)
0193     printf("//[intersect_leaf typecode %d name %s gtransformIdx %d \n", typecode, CSG::Name(typecode), gtransformIdx );
0194 #endif
0195
0196 #if !defined(PRODUCTION) && defined(DEBUG)
0197     //printf("//[intersect_leaf typecode %d name %s gtransformIdx %d \n", typecode, CSG::Name(typecode), gtransformIdx );
0198     //printf("//intersect_leaf ray_origin (%10.4f,%10.4f,%10.4f) \n",  ray_origin.x, ray_origin.y, ray_origin.z );
0199     //printf("//intersect_leaf ray_direction (%10.4f,%10.4f,%10.4f) \n",  ray_direction.x, ray_direction.y, ray_direction.z );
0200     /*
0201     if(q)
0202     {
0203         printf("//intersect_leaf q.q0.f (%10.4f,%10.4f,%10.4f,%10.4f)  \n", q->q0.f.x,q->q0.f.y,q->q0.f.z,q->q0.f.w  );
0204         printf("//intersect_leaf q.q1.f (%10.4f,%10.4f,%10.4f,%10.4f)  \n", q->q1.f.x,q->q1.f.y,q->q1.f.z,q->q1.f.w  );
0205         printf("//intersect_leaf q.q2.f (%10.4f,%10.4f,%10.4f,%10.4f)  \n", q->q2.f.x,q->q2.f.y,q->q2.f.z,q->q2.f.w  );
0206         printf("//intersect_leaf q.q3.f (%10.4f,%10.4f,%10.4f,%10.4f)  \n", q->q3.f.x,q->q3.f.y,q->q3.f.z,q->q3.f.w  );
0207         printf("//intersect_leaf origin (%10.4f,%10.4f,%10.4f) \n",  origin.x, origin.y, origin.z );
0208         printf("//intersect_leaf direction (%10.4f,%10.4f,%10.4f) \n",  direction.x, direction.y, direction.z );
0209     }
0210     */
0211 #endif
0212
0213     switch(typecode)
0214     {
0215         case CSG_SPHERE:           intersect_leaf_sphere(           valid_isect, isect, node->q0,               t_min, origin, direction ) ; break ;
0216         case CSG_ZSPHERE:          intersect_leaf_zsphere(          valid_isect, isect, node->q0, node->q1,     t_min, origin, direction ) ; break ;
0217         case CSG_CYLINDER:         intersect_leaf_cylinder(         valid_isect, isect, node->q0, node->q1,     t_min, origin, direction ) ; break ;
0218         case CSG_BOX3:             intersect_leaf_box3(             valid_isect, isect, node->q0,               t_min, origin, direction ) ; break ;
0219         case CSG_CONE:             intersect_leaf_newcone(          valid_isect, isect, node->q0,               t_min, origin, direction ) ; break ;
0220         case CSG_CONVEXPOLYHEDRON: intersect_leaf_convexpolyhedron( valid_isect, isect, node, plan,             t_min, origin, direction ) ; break ;
0221         case CSG_HYPERBOLOID:      intersect_leaf_hyperboloid(      valid_isect, isect, node->q0,               t_min, origin, direction ) ; break ;
0222         case CSG_PHICUT:           intersect_leaf_phicut_simple(    valid_isect, isect, node->q0,               t_min, origin, direction ) ; break ;
0223         case CSG_HALFSPACE:        intersect_leaf_halfspace(        valid_isect, isect, node->q0,               t_min, origin, direction ) ; break ;
0224 #if !defined(PRODUCTION) && defined(CSG_EXTRA)
0225         case CSG_PLANE:            intersect_leaf_plane(            valid_isect, isect, node->q0,               t_min, origin, direction ) ; break ;
0226         case CSG_SLAB:             intersect_leaf_slab(             valid_isect, isect, node->q0, node->q1,     t_min, origin, direction ) ; break ;
0227         case CSG_OLDCYLINDER:      intersect_leaf_oldcylinder(      valid_isect, isect, node->q0, node->q1,     t_min, origin, direction ) ; break ;
0228         case CSG_THETACUT:         intersect_leaf_thetacut(         valid_isect, isect, node->q0, node->q1,     t_min, origin, direction ) ; break ;
0229         case CSG_OLDCONE:          intersect_leaf_oldcone(          valid_isect, isect, node->q0,               t_min, origin, direction ) ; break ;
0230         case CSG_INFCYLINDER:      intersect_leaf_infcylinder(      valid_isect, isect, node->q0, node->q1,     t_min, origin, direction ) ; break ;
0231         case CSG_DISC:             intersect_leaf_disc(             valid_isect, isect, node->q0, node->q1,     t_min, origin, direction ) ; break ;
0232 #endif
0233     }
0234     // NB: changing typecode->imp mapping is a handy way to use old imp with current geometry
0235
0236
0237
0238 #if defined(DEBUG_PIDXYZ)
0239     // HMM MAGIC ACTIVE HERE TOO
0240     //if(dumpxyz) printf("//[intersect_leaf.MAGIC typecode %d valid_isect %d isect (%10.4f %10.4f %10.4f %10.4f) complement %d \n",  typecode, valid_isect, isect.x, isect.y, isect.z, isect.w, complement );
0241     //if(dumpxyz) printf("//[intersect_leaf.MAGIC typecode %d                isect (%10.4f %10.4f %10.4f %10.4f) complement %d \n",  typecode,              isect.x, isect.y, isect.z, isect.w, complement );
0242     //if(dumpxyz) printf("//[intersect_leaf.MAGIC typecode %d \n",  typecode );
0243     //if(dumpxyz) printf("//[intersect_leaf.MAGIC \n");
0244
0245     /**
0246     when applying the MAGIC here just a string will do
0247     **/
0248 #endif
0249
0250
0251 //#define WITH_HEISENBUG 1
0252 #if !defined(WITH_HEISENBUG)
0253
0254    if(valid_isect && q ) q->left_multiply_inplace( isect, 0.f ) ;
0255     // normals transform differently : with inverse-transform-transposed
0256     // so left_multiply the normal by the inverse-transform rather than the right_multiply
0257     // done above to get the inverse transformed origin and direction
0258
0259
0260     /// BIZARRO : RE-EXPRESSING THE MISS-COMPLEMENT-SIGNALLING IMPL FROM THE ABOVE TO A
0261     /// TERSE FORM BELOW (WHICH SHOULD DO EFFECTIVELY THE SAME THING)
0262     /// AVOIDS THE HEISENBUG : NO NEED FOR MAGIC PRINTF IN intersect_leaf
0263
0264     if(complement)
0265     {
0266 #if defined(DEBUG_PIDXYZ)
0267         //if(dumpxyz) printf("// intersect_leaf complement %d valid_isect %d \n", complement, valid_isect );
0268 #endif
0269
0270         // flip normal for hit, signal complement for miss
0271         isect.x = valid_isect ? -isect.x : -0.f ;    // miss needs to signal complement with -0.f signbit
0272         isect.y = valid_isect ? -isect.y : isect.y ; // miss unbounded exit signalled in isect.y for intersect_tree
0273         isect.z = valid_isect ? -isect.z : isect.z ;
0274     }
0275
0276 #else
0277     /// CAUTION : SOMETHING ABOUT THE BELOW MISS-COMPLEMENT-SIGNALLING
0278     /// CODE CAUSES OPTIX 7.5 AND 8.0 HEISENBUG WITH CUDA 12.4 AS REVEALED
0279     /// WITH RTX 5000 Ada GENERATION GPU
0280
0281      if(valid_isect)
0282      {
0283          if(q) q->left_multiply_inplace( isect, 0.f );
0284
0285          if(complement)  // flip normal for complement
0286          {
0287             isect.x = -isect.x ;
0288             isect.y = -isect.y ;
0289             isect.z = -isect.z ;
0290         }
0291     }
0292     else
0293     {
0294         // even for miss need to signal the complement with a -0.f in isect.x
0295         if(complement) isect.x = -0.f ;
0296         // note that isect.y is also flipped for unbounded exit : for consumption by intersect_tree
0297     }
0298 #endif
0299
0300
0301     // NOTICE THAT "VALID_ISECT" IS A BIT MIS-NAMED : AS FALSE JUST MEANS A GEOMETRY MISS
0302
0303 #if !defined(PRODUCTION) && defined(DEBUG_RECORD)
0304     printf("//]intersect_leaf typecode %d name %s valid_isect %d isect (%10.4f %10.4f %10.4f %10.4f)   \n", typecode, CSG::Name(typecode), valid_isect, isect.x, isect.y, isect.z, isect.w);
0305 #endif
0306
0307 #if defined(DEBUG_PIDXYZ)
0308     // BIZARRELY WITH OptiX 7.5.0 CUDA 12.4 "RTX 5000 Ada Generation" : commenting the below line breaks boolean intersects
0309     // NOPE SAME WITH OptiX 8.0.0 CUDA 12.4 "RTX 5000 Ada Generation"
0310
0311     //if(dumpxyz) printf("%d\n", valid_isect );  // HUH : NEED THIS LINE WITH OPTIX 7.5 CUDA 12.4 RTX 5000 ADA
0312     //if(dumpxyz) printf("//]intersect_leaf typecode %d valid_isect %d isect (%10.4f %10.4f %10.4f %10.4f) complement %d \n",  typecode, valid_isect, isect.x, isect.y, isect.z, isect.w, complement );
0313
0314     //if(dumpxyz) printf("//hello\n");
0315     //if(dumpxyz) printf("//]intersect_leaf typecode %d \n", typecode );
0316     //if(dumpxyz) printf("//]intersect_leaf isect (%10.4f %10.4f %10.4f %10.4f) \n", isect.x, isect.y, isect.z, isect.w );
0317     //if(dumpxyz) printf("//]intersect_leaf complement %d \n", complement );
0318
0319     /**
0320     Seems have to potentially dump valid_isect here for the restorative sorcery to work
0321     **/
0322
0323 #endif
0324 }
0325