fix line source

fangq · Jun 4, 2022 · 8af3631 · 8af3631
1 parent 149b1ef
commit 8af3631
Showing 1 changed file with 6 additions and 5 deletions.
diff --git a/src/mcx_core.cu b/src/mcx_core.cu
@@ -1350,11 +1350,12 @@ __device__ inline int launchnewphoton(MCXpos* p, MCXdir* v, Stokes* s, MCXtime*
                                            p->w);
 
                     if (gcfg->srctype == MCX_SRC_LINE) {
-                        float s, q;
-                        r = 1.f - 2.f * rand_uniform01(t);
-                        s = 1.f - 2.f * rand_uniform01(t);
-                        q = sqrtf(1.f - v->x * v->x - v->y * v->y) * (rand_uniform01(t) > 0.5f ? 1.f : -1.f);
-                        *((float4*)v) = float4(v->y * q - v->z * s, v->z * r - v->x * q, v->x * s - v->y * r, v->nscat);
+                        float sphi, cphi;
+                        r = rsqrtf(gcfg->srcparam1.x * gcfg->srcparam1.x + gcfg->srcparam1.y * gcfg->srcparam1.y + gcfg->srcparam1.z * gcfg->srcparam1.z);
+                        *((float4*)v) = float4(gcfg->srcparam1.x * r, gcfg->srcparam1.y * r, gcfg->srcparam1.z * r, v->nscat);
+                        r = TWO_PI * rand_uniform01(t); // phi
+                        sincosf(r, &sphi, &cphi); // y=sin(phi), x=cos(phi)
+                        rotatevector(v, 1.f, 0.f, sphi, cphi);
                     }
 
                     *rv = float3(rv->x + (gcfg->srcparam1.x) * 0.5f,