1 /// Helper functions for working with geometrical shapes
2 module ggplotd.geometry;
3 
4 version( unittest )
5 {
6     import dunit.toolkit;
7 }
8 
9 /// Vertex in 3D space
10 struct Vertex3D
11 {
12     /// x value
13     double x;
14     /// y value
15     double y;
16     /// z value
17     double z;
18 
19     /// Copy constructor
20     this( in Vertex3D v )
21     {
22         x = v.x; y = v.y; z = v.z;
23     }
24 
25     /// Constructor taking x, y and z value
26     this( in double _x, in double _y, in double _z )
27     {
28         x = _x; y = _y; z = _z;
29     }
30 
31     /// Vertex substraction or addition
32     Vertex3D opBinary(string s)( in Vertex3D v2 ) const if (s == "-" || s == "+" )
33     {
34         import std.format : format;
35         mixin( format("return Vertex3D( x %s v2.x, y %s v2.y, z %s v2.z );",s,s,s) );
36     }
37 }
38 unittest
39 {
40     immutable v1 = Vertex3D(1,2,3);
41     immutable v2 = Vertex3D(3,2,1);
42     auto v =  v1-v2;
43     assertEqual( v, Vertex3D(-2,0,2) );
44     v =  v1+v2;
45     assertEqual( v, Vertex3D(4,4,4) );
46 }
47 
48 private Vertex3D crossProduct( in Vertex3D v1, in Vertex3D v2 )
49 {
50     return Vertex3D( v1.y*v2.z-v1.z*v2.y, 
51         v1.z*v2.x-v1.x*v2.z, 
52         v1.x*v2.y-v1.y*v2.x );
53 }
54 
55 unittest
56 {
57     auto v = Vertex3D(1,2,3).crossProduct( Vertex3D( 3,2,1 ) );
58     assertEqual( v, Vertex3D(-4,8,-4) );
59 }
60 
61 /// Calculate gradientVector based on a Triangle. The vertices of the triangle are assumed to be sorted by height.
62 auto gradientVector( T )( in T triangle )
63 {
64     import std.math : pow;
65     assert( triangle[0].z <= triangle[1].z && triangle[1].z <= triangle[2].z,
66         "gradientVector expects the triangle vertices to be sorted by height" );
67     auto gVector = [ Vertex3D(triangle[0]), Vertex3D(triangle[2]) ];
68 
69     if (triangle[0].z == triangle[2].z) {
70         return gVector;
71     }
72 
73     auto e1 = triangle[2]-triangle[0];
74     auto e2 = triangle[1]-triangle[0];
75     auto normal = e1.crossProduct( e2 );
76 
77     Vertex3D v;
78     if (normal.x == 0)
79     {
80         v = Vertex3D( normal.x/normal.y, 1, 
81             -(pow(normal.y,2)+pow(normal.x,2))/(normal.y*normal.z)
82             );
83     } else 
84     {
85         v = Vertex3D( 1, normal.y/normal.x, 
86             -(pow(normal.y,2)+pow(normal.x,2))/(normal.x*normal.z)
87             );
88     }
89 
90     immutable scalar = (gVector[1].z-gVector[0].z)/v.z;
91 
92     gVector[1].x = gVector[0].x + scalar*v.x;
93     gVector[1].y = gVector[0].y + scalar*v.y;
94     gVector[1].z = gVector[0].z + scalar*v.z;
95 
96     return gVector; 
97 }
98 
99 unittest
100 {
101     auto v = gradientVector( [ Vertex3D( 0, 0, 1 ),
102         Vertex3D( 0, 1, 1 ),
103         Vertex3D( 1, 1, 1 )] );
104     assertEqual( v[0].z, v[1].z );
105 
106     v = gradientVector( [ Vertex3D( 0, 0, 0 ),
107         Vertex3D( 0, 1, 0.5 ),
108         Vertex3D( 1, 1, 1 )] );
109     assertEqual( v, [Vertex3D(0, 0, 0), Vertex3D(1, 1, 1)] );
110 
111     v = gradientVector( [ Vertex3D( 0, 0, 0 ),
112         Vertex3D( 0, 1, 0.1 ),
113         Vertex3D( 1, 1, 1 )] );
114     assertEqual( v[0], Vertex3D(0, 0, 0) );
115     assertApprox( v[1].x, 1.09756, 1e-5 );
116     assertApprox( v[1].y, 0.121951, 1e-5 );
117     assertApprox( v[1].z, 1.0 );
118 }
119