Hello Shailesh,
Your patch needed to be manually merged, because it had merge conflicts
with the latest version at SVN HEAD.
Other than that you can't use just any value for ID_HRT. You need to use
the value declared in "include/rt/defines.h". That produces visible output.
But there are still problems with the intersections. Depending on the view
direction the output can be completely out of whack. It is probably some
issue either with the orientation matrix/computations or with the equation
solver.
Regards,
---------- Forwarded message ----------
From: Shailesh Tripathi <shailesh.tripathi.ec...@itbhu.ac.in>
Date: Thu, Apr 13, 2017 at 12:12 AM
Subject: Re: request for some help
To: Vasco Alexandre da Silva Costa <vasco.co...@gmail.com>
...
I am attaching my patch containing the latest changes.
The commands used:
>in hrt.s hrt 0 0 0 1 0 0 0 1 0 0 0 1 14
>rt -z1
Also, please tell me is there a log file generated with the rt log, as I am
not being able to copy my log to ask it on the mailing list.
--
Vasco Alexandre da Silva Costa
PhD in Computer Engineering (Computer Graphics)
Instituto Superior Técnico/University of Lisbon, Portugal
Index: src/librt/primitives/primitive_util.c
===================================================================
--- src/librt/primitives/primitive_util.c (revision 69522)
+++ src/librt/primitives/primitive_util.c (working copy)
@@ -627,6 +627,7 @@
"rec_shot.cl",
"tgc_shot.cl",
"tor_shot.cl",
+ "hrt_shot.cl",
"rt.cl",
};
@@ -689,6 +690,7 @@
case ID_BOT: size = clt_bot_pack(pool, stp); break;
case ID_EPA: size = clt_epa_pack(pool, stp); break;
case ID_ETO: size = clt_eto_pack(pool, stp); break;
+ case ID_HRT: size = clt_hrt_pack(pool, stp); break;
default: size = 0; break;
}
return size;
Index: src/librt/primitives/hrt/hrt_shot.cl
===================================================================
--- src/librt/primitives/hrt/hrt_shot.cl (revision 0)
+++ src/librt/primitives/hrt/hrt_shot.cl (revision 0)
@@ -0,0 +1,238 @@
+#include "common.cl"
+
+
+struct hrt_specific {
+ double hrt_V[3]; /* Vector to center of heart */
+ double hrt_SoR[16]; /* Scale(Rot(vect)) */
+};
+
+
+int
+hrt_shot(RESULT_TYPE *res, const double3 r_pt, const double3 r_dir, const uint idx, global const struct hrt_specific *hrt)
+{
+ double3 dprime; /* D' : The new shot direction */
+ double3 pprime; /* P' : The new shot point */
+ double3 trans; /* Translated shot vector */
+ double3 norm_pprime; /* P' with normalized distance from heart */
+ double Xsqr[3]; /* X^2 */
+ double Ysqr[3]; /* 9/4*Y^2 */
+ double Zsqr[3]; /* Z^2 - 1 */
+ double Acube[7], S[7]; /* The sextic equation (of power 6) and A^3 */
+ double Zcube[4], A[3]; /* Z^3 and X^2 + 9/4 * Y^2 + Z^2 - 1 */
+ double X2_Y2[3], Z3_X2_Y2[6]; /* X^2 + 9/80*Y^2 and Z^3*(X^2 + 9/80*Y^2) */
+ bn_complex_t complex[6]; /* The complex roots */
+ double real[6]; /* The real roots */
+ int i;
+ int j;
+ int npts;
+
+ /* Translate the ray point */
+ trans = r_pt - vload3(0, hrt->hrt_V);
+
+ /* Scale and Rotate point to get P' */
+ pprime.x = (hrt->hrt_SoR[0]*trans.x + hrt->hrt_SoR[1]*trans.y + hrt->hrt_SoR[2]*trans.z) * 1.0/(hrt->hrt_SoR[15]);
+ pprime.y = (hrt->hrt_SoR[4]*trans.x + hrt->hrt_SoR[5]*trans.y + hrt->hrt_SoR[6]*trans.z) * 1.0/(hrt->hrt_SoR[15]);
+ pprime.z = (hrt->hrt_SoR[8]*trans.x + hrt->hrt_SoR[9]*trans.y + hrt->hrt_SoR[10]*trans.z) * 1.0/(hrt->hrt_SoR[15]);
+ /* Translate ray direction vector */
+ dprime = MAT4X3VEC(hrt->hrt_SoR, r_dir);
+ dprime = normalize(dprime);
+
+ /* Normalize distance from the heart. Substitutes a corrected ray
+ * point, which contains a translation along the ray direction to
+ * the closest approach to vertex of the heart.Translating the ray
+ * along the direction of the ray to the closest point near the
+ * heart's center vertex. Thus, the New ray origin is normalized.
+ */
+ norm_pprime = dprime * dot(pprime, dprime);
+ norm_pprime = pprime - norm_pprime;
+
+ /**
+ * Generate the sextic equation S(t) = 0 to be passed through the root finder.
+ */
+ /* X**2 */
+ Xsqr[0] = dprime.x * dprime.x;
+ Xsqr[1] = 2 * dprime.x * pprime.x;
+ Xsqr[2] = pprime.x * pprime.x;
+
+ /* 9/4 * Y**2*/
+ Ysqr[0] = 9/4 * dprime.y * dprime.y;
+ Ysqr[1] = 9/2 * dprime.y * pprime.y;
+ Ysqr[2] = 9/4 * (pprime.y * pprime.y);
+
+ /* Z**2 - 1 */
+ Zsqr[0] = dprime.z * dprime.z;
+ Zsqr[1] = 2 * dprime.z * pprime.z;
+ Zsqr[2] = pprime.z * pprime.z - 1.0 ;
+
+ /* A = X^2 + 9/4 * Y^2 + Z^2 - 1 */
+ A[0] = Xsqr[0] + Ysqr[0] + Zsqr[0];
+ A[1] = Xsqr[1] + Ysqr[1] + Zsqr[1];
+ A[2] = Xsqr[2] + Ysqr[2] + Zsqr[2];
+
+ /* Z**3 */
+ Zcube[0] = dprime.z * Zsqr[0];
+ Zcube[1] = 1.5 * dprime.z * Zsqr[1];
+ Zcube[2] = 1.5 * pprime.z * Zsqr[1];
+ Zcube[3] = pprime.z * ( Zsqr[2] + 1.0 );
+
+ Acube[0] = A[0] * A[0] * A[0];
+ Acube[1] = 3.0 * A[0] * A[0] * A[1];
+ Acube[2] = 3.0 * (A[0] * A[0] * A[2] + A[0] * A[1] * A[1]);
+ Acube[3] = 6.0 * A[0] * A[1] * A[2] + A[1] * A[1] * A[1];
+ Acube[4] = 3.0 * (A[0] * A[2] * A[2] + A[1] * A[1] * A[2]);
+ Acube[5] = 3.0 * A[1] * A[2] * A[2];
+ Acube[6] = A[2] * A[2] * A[2];
+
+ /* X**2 + 9/80 Y**2 */
+ X2_Y2[0] = Xsqr[0] + Ysqr[0] / 20 ;
+ X2_Y2[1] = Xsqr[1] + Ysqr[1] / 20 ;
+ X2_Y2[2] = Xsqr[2] + Ysqr[2] / 20 ;
+
+ /* Z**3 * (X**2 + 9/80 * Y**2) */
+ Z3_X2_Y2[0] = Zcube[0] * X2_Y2[0];
+ Z3_X2_Y2[1] = X2_Y2[0] * Zcube[1];
+ Z3_X2_Y2[2] = X2_Y2[0] * Zcube[2] + X2_Y2[1] * Zcube[0] + X2_Y2[1] * Zcube[1] + X2_Y2[2] * Zcube[0];
+ Z3_X2_Y2[3] = X2_Y2[0] * Zcube[3] + X2_Y2[1] * Zcube[2] + X2_Y2[2] * Zcube[1];
+ Z3_X2_Y2[4] = X2_Y2[1] * Zcube[3] + X2_Y2[2] * Zcube[2];
+ Z3_X2_Y2[5] = X2_Y2[2] * Zcube[3];
+
+ S[0] = Acube[0];
+ S[1] = Acube[1] - Z3_X2_Y2[0];
+ S[2] = Acube[2] - Z3_X2_Y2[1];
+ S[3] = Acube[3] - Z3_X2_Y2[2];
+ S[4] = Acube[4] - Z3_X2_Y2[3];
+ S[5] = Acube[5] - Z3_X2_Y2[4];
+ S[6] = Acube[6] - Z3_X2_Y2[5];
+
+ /* It is known that the equation is sextic (of order 6). Therefore, if the
+ * root finder returns other than 6 roots, return an error.
+ */
+ if ((i = rt_poly_roots(S, 6, complex)) != 6) {
+ return 0; /* MISS */
+ }
+
+ /* Only real roots indicate an intersection in real space.
+ *
+ * Look at each root returned; if the imaginary part is zero or
+ * sufficiently close, then use the real part as one value of 't'
+ * for the intersections
+ */
+ for (j = 0, i = 0; j < 6; j++) {
+ if (NEAR_ZERO(complex[j].im, RT_PCOEF_TOL))
+ real[i++] = complex[j].re;
+ }
+ /* Here, 'i' is number of points found */
+ switch (i) {
+ case 0:
+ return 0; /* No hit */
+
+ default:
+ // bu_log("rt_hrt_shot: reduced 6 to %d roots\n", i);
+ // bn_pr_roots(stp->st_name, complex, 6);
+ return 0; /* No hit */
+
+ case 2:
+ {
+ /* Sort most distant to least distant. */
+ double u;
+ if ((u=real[0]) < real[1]) {
+ /* bubble larger towards [0] */
+ real[0] = real[1];
+ real[1] = u;
+ }
+ }
+ break;
+ case 4:
+ {
+ short n;
+ short lim;
+
+ /* Inline rt_pt_sort(). Sorts real[] into descending order. */
+ for (lim = i-1; lim > 0; lim--) {
+ for (n = 0; n < lim; n++) {
+ double u;
+ if ((u=real[n]) < real[n+1]) {
+ /* bubble larger towards [0] */
+ real[n] = real[n+1];
+ real[n+1] = u;
+ }
+ }
+ }
+ }
+ break;
+ case 6:
+ {
+ short num;
+ short limit;
+
+ /* Inline rt_pt_sort(). Sorts real[] into descending order. */
+ for (limit = i-1; limit > 0; limit--) {
+ for (num = 0; num < limit; num++) {
+ double u;
+ if ((u=real[num]) < real[num+1]) {
+ /* bubble larger towards [0] */
+ real[num] = real[num+1];
+ real[num+1] = u;
+ }
+ }
+ }
+ }
+ break;
+ }
+
+ npts = i;
+ /* Now, t[0] > t[npts-1] */
+ /* real[1] is entry point, and real[0] is farthest exit point */
+ /* real[3] is entry point, and real[2] is exit point */
+ /* real[5] is entry point, and real[4] is exit point */
+ for (i=npts-1; i>0; i -= 2) {
+ struct hit hits[2];
+
+ hits[0].hit_dist = real[i];
+ hits[0].hit_surfno = 0;
+ /* Set aside vector for rt_tor_norm() later */
+ hits[0].hit_vpriv = pprime + real[i] * dprime;
+
+ hits[1].hit_dist = real[i-1];
+ hits[1].hit_surfno = 0;
+ /* Set aside vector for rt_tor_norm() later */
+ hits[1].hit_vpriv = pprime + real[i-1] * dprime;
+
+ do_segp(res, idx, &hits[0], &hits[1]);
+ }
+
+ return npts;
+}
+
+
+
+void
+hrt_norm(register struct hit *hitp, const double3 r_pt, const double3 r_dir, global const struct hrt_specific *hrt)
+{
+
+ double w, fx, fy, fz;
+ double3 work;
+
+ hitp->hit_point = r_pt + r_dir * hitp->hit_dist;
+ w = hitp->hit_vpriv.x * hitp->hit_vpriv.x
+ + 9.0/4.0 * hitp->hit_vpriv.y * hitp->hit_vpriv.y
+ + hitp->hit_vpriv.z * hitp->hit_vpriv.z - 1.0;
+ fx = (w * w - 1/3 * hitp->hit_vpriv.z * hitp->hit_vpriv.z * hitp->hit_vpriv.z) * hitp->hit_vpriv.x;
+ fy = hitp->hit_vpriv.y * (12/27 * w * w - 80/3 * hitp->hit_vpriv.z * hitp->hit_vpriv.z * hitp->hit_vpriv.z);
+ fz = (w * w - 0.5 * hitp->hit_vpriv.z * (hitp->hit_vpriv.x * hitp->hit_vpriv.x + 9/80 * hitp->hit_vpriv.y * hitp->hit_vpriv.y)) * hitp->hit_vpriv.z;
+ work = (double3){
+ fx,
+ fy,
+ fz};
+ hitp->hit_normal = work;
+}
+
+/*
+ * Local Variables:
+ * mode: C
+ * tab-width: 8
+ * indent-tabs-mode: t
+ * c-file-style: "stroustrup"
+ * End:
+ * ex: shiftwidth=4 tabstop=8
+ */
Property changes on: src/librt/primitives/hrt/hrt_shot.cl
___________________________________________________________________
Added: svn:executable
+ *
Index: src/librt/primitives/hrt/hrt.c
===================================================================
--- src/librt/primitives/hrt/hrt.c (revision 69522)
+++ src/librt/primitives/hrt/hrt.c (working copy)
@@ -147,7 +147,34 @@
mat_t hrt_invR; /* transposed rotation matrix */
};
+#ifdef USE_OPENCL
+/* largest data members first */
+struct clt_hrt_specific {
+ /* REMAINING ELEMENTS PROVIDED BY IMPORT, UNUSED BY EXPORT */
+ cl_double hrt_V[3]; /* Vector to center of heart */
+ cl_double hrt_SoR[16]; /* Scale(Rot(vect)) */
+};
+
+
+size_t
+clt_hrt_pack(struct bu_pool *pool, struct soltab *stp)
+{
+ struct hrt_specific *hrt =
+ (struct hrt_specific *)stp->st_specific;
+ struct clt_hrt_specific *args;
+
+ const size_t size = sizeof(*args);
+ args = (struct clt_hrt_specific*)bu_pool_alloc(pool, 1, size);
+
+ VMOVE(args->hrt_V, hrt->hrt_V);
+
+ MAT_COPY(args->hrt_SoR, hrt->hrt_SoR);
+ return size;
+}
+
+#endif /* USE_OPENCL */
+
/**
* Compute the bounding RPP for a heart.
*/
Index: src/librt/primitives/common.cl
===================================================================
--- src/librt/primitives/common.cl (revision 69522)
+++ src/librt/primitives/common.cl (working copy)
@@ -24,6 +24,8 @@
#define RT_PCOEF_TOL (1.0e-10)
#define RT_DOT_TOL (0.001)
+#define BN_MAX_POLY_DEGREE 6 /* Maximum Poly Order */
+
#define MAX_FASTF (1.0e73) /* Very close to the largest number */
#define SQRT_MAX_FASTF (1.0e36) /* This squared just avoids overflow */
#define SMALL_FASTF (1.0e-77)
@@ -42,7 +44,6 @@
struct { double re, im; };
} bn_complex_t;
-
struct hit {
double3 hit_point;
double3 hit_normal;
@@ -125,6 +126,7 @@
struct sph_specific;
struct tgc_specific;
struct tor_specific;
+struct hrt_specific;
extern int arb_shot(RESULT_TYPE *res, const double3 r_pt, const double3 r_dir, const uint idx, global const struct arb_specific *arb);
extern int bot_shot(RESULT_TYPE *res, const double3 r_pt, double3 r_dir, const uint idx, global const uchar *args);
@@ -137,6 +139,7 @@
extern int sph_shot(RESULT_TYPE *res, const double3 r_pt, const double3 r_dir, const uint idx, global const struct sph_specific *sph);
extern int tgc_shot(RESULT_TYPE *res, const double3 r_pt, const double3 r_dir, const uint idx, global const struct tgc_specific *tgc);
extern int tor_shot(RESULT_TYPE *res, const double3 r_pt, const double3 r_dir, const uint idx, global const struct tor_specific *tor);
+extern int hrt_shot(RESULT_TYPE *res, const double3 r_pt, const double3 r_dir, const uint idx, global const struct hrt_specific *hrt);
extern void arb_norm(struct hit *hitp, const double3 r_pt, const double3 r_dir, global const struct arb_specific *arb);
extern void bot_norm(struct hit *hitp, const double3 r_pt, const double3 r_dir, global const uchar *args);
@@ -149,6 +152,7 @@
extern void sph_norm(struct hit *hitp, const double3 r_pt, const double3 r_dir, global const struct sph_specific *sph);
extern void tgc_norm(struct hit *hitp, const double3 r_pt, const double3 r_dir, global const struct tgc_specific *tgc);
extern void tor_norm(struct hit *hitp, const double3 r_pt, const double3 r_dir, global const struct tor_specific *tor);
+extern void hrt_norm(struct hit *hitp, const double3 r_pt, const double3 r_dir, global const struct hrt_specific *hrt);
#endif /* COMMON_CL */
Index: src/librt/primitives/rt.cl
===================================================================
--- src/librt/primitives/rt.cl (revision 69522)
+++ src/librt/primitives/rt.cl (working copy)
@@ -94,6 +94,7 @@
#define ID_EHY 20 /**< @brief Elliptical Hyperboloid */
#define ID_ETO 21 /**< @brief Elliptical Torus */
#define ID_BOT 30 /**< @brief Bag o' triangles */
+#define ID_HRT 32 /**< @brief heart */
inline int shot(RESULT_TYPE *res, const double3 r_pt, const double3 r_dir, const uint idx, const int id, global const void *args)
@@ -111,6 +112,7 @@
case ID_BOT: return bot_shot(res, r_pt, r_dir, idx, args);
case ID_EPA: return epa_shot(res, r_pt, r_dir, idx, args);
case ID_ETO: return eto_shot(res, r_pt, r_dir, idx, args);
+ case ID_HRT: return hrt_shot(res, r_pt, r_dir, idx, args);
default: return 0;
};
}
@@ -130,6 +132,7 @@
case ID_BOT: bot_norm(hitp, r_pt, r_dir, args); break;
case ID_EPA: epa_norm(hitp, r_pt, r_dir, args); break;
case ID_ETO: eto_norm(hitp, r_pt, r_dir, args); break;
+ case ID_HRT: hrt_norm(hitp, r_pt, r_dir, args); break;
default: break;
};
}
Index: src/librt/librt_private.h
===================================================================
--- src/librt/librt_private.h (revision 69522)
+++ src/librt/librt_private.h (working copy)
@@ -229,6 +229,7 @@
CLT_DECLARE_INTERFACE(ehy);
CLT_DECLARE_INTERFACE(bot);
CLT_DECLARE_INTERFACE(eto);
+CLT_DECLARE_INTERFACE(hrt);
extern size_t clt_bot_pack(struct bu_pool *pool, struct soltab *stp);
#endif
Index: src/librt/CMakeLists.txt
===================================================================
--- src/librt/CMakeLists.txt (revision 69522)
+++ src/librt/CMakeLists.txt (working copy)
@@ -279,6 +279,7 @@
primitives/sph/sph_shot.cl
primitives/tgc/tgc_shot.cl
primitives/tor/tor_shot.cl
+ primitives/hrt/hrt_shot.cl
primitives/xxx/xxx.h
search.h
subd_test_bot.asc
@@ -309,6 +310,7 @@
primitives/rec/rec_shot.cl
primitives/tgc/tgc_shot.cl
primitives/tor/tor_shot.cl
+ primitives/hrt/hrt_shot.cl
)
BRLCAD_ADDDATA(nurb_example.c sample_applications)
------------------------------------------------------------------------------
Check out the vibrant tech community on one of the world's most
engaging tech sites, Slashdot.org! http://sdm.link/slashdot
_______________________________________________
BRL-CAD Developer mailing list
brlcad-devel@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/brlcad-devel
