1 /**
2 Flexible vector template with some math utils. 
3 
4 By default defines glsl style aliases 
5 */
6 module sily.vector;
7 
8 import std.math;
9 import std.numeric;
10 import std.conv;
11 import std.traits;
12 import std.typetuple;
13 import std.algorithm;
14 import std.stdio;
15 import std.string;
16 import std.format;
17 
18 import sily.meta.swizzle;
19 import sily.math;
20 import sily.array;
21 import sily.color;
22 import sily.matrix;
23 import sily.quat;
24 
25 /// GLSL style alias
26 alias vec(T, size_t N) = Vector!(T, N);
27 /// Ditto
28 alias vec(size_t N) = Vector!(float, N);
29 /// Ditto
30 alias vec2 = vec!2;
31 /// Ditto
32 alias vec3 = vec!3;
33 /// Ditto
34 alias vec4 = vec!4;
35 
36 /// Ditto
37 alias dvec(size_t N) = Vector!(double, N);
38 /// Ditto
39 alias dvec2 = dvec!2;
40 /// Ditto
41 alias dvec3 = dvec!3;
42 /// Ditto
43 alias dvec4 = dvec!4;
44 
45 /// Ditto
46 alias ivec(size_t N) = Vector!(int, N);
47 /// Ditto
48 alias ivec2 = ivec!2;
49 /// Ditto
50 alias ivec3 = ivec!3;
51 /// Ditto
52 alias ivec4 = ivec!4;
53 
54 /// Ditto
55 alias uvec(size_t N) = Vector!(uint, N);
56 /// Ditto
57 alias uvec2 = uvec!2;
58 /// Ditto
59 alias uvec3 = uvec!3;
60 /// Ditto
61 alias uvec4 = uvec!4;
62 
63 /++
64 Vector structure with data accesible with `[N]` or swizzling.
65 
66 All operations on Vector (*, +, /, -, etc...) are scalar.
67 
68 Allows casting to sily.color or sily.matrix.
69 +/
70 struct Vector(T, size_t N) if (isNumeric!T && N > 0)  {
71     /// Vector data
72     public T[N] data = [ 0 ];
73 
74     /// Alias to allow easy `data` access
75     alias data this;
76     /// Alias to data type (e.g. float, int)
77     alias dataType = T;
78     /** 
79     Alias to vector type. Can be used to contruct vectors
80     of same type
81     ---
82     auto rvec7 = Vector!(real, 7)(10);
83     auto rvec7s = rvec7.VecType(20);
84     ---
85     */
86     alias VecType = Vector!(T, N);
87     /// Alias to vector size
88     enum size_t size = N;
89     
90     /**
91     Constructs Vector from components. If no components present
92     vector will be filled with 0
93     Example:
94     ---
95     // Vector can be constructed manually or with aliases
96     auto v2 = ivec2(10, 20);
97     // Also vector can be given only one value,
98     // in that case it'll be filled with that value
99     auto v5 = ivec4(13);
100     auto v6 = vec4(0.3f);
101     // Vector values can be accessed with array slicing,
102     // by using color symbols or swizzling
103     float v6x = v6.x;
104     float v6z = v6.z;
105     float[] v6yzx = v6.yzx;
106     float v6y = v6[1];
107     // Valid vector accessors are:
108     // Vector2 - [x, y], [w, h], [u, v]
109     // Vector3 - [x, y, z], [w, h, d], [u, v, t], [r, g, b]
110     // Vector4 - [x, y, z, w], [r, g, b, a]
111     // Other sizes must be accessed with index
112     ---
113     */
114     this(in T val) {
115         foreach (i; 0 .. size) { data[i] = val; }
116     }
117     /// Ditto
118     this(in T[N] vals...) {
119         data = vals;
120     }
121 
122     static if (N == 3 && isFloatingPoint!T) {
123         /// Construct euler vector from quaternion
124         this(in T[4] vals...) {
125             T x0 = 2.0 * (vals[3] * vals[0] + vals[1] * vals[2]);
126             T x1 = 1.0 - 2.0 * (vals[0] * vals[0] + vals[1] * vals[1]);
127             
128             T y0 = 2.0 * (vals[3] * vals[1] - vals[2] * vals[0]);
129             y0 = y0.clamp(-1.0, 1.0);
130 
131             T z0 = 2.0 * (vals[3] * vals[2] + vals[0] * vals[1]);
132             T z1 = 1.0 - 2.0 * (vals[1] * vals[1] + vals[2] * vals[2]);
133 
134             data = [
135                 atan2(x0, x1),
136                 asin(y0),
137                 atan2(z0, z1)
138             ];
139         }
140     }
141 
142     /* -------------------------------------------------------------------------- */
143     /*                         UNARY OPERATIONS OVERRIDES                         */
144     /* -------------------------------------------------------------------------- */
145     
146     /// opBinary x [+, -, *, /, %] y
147     VecType opBinary(string op, R)(in Vector!(R, N) b) const if ( isNumeric!R ) {
148         // assert(/* this !is null && */ b !is null, "\nOP::ERROR nullptr Vector!" ~ size.to!string ~ ".");
149         VecType ret = VecType();
150         foreach (i; 0 .. size) { mixin( "ret[i] = data[i] " ~ op ~ " b.data[i];" ); }
151         return ret;
152     }
153 
154     /// Ditto
155     VecType opBinaryRight(string op, R)(in Vector!(R, N) b) const if ( isNumeric!R ) {
156         // assert(/* this !is null && */ b !is null, "\nOP::ERROR nullptr Vector!" ~ size.to!string ~ ".");
157         VecType ret = VecType();
158         foreach (i; 0 .. size) { mixin( "ret[i] = b.data[i] " ~ op ~ " data[i];" ); }
159         return ret;
160     }
161 
162     /// Ditto
163     VecType opBinary(string op, R)(in R b) const if ( isNumeric!R ) {
164         // assert(this !is null, "\nOP::ERROR nullptr Vector!" ~ size.to!string ~ ".");
165         VecType ret = VecType();
166         foreach (i; 0 .. size) { mixin( "ret[i] = data[i] " ~ op ~ " b;" ); }
167         return ret;
168     }
169 
170     /// Ditto
171     VecType opBinaryRight(string op, R)(in R b) const if ( isNumeric!R ) {
172         // assert(this !is null, "\nOP::ERROR nullptr Vector!" ~ size.to!string ~ ".");
173         VecType ret;
174         foreach (i; 0 .. size) { mixin( "ret[i] = b " ~ op ~ " data[i];" ); }
175         return ret;
176     }
177 
178     /// opEquals x == y
179     bool opEquals(R)(in Vector!(R, size) b) const if ( isNumeric!R ) {
180         // assert(/* this !is null && */ b !is null, "\nOP::ERROR nullptr Vector!" ~ size.to!string ~ ".");
181         bool eq = true;
182         foreach (i; 0 .. size) { 
183             eq = eq && data[i] == b.data[i]; 
184             if (!eq) break;
185         }
186         return eq;
187     }
188 
189     /// opCmp x [< > <= >=] y
190     int opCmp(R)(in Vector!(R, N) b) const if ( isNumeric!R ) {
191         // assert(/* this !is null && */ b !is null, "\nOP::ERROR nullptr Vector!" ~ size.to!string ~ ".");
192         double al = cast(double) length();
193         double bl = cast(double) b.length();
194         if (al == bl) return 0;
195         if (al < bl) return -1;
196         return 1;
197     }
198 
199     /// opUnary [-, +, --, ++] x
200     VecType opUnary(string op)() if (op == "-" || op == "+") {
201         // assert(this !is null, "\nOP::ERROR nullptr Vector!" ~ size.to!string ~ ".");
202         VecType ret;
203         if (op == "-")
204             foreach (i; 0 .. size) { ret[i] = -data[i]; }
205         if (op == "+")
206             foreach (i; 0 .. size) { ret[i] = data[i]; }
207         return ret;
208     }
209     
210     /// Invert vector
211     VecType opUnary(string op)() if (op == "~" && isFloatingPoint!T) {
212         // assert(this !is null, "\nOP::ERROR nullptr Vector!" ~ size.to!string ~ ".");
213         VecType ret;
214         foreach (i; 0 .. size) { ret[i] = 1.0 / data[i]; }
215         return ret;
216     }
217     
218     /// Ditto
219     dvec!N opUnary(string op)() if (op == "~" && !isFloatingPoint!T) {
220         // assert(this !is null, "\nOP::ERROR nullptr Vector!" ~ size.to!string ~ ".");
221         dvec!N ret;
222         foreach (i; 0 .. size) { ret[i] = 1.0 / data[i]; }
223         return ret;
224     }
225 
226     /// opOpAssign x [+, -, *, /, %]= y
227     VecType opOpAssign(string op, R)( in Vector!(R, N) b ) if ( isNumeric!R ) { 
228         // assert(/* this !is null && */ b !is null, "\nOP::ERROR nullptr Vector!" ~ size.to!string ~ ".");
229         foreach (i; 0 .. size) { mixin( "data[i] = data[i] " ~ op ~ " b.data[i];" ); }
230         return this;
231     }
232     
233     /// Ditto
234     VecType opOpAssign(string op, R)( in R b ) if ( isNumeric!R ) { 
235         // assert(this !is null, "\nOP::ERROR nullptr Vector!" ~ size.to!string ~ ".");
236         foreach (i; 0 .. size) { mixin( "data[i] = data[i] " ~ op ~ " b;" ); }
237         return this;
238     }
239 
240     /// opCast cast(x) y
241     R opCast(R)() const if (isVector!(R, N) && R.size == N){
242         R ret;
243         foreach (i; 0 ..  N) {
244             ret[i] = cast(R.dataType) data[i];
245         }
246         return ret;
247     }
248 
249     /// Ditto
250     R opCast(R)() const if (is(R == quat) && N == 4){
251         R ret;
252         foreach (i; 0 ..  N) {
253             ret[i] = cast(float) data[i];
254         }
255         return ret;
256     }
257 
258     /// Cast to matrix (column/row matrix)
259     R opCast(R)() const if (isMatrix!(R, N, 1)) {
260         R ret;
261         foreach (i; 0 ..  N) {
262             ret[i][0] = cast(R.dataType) data[i];
263         }
264         return ret;
265     }
266 
267     /// Ditto
268     R opCast(R)() const if (isMatrix!(R, 1, N)) {
269         R ret;
270         foreach (i; 0 ..  N) {
271             ret[0][i] = cast(R.dataType) data[i];
272         }
273         return ret;
274     }
275 
276     /// Cast to color
277     Color opCast(R)() const if (is(R == Color) && (N == 3 || N == 4) && isFloatingPoint!T){
278         Color ret;
279         foreach (i; 0 ..  N) {
280             ret[i] = cast(float) data[i];
281         }
282         if (N == 3) ret[3] = 1.0f;
283         return ret;
284     }
285 
286     /// Ditto
287     Color opCast(R)() const if (is(R == Color) && (N == 3 || N == 4) && !isFloatingPoint!T){
288         Color ret;
289         foreach (i; 0 ..  N) {
290             ret[i] = cast(float) (data[i] / 255.0f);
291         }
292         if (N == 3) ret[3] = 1.0f;
293         return ret;
294     }
295 
296     /// Cast to bool
297     bool opCast(T)() const if (is(T == bool)) {
298         return !lengthSquared.isClose(0, float.epsilon);
299     }
300 
301     /// Returns hash 
302     size_t toHash() const @safe nothrow {
303         return typeid(data).getHash(&data);
304     }
305     
306     static if (N == 2 || N == 3 || N == 4) {
307         static if (N == 2) private enum AS = "x y|w h|u v"; 
308         else
309         static if (N == 3) private enum AS = "x y z|w h d|u v t|r g b"; 
310         else
311         static if (N == 4) private enum AS = "x y z w|r g b a";
312         /// Mixes in swizzle
313         mixin accessByString!(T, N, "data", AS);
314     }
315 
316     /// Returns copy of vector
317     public VecType copyof() {
318         return VecType(data);
319     }
320 
321     /// Returns string representation of vector: `[1.00, 1.00,... , 1.00]`
322     public string toString() const {
323         import std.conv : to;
324         string s;
325         s ~= "[";
326         foreach (i; 0 .. size) {
327             s ~= isFloatingPoint!T ? format("%.2f", data[i]) : format("%d", data[i]);
328             if (i != size - 1) s ~= ", ";
329         }
330         s ~= "]";
331         return s;
332     }
333 
334     /// Returns pointer to data
335     public T* ptr() return {
336         return data.ptr;
337     }
338 
339     /* -------------------------------------------------------------------------- */
340     /*                         STATIC GETTERS AND SETTERS                         */
341     /* -------------------------------------------------------------------------- */
342     
343     /// Constructs predefined vector
344     static alias zero  = () => VecType(0);
345     /// Ditto
346     static alias one   = () => VecType(1);
347 
348     static if(isFloatingPoint!T) {
349         /// Ditto
350         static alias inf   = () => VecType(float.infinity);
351     }
352 
353     static if(N == 2) {
354         /// Ditto
355         static alias left  = () => VecType(-1, 0);
356         /// Ditto
357         static alias right = () => VecType(1, 0);
358         /// Ditto
359         static alias up    = () => VecType(0, -1);
360         /// Ditto
361         static alias down  = () => VecType(0, 1);
362     }
363 
364     static if(N == 3) {
365         static alias forward = () => VecType(0, 0, -1);
366         /// Ditto
367         static alias back    = () => VecType(0, 0, 1);
368         /// Ditto
369         static alias left    = () => VecType(-1, 0, 0);
370         /// Ditto
371         static alias right   = () => VecType(1, 0, 0);
372         /// Ditto
373         static alias up      = () => VecType(0, 1, 0);
374         /// Ditto
375         static alias down    = () => VecType(0, -1, 0);
376     }
377 
378     /* -------------------------------------------------------------------------- */
379     /*                                    MATH                                    */
380     /* -------------------------------------------------------------------------- */
381     
382     /// Returns vector length
383     public double length() const {
384         return sqrt(cast(double) lengthSquared);
385     }
386 
387     /// Returns squared vector length
388     public T lengthSquared() const {
389         T l = 0;
390         foreach (i; 0 .. size) { l += data[i] * data[i]; }
391         return l;
392     }
393 
394     /** 
395     Returns squared distance from vector to `b`
396     Params:
397       b = Vector to calculate distance to
398     Returns: Distance
399     */
400     public T distanceSquaredTo(VecType b) {
401         T dist = 0;
402         foreach (i; 0 .. size) { dist += (data[i] - b.data[i]) * (data[i] - b.data[i]); }
403         return dist;
404     }
405 
406     /** 
407     Calculates distance to vector `b`
408     Params:
409       b = Vector
410     Returns: Distance
411     */
412     public double distanceTo(VecType b) {
413         return sqrt(cast(double) distanceSquaredTo(b));
414     }
415 
416     /** 
417     Is vector approximately close to `v`
418     Params:
419       v = Vector to compare
420     Returns: 
421     */
422     public bool isClose(VecType v) {
423         bool eq = true;
424         foreach (i; 0 .. size) { eq = eq && data[i].isClose(v[i], float.epsilon); }
425         return eq;
426     }
427 
428     /// Normalises vector
429     public VecType normalized() {
430         T q = lengthSquared;
431         if (q == 0 || q.isClose(0, float.epsilon)) return this;
432         VecType ret;
433         double l = sqrt(cast(double) lengthSquared);
434         foreach (i; 0 .. size) { ret[i] = cast(T) (data[i] / l); }
435         return ret;
436     }
437     /// Ditto
438     alias normalised = normalized;
439     
440     /// Normalises vector in place
441     public VecType normalize() {
442         T q = lengthSquared;
443         if (q == 0 || q.isClose(0, float.epsilon)) return this;
444         double l = sqrt(cast(double) lengthSquared);
445         foreach (i; 0 .. size) { data[i] = cast(T) (data[i] / l); }
446         return this;
447     }
448     /// Ditto
449     alias normalise = normalize;
450 
451     /// Returns true if vector is normalised
452     public bool isNormalized() {
453         return lengthSquared.isClose(1, float.epsilon);
454     }
455     /// Ditto
456     alias isNormalised = isNormalized;
457 
458     /** 
459     Performs dot product
460     Params:
461       b = Vector
462     Returns: dot product
463     */
464     public double dot(VecType b) {
465         double d = 0;
466         foreach (i; 0 .. size) { d += cast(double) data[i] * cast(double) b.data[i]; }
467         return d;
468     }
469 
470     /// Signs current vector
471     public VecType sign() {
472         VecType ret;
473         foreach (i; 0 .. size) { ret[i] = data[i].sgn(); }
474         return ret;
475     }
476     
477     // Opearations that only make sense on floats
478     static if (isFloatingPoint!T) {
479         /// Floors vector values
480         public VecType floor() {
481             VecType ret;
482             foreach (i; 0 .. size) { ret[i] = data[i].floor(); }
483             return ret;
484         }
485 
486         /// Ceils vector values
487         public VecType ceil() {
488             VecType ret;
489             foreach (i; 0 .. size) { ret[i] = data[i].ceil(); }
490             return ret;
491         }
492 
493         /// Rounds vector values
494         public VecType round() {
495             VecType ret;
496             foreach (i; 0 .. size) { ret[i] = data[i].round(); }
497             return ret;
498         }
499 
500         /** 
501         Linear interpolates vector
502         Params:
503           to = Vector to interpolate to
504           weight = Interpolation weight in range [0.0, 1.0]
505         */
506         public VecType lerp(VecType to, double weight) {
507             VecType ret;
508             foreach (i; 0 .. size) { ret[i] = data[i] + (weight * (to.data[i] - data[i])); }
509             return ret;
510         }
511 
512     }
513 
514     /// Abs vector values
515     public VecType abs() {
516         VecType ret;
517         foreach (i; 0 .. size) { ret[i] = data[i].abs(); }
518         return ret;
519     }
520 
521     /** 
522     Clamps vector values to min
523     Params:
524       b = Minimal Vector
525     */
526     public VecType min(VecType b) {
527         VecType ret;
528         foreach (i; 0 .. size) { ret[i] = data[i].min(b.data[i]); }
529         return ret;
530     }
531 
532     /** 
533     Clamps vector values to max
534     Params:
535       b = Maximal Vector
536     */
537     public VecType max(VecType b) {
538         VecType ret;
539         foreach (i; 0 .. size) { ret[i] = data[i].max(b.data[i]); }
540         return ret;
541     }
542     
543     /** 
544     Clamps vector values
545     Params:
546       b = Minimal Vector
547       b = Maximal Vector
548     */
549     public VecType clamp(VecType p_min, VecType p_max) {
550         VecType ret;
551         foreach (i; 0 .. size) { ret[i] = data[i].clamp(p_min.data[i], p_max.data[i]); }
552         return ret;
553     }
554     
555     /** 
556     Clamps vector values in place
557     Params:
558       b = Minimal Vector
559       b = Maximal Vector
560     */
561     public VecType clamped(VecType p_min, VecType p_max) {
562         foreach (i; 0 .. size) { data[i] = data[i].clamp(p_min.data[i], p_max.data[i]); }
563         return this;
564     }
565     
566     /** 
567     Clamps vector values
568     Params:
569       b = Minimal value
570       b = Maximal value
571     */
572     public VecType clamp(T p_min, T p_max) {
573         VecType ret;
574         foreach (i; 0 .. size) { ret[i] = data[i].clamp(p_min, p_max); }
575         return ret;
576     }
577     
578     /** 
579     Clamps vector values in place
580     Params:
581       b = Minimal value
582       b = Maximal value
583     */
584     public VecType clamped(T p_min, T p_max) {
585         foreach (i; 0 .. size) { data[i] = data[i].clamp(p_min, p_max); }
586         return this;
587     }
588 
589     /** 
590     Snaps vector values
591     Params:
592       p_step = Vector to snap to
593     */
594     public VecType snap(VecType p_step) {
595         VecType ret;
596         foreach (i; 0 .. size) { 
597             ret[i] = data[i].snap(p_step[i]);
598         }
599         return ret;
600     }
601 
602     /** 
603     Snaps vector values in place
604     Params:
605       p_step = Vector to snap to
606     */
607     public VecType snapped(VecType p_step) {
608         foreach (i; 0 .. size) { 
609             data[i] = data[i].snap(p_step[i]);
610         }
611         return this;
612     }
613 
614     /** 
615     Snaps vector values
616     Params:
617       p_step = value to snap to
618     */
619     public VecType snap(T p_step) {
620         VecType ret;
621         foreach (i; 0 .. size) { 
622             ret[i] = data[i].snap(p_step);
623         }
624         return ret;
625     }
626 
627     /** 
628     Snaps vector values in place
629     Params:
630       p_step = value to snap to
631     */
632     public VecType snapped(T p_step) {
633         foreach (i; 0 .. size) { 
634             data[i] = data[i].snap(p_step);
635         }
636         return this;
637     }
638 
639     /// Calculates reflected vector to normal
640     public VecType reflect(VecType normal) {
641         double d = dot(normal);
642         VecType ret;
643         foreach (i; 0..size) {
644             ret[i] = cast(T) (data[i] - (2.0 * normal.data[i]) * d); 
645         }
646         return ret;
647     }
648 
649     /++
650     Calculates direction of refracted ray where this is incoming ray,
651     n is normal vector and r is refractive ratio
652     +/
653     public VecType refract(VecType n, double r) {
654         double dt = dot(n);
655         double d = 1.0 - r * r * (1.0 - dt * dt);
656         VecType ret = 0;
657         if (d >= 0) {
658             foreach (i; 0..size) {
659                 ret[i] = cast(T) (r * data[i] - (r * dt + d) * n.data[i]); 
660             }
661         }
662         return ret;
663     }
664 
665 
666     /// Moves vector toward target
667     public VecType moveTowards(VecType to, double maxDist) {
668         double[size] d;
669         double val = 0;
670         
671         foreach (i; 0..size) {
672             d[i] = to[i] - data[i];
673             val += d[i] * d[i];
674         }
675 
676         if (val == 0 || (maxDist >= 0 && val <= maxDist * maxDist)) return to;
677         
678         VecType ret;
679         double dist = sqrt(val);
680 
681         foreach (i; 0..size) {
682             ret[i] = cast(T) (data[i] + d[i] / dist * maxDist);
683         }
684 
685         return ret;
686     }
687         
688     /// Limits length of vector and returns resulting vector
689     public VecType limitLength(double p_min, double p_max) {
690         VecType ret = VecType(data);
691         double len = length;
692 
693         if (len > 0.0) {
694             len = sqrt(len);
695 
696             if (len < p_min) {
697                 double scale = p_min / len;
698                 foreach (i; 0..size) {
699                     ret[i] = cast(T) (ret[i] * scale);
700                 }
701             }
702 
703             if (len > p_max) {
704                 double scale = p_max / len;
705                 foreach (i; 0..size) {
706                     ret[i] = cast(T) (ret[i] * scale);
707                 }
708             }
709         }
710 
711         return ret;
712     }
713 
714     static if (N == 2) {
715         /// Calculates angle between this vector and v2 from [0, 0]
716         public double angle(VecType v2) {
717             static if (isFloatingPoint!T) {
718                 return atan2(v2.y - data[1], v2.x - data[0]);
719             } else {
720                 return atan2(cast(double) v2.y - data[1], cast(double) v2.x - data[0]);
721             }
722         }
723         
724         /// Calculates angle between line this->to and X ordinate
725         public double angleTo(VecType to) {
726             return acos(dot(to).clamp(-1.0, 1.0));
727         }
728 
729         /// Calculates cross product of this and b, assumes that Z = 0
730         public double cross(VecType b) {
731             return data[0] * b.y - data[1] * b.x;
732         }
733 
734         /// Rotates vector by angle
735         public VecType rotate(double angle) {
736             double sina = sin(angle);
737             double cosa = cos(angle);
738             return VecType( cast(T) (data[0] * cosa - data[1] * sina), cast(T) (data[0] * sina + data[1] * cosa) );
739         }
740     }
741 
742     static if (N == 3) {
743         /// Calculates angle between this vector and v2 from [0, 0]
744         public double angle(VecType v2) {
745             VecType cr = cross(v2);
746             double len = cr.length();
747             double dot = dot(v2);
748             return atan2(len, dot);
749         }
750 
751 
752         /// Calculates cross product of this and b
753         public VecType cross(VecType b) {
754             VecType cr = VecType();            
755             T ax = data[0],
756               ay = data[1],
757               az = data[2];
758             T bx = b.data[0],
759               by = b.data[1],
760               bz = b.data[2];
761             cr.data[0] = ay * bz - az * by;
762             cr.data[1] = az * bx - ax * bz;
763             cr.data[2] = ax * by - ay * bx;
764             return cr;
765         }
766         
767         /// Returns vector perpendicular to this
768         public VecType perpendicular() {
769             double p_min = cast(double) data[0].abs;
770             VecType cardinal = VecType(1, 0, 0);
771             
772             if (data[1].abs < p_min) {
773                 p_min = cast(double) data[1].abs;
774                 cardinal = VecType(0, 1, 0);
775             }
776 
777             if (data[2].abs < p_min) {
778                 cardinal = VecType(0, 0, 1);
779             }
780 
781             return VecType(
782                 data[1] * cardinal.z - data[2] * cardinal.y,
783                 data[2] * cardinal.x - data[0] * cardinal.z,
784                 data[0] * cardinal.y - data[1] * cardinal.x
785             );
786         }
787 
788         /// Projects vector3 from screen space into object space
789         public VecType unproject(mat4 proj, mat4 view) {
790             mat4 vpnorm = (view * proj);
791             mat4 vpinv = ~vpnorm;
792             quat qt = [data[0], data[1], data[2], 1.0f];
793             quat qtrs = qt * vpinv;
794 
795             return VecType( cast(T) qtrs.x, cast(T) qtrs.y, cast(T) qtrs.z ) / cast(T) qtrs.w;
796         }
797 
798         // todo project(mat4 proj, mat4 view) {}
799     }
800 
801     static if (N == 4) {
802         /// Calculates cross product of this and b, 
803         /// W is multiplication of a.w and b.w
804         public VecType cross(VecType b) {
805             VecType cr = VecType();            
806             T ax = data[0],
807               ay = data[1],
808               az = data[2],
809               aw = data[3];
810             T bx = b.data[0],
811               by = b.data[1],
812               bz = b.data[2],
813               bw = b.data[3];
814             cr.data[0] = ay * bz - az * by;
815             cr.data[1] = az * bx - ax * bz;
816             cr.data[2] = ax * by - ay * bx;
817             cr.data[3] = aw * bw;
818             return cr;
819 
820         }
821     }
822 }
823 
824 /// Is V a vector with any size and any type
825 template isVector(V) {
826     enum isVector = is(Unqual!V == Vector!(U, n), size_t n, U);
827 }
828 
829 /// Is V a vector with size n and any type
830 template isVector(V, size_t n) {
831     static if(is(Unqual!V == Vector!(U, n2), size_t n2, U)) {
832         enum isVector = n == n2;
833     } else {
834         enum isVector = false;
835     }
836 }
837 
838 /// Is V a vector with any size and type T
839 template isVector(V, T) {
840     static if(is(Unqual!V == Vector!(U, n), size_t n, U)) {
841         enum isVector = is(T == U);
842     } else {
843         enum isVector = false;
844     }
845 }
846 
847 /// Is V a vector with size N and type T
848 template isVector(V, T, size_t N) {
849     enum isVector = is(Unqual!V == Vector!(T, N));
850 }
851