1 module ggplotd.geom;
2 
3 import std.range : front, popFront, empty;
4 
5 import cairo = cairo.cairo;
6 
7 import ggplotd.bounds;
8 import ggplotd.aes;
9 import ggplotd.colour : ColourID, ColourMap;
10 
11 version (unittest)
12 {
13     import dunit.toolkit;
14 }
15 
16 ///
17 struct Geom
18 {
19     this(T)( in T tup ) //if (is(T==Tuple))
20     {
21         mask = tup.mask;
22     }
23 
24     alias drawFunction = cairo.Context delegate(cairo.Context context, 
25         ColourMap colourMap);
26     drawFunction draw; ///
27     ColourID[] colours; ///
28     AdaptiveBounds bounds; ///
29 
30     bool mask = true; /// Whether to mask/prevent drawing outside plotting area
31 
32     import std.typecons : Tuple;
33 
34     Tuple!(double, string)[] xTickLabels; ///
35     Tuple!(double, string)[] yTickLabels; ///
36 }
37 
38 ///
39 auto geomPoint(AES)(AES aes)
40 {
41     alias CoordX = typeof(NumericLabel!(typeof(AES.x))(AES.x));
42     alias CoordY = typeof(NumericLabel!(typeof(AES.y))(AES.y));
43     alias CoordType = typeof(merge(aes, Aes!(CoordX, "x", CoordY,
44         "y")(CoordX(AES.x), CoordY(AES.y))));
45 
46     struct GeomRange(T)
47     {
48         this(T aes)
49         {
50             _aes = merge(aes, Aes!(CoordX, "x", CoordY, "y")(CoordX(aes.x), CoordY(aes.y)));
51         }
52 
53         @property auto front()
54         {
55             immutable tup = _aes.front;
56             auto f = delegate(cairo.Context context, ColourMap colourMap ) 
57             {
58                 auto devP = context.userToDevice(cairo.Point!double(tup.x[0], tup.y[0]));
59                 context.save();
60                 context.identityMatrix;
61                 context.rectangle(devP.x - 0.5 * tup.size, devP.y - 0.5 * tup.size, tup.size, tup.size);
62                 context.restore();
63 
64                 auto col = colourMap(ColourID(tup.colour));
65                 import cairo.cairo : RGBA;
66 
67                 context.identityMatrix();
68 
69                 context.setSourceRGBA(RGBA(col.red, col.green, col.blue, tup.alpha));
70                 context.fill();
71 
72                 return context;
73             };
74 
75             AdaptiveBounds bounds;
76             bounds.adapt(Point(tup.x[0], tup.y[0]));
77             auto geom = Geom( tup );
78             geom.draw = f;
79             geom.colours ~= ColourID(tup.colour);
80             geom.bounds = bounds;
81             return geom;
82         }
83 
84         void popFront()
85         {
86             _aes.popFront();
87         }
88 
89         @property bool empty()
90         {
91             return _aes.empty;
92         }
93 
94     private:
95         CoordType _aes;
96     }
97 
98     return GeomRange!AES(aes);
99 }
100 
101 ///
102 unittest
103 {
104     auto aes = Aes!(double[], "x", double[], "y")([1.0], [2.0]);
105     auto gl = geomPoint(aes);
106     assertEqual(gl.front.colours[0][1], "black");
107     gl.popFront;
108     assert(gl.empty);
109 }
110 
111 ///
112 auto geomLine(AES)(AES aes)
113 {
114     import std.algorithm : map;
115     import std.range : array, zip;
116 
117     struct GeomRange(T)
118     {
119         this(T aes)
120         {
121             groupedAes = aes.group;
122         }
123 
124         @property auto front()
125         {
126             auto xs = NumericLabel!(typeof(groupedAes.front.front.x)[])(
127                 groupedAes.front.map!((t) => t.x).array);
128             auto ys = NumericLabel!(typeof(groupedAes.front.front.y)[])(
129                 groupedAes.front.map!((t) => t.y).array);
130             auto coords = zip(xs, ys);
131 
132             immutable flags = groupedAes.front.front;
133             auto f = delegate(cairo.Context context, ColourMap colourMap ) {
134                 auto fr = coords.front;
135                 context.moveTo(fr[0][0], fr[1][0]);
136                 coords.popFront;
137                 foreach (tup; coords)
138                 {
139                     context.lineTo(tup[0][0], tup[1][0]);
140                 }
141 
142                 auto col = colourMap(ColourID(flags.colour));
143                 import cairo.cairo : RGBA;
144 
145                 context.identityMatrix();
146                 if (flags.fill>0)
147                 {
148                     context.setSourceRGBA(RGBA(col.red, col.green, col.blue, flags.fill));
149                     context.fillPreserve();
150                 }
151                 context.setSourceRGBA(RGBA(col.red, col.green, col.blue, flags.alpha));
152                 context.stroke();
153 
154                 return context;
155             };
156 
157             AdaptiveBounds bounds;
158             coords = zip(xs, ys);
159             auto geom = Geom(groupedAes.front.front);
160             foreach (tup; coords)
161             {
162                 bounds.adapt(Point(tup[0][0], tup[1][0]));
163                 if (!xs.numeric)
164                     geom.xTickLabels ~= tup[0];
165                 if (!ys.numeric)
166                     geom.yTickLabels ~= tup[1];
167             }
168             geom.draw = f;
169             geom.colours ~= ColourID(groupedAes.front.front.colour);
170             geom.bounds = bounds;
171             return geom;
172         }
173 
174         void popFront()
175         {
176             groupedAes.popFront;
177         }
178 
179         @property bool empty()
180         {
181             return groupedAes.empty;
182         }
183 
184     private:
185         typeof(group(T.init)) groupedAes;
186     }
187 
188     return GeomRange!AES(aes);
189 }
190 
191 ///
192 unittest
193 {
194     auto aes = Aes!(double[], "x", double[], "y", string[], "colour")([1.0,
195         2.0, 1.1, 3.0], [3.0, 1.5, 1.1, 1.8], ["a", "b", "a", "b"]);
196 
197     auto gl = geomLine(aes);
198 
199     import std.range : empty;
200 
201     assert(gl.front.xTickLabels.empty);
202     assert(gl.front.yTickLabels.empty);
203 
204     assertEqual(gl.front.colours[0][1], "a");
205     assertEqual(gl.front.bounds.min_x, 1.0);
206     assertEqual(gl.front.bounds.max_x, 1.1);
207     gl.popFront;
208     assertEqual(gl.front.colours[0][1], "b");
209     assertEqual(gl.front.bounds.max_x, 3.0);
210     gl.popFront;
211     assert(gl.empty);
212 }
213 
214 unittest
215 {
216     auto aes = Aes!(string[], "x", string[], "y", string[], "colour")(["a",
217         "b", "c", "b"], ["a", "b", "b", "a"], ["b", "b", "b", "b"]);
218 
219     auto gl = geomLine(aes);
220     assertEqual(gl.front.xTickLabels.length, 4);
221     assertEqual(gl.front.yTickLabels.length, 4);
222 }
223 
224 unittest
225 {
226     auto aes = Aes!(string[], "x", string[], "y", string[], "colour")(["a",
227         "b", "c", "b"], ["a", "b", "b", "a"], ["b", "b", "b", "b"]);
228 
229     auto gl = geomLine(aes);
230     auto aes2 = Aes!(string[], "x", string[], "y", double[], "colour")(["a",
231         "b", "c", "b"], ["a", "b", "b", "a"], [0, 1, 0, 0.1]);
232 
233     auto gl2 = geomLine(aes2);
234 
235     import std.range : chain, walkLength;
236 
237     assertEqual(gl.chain(gl2).walkLength, 4);
238 }
239 
240 // Bin a range of data
241 private auto bin(R)(R xs, size_t noBins = 10)
242 {
243     struct Bin
244     {
245         double[] range;
246         size_t count;
247     }
248 
249     import std.typecons : Tuple;
250     import std.algorithm : group;
251 
252     struct BinRange(Range)
253     {
254         this(Range xs, size_t noBins)
255         {
256             import std.math : floor;
257             import std.algorithm : min, max, reduce, sort, map;
258             import std.array : array;
259             import std.range : walkLength;
260 
261             assert(xs.walkLength > 0);
262 
263             // Find the min and max values
264             auto minmax = xs.reduce!((a, b) => min(a, b), (a, b) => max(a, b));
265             _width = (minmax[1] - minmax[0]) / (noBins - 1);
266             _noBins = noBins;
267             // If min == max we need to set a custom width
268             if (_width == 0)
269                 _width = 0.1;
270             _min = minmax[0] - 0.5 * _width;
271 
272             // Count the number of data points that fall in a
273             // bin. This is done by scaling them into whole numbers
274             counts = xs.map!((a) => floor((a - _min) / _width)).array.sort().array.group();
275 
276             // Initialize our bins
277             if (counts.front[0] == _binID)
278             {
279                 _cnt = counts.front[1];
280                 counts.popFront;
281             }
282         }
283 
284         /// Return a bin describing the range and number of data points (count) that fall within that range.
285         @property auto front()
286         {
287             return Bin([_min, _min + _width], _cnt);
288         }
289 
290         void popFront()
291         {
292             _min += _width;
293             _cnt = 0;
294             ++_binID;
295             if (!counts.empty && counts.front[0] == _binID)
296             {
297                 _cnt = counts.front[1];
298                 counts.popFront;
299             }
300         }
301 
302         @property bool empty()
303         {
304             return _binID >= _noBins;
305         }
306 
307     private:
308         double _min;
309         double _width;
310         size_t _noBins;
311         size_t _binID = 0;
312         typeof(group(Range.init)) counts;
313         size_t _cnt = 0;
314     }
315 
316     return BinRange!R(xs, noBins);
317 }
318 
319 unittest
320 {
321     import std.array : array;
322     import std.range : back, walkLength;
323 
324     auto binR = bin!(double[])([0.5, 0.01, 0.0, 0.9, 1.0, 0.99], 11);
325     assertEqual(binR.walkLength, 11);
326     assertEqual(binR.front.range, [-0.05, 0.05]);
327     assertEqual(binR.front.count, 2);
328     assertLessThan(binR.array.back.range[0], 1);
329     assertGreaterThan(binR.array.back.range[1], 1);
330     assertEqual(binR.array.back.count, 2);
331 
332     binR = bin!(double[])([0.01], 11);
333     assertEqual(binR.walkLength, 11);
334     assertEqual(binR.front.count, 1);
335 
336     binR = bin!(double[])([-0.01, 0, 0, 0, 0.01], 11);
337     assertEqual(binR.walkLength, 11);
338     assertLessThan(binR.front.range[0], -0.01);
339     assertGreaterThan(binR.front.range[1], -0.01);
340     assertEqual(binR.front.count, 1);
341     assertLessThan(binR.array.back.range[0], 0.01);
342     assertGreaterThan(binR.array.back.range[1], 0.01);
343     assertEqual(binR.array.back.count, 1);
344     assertEqual(binR.array[5].count, 3);
345     assertLessThan(binR.array[5].range[0], 0.0);
346     assertGreaterThan(binR.array[5].range[1], 0.0);
347 }
348 
349 
350 /// Draw histograms based on the x coordinates of the data (aes)
351 auto geomHist(AES)(AES aes)
352 {
353     import std.algorithm : map;
354     import std.array : Appender, array;
355     import std.range : repeat;
356     import std.typecons : Tuple;
357 
358     // New appender to hold lines for drawing histogram
359     auto appender = Appender!(Geom[])([]);
360 
361     foreach (grouped; group(aes)) // Split data by colour/id
362     {
363         auto bins = grouped.map!((t) => t.x) // Extract the x coordinates
364             .array.bin(11); // Bin the data
365 
366         foreach (bin; bins)
367         {
368             // Specifying the boxes for the histogram. The merge is used to keep the colour etc. information
369             // contained in the original aes passed to geomHist.
370             appender.put(
371                 geomLine( [
372                     grouped.front.merge(Tuple!(double, "x", double, "y" )( 
373                             bin.range[0], 0.0 )),
374                     grouped.front.merge(Tuple!(double, "x", double, "y" )( 
375                             bin.range[0], bin.count )),
376                     grouped.front.merge(Tuple!(double, "x", double, "y" )( 
377                             bin.range[1], bin.count )),
378                     grouped.front.merge(Tuple!(double, "x", double, "y" )( 
379                             bin.range[1], 0.0 )),
380                 ] )
381             );
382         }
383     }
384 
385     // Return the different lines 
386     return appender.data;
387 }
388 
389 /// Draw axis, first and last location are start/finish
390 /// others are ticks (perpendicular)
391 auto geomAxis(AES)(AES aes, double tickLength, string label)
392 {
393     import std.algorithm : find;
394     import std.array : array;
395     import std.range : chain, empty, repeat;
396     import std.math : sqrt, pow;
397 
398     double[] xs;
399     double[] ys;
400 
401     double[] lxs;
402     double[] lys;
403     double[] langles;
404     string[] lbls;
405 
406     auto colour = aes.front.colour;
407     auto toDir = aes.find!("a.x != b.x || a.y != b.y")(aes.front).front; 
408     auto direction = [toDir.x - aes.front.x, toDir.y - aes.front.y];
409     auto dirLength = sqrt(pow(direction[0], 2) + pow(direction[1], 2));
410     direction[0] *= tickLength / dirLength;
411     direction[1] *= tickLength / dirLength;
412  
413     while (!aes.empty)
414     {
415         auto tick = aes.front;
416         xs ~= tick.x;
417         ys ~= tick.y;
418 
419         aes.popFront;
420 
421         // Draw ticks perpendicular to main axis;
422         if (xs.length > 1 && !aes.empty)
423         {
424             xs ~= [tick.x + direction[1], tick.x];
425             ys ~= [tick.y + direction[0], tick.y];
426 
427             lxs ~= tick.x - 1.5*direction[1];
428             lys ~= tick.y - 1.5*direction[0];
429             lbls ~= tick.label;
430             langles ~= tick.angle;
431         }
432     }
433 
434     // Main label
435     auto xm = xs[0] + 0.5*(xs[$-1]-xs[0]) - 4.0*direction[1];
436     auto ym = ys[0] + 0.5*(ys[$-1]-ys[0]) - 4.0*direction[0];
437     auto aesM = Aes!(double[], "x", double[], "y", string[], "label", 
438         double[], "angle", bool[], "mask")( [xm], [ym], [label], 
439             langles, [false]);
440 
441     return geomLine(Aes!(typeof(xs), "x", typeof(ys), "y", bool[], "mask")(
442         xs, ys, false.repeat(xs.length).array)).chain(
443         geomLabel(Aes!(double[], "x", double[], "y", string[], "label",
444         double[], "angle", bool[], "mask")(lxs, lys, lbls, langles, 
445             false.repeat(lxs.length).array)))
446             .chain( geomLabel(aesM) );
447 }
448 
449 /// Draw Label at given x and y position
450 auto geomLabel(AES)(AES aes)
451 {
452     alias CoordX = typeof(NumericLabel!(typeof(AES.x))(AES.x));
453     alias CoordY = typeof(NumericLabel!(typeof(AES.y))(AES.y));
454     alias CoordType = typeof(merge(aes, Aes!(CoordX, "x", CoordY,
455         "y")(CoordX(AES.x), CoordY(AES.y))));
456 
457     struct GeomRange(T)
458     {
459         size_t size = 6;
460         this(T aes)
461         {
462             _aes = merge(aes, Aes!(CoordX, "x", CoordY, "y")(CoordX(aes.x), CoordY(aes.y)));
463         }
464 
465         @property auto front()
466         {
467             immutable tup = _aes.front;
468             auto f = delegate(cairo.Context context, ColourMap colourMap) {
469                 context.setFontSize(14.0);
470                 context.moveTo(tup.x[0], tup.y[0]);
471                 context.save();
472                 context.identityMatrix;
473                 context.rotate(tup.angle);
474                 auto extents = context.textExtents(tup.label);
475                 auto textSize = cairo.Point!double(0.5 * extents.width, 0.5 * extents.height);
476                 context.relMoveTo(-textSize.x, textSize.y);
477 
478                 auto col = colourMap(ColourID(tup.colour));
479                 import cairo.cairo : RGBA;
480 
481                 context.setSourceRGBA(RGBA(col.red, col.green, col.blue, tup.alpha));
482  
483                 context.showText(tup.label);
484                 context.restore();
485                 return context;
486             };
487 
488             AdaptiveBounds bounds;
489             bounds.adapt(Point(tup.x[0], tup.y[0]));
490 
491             auto geom = Geom( tup );
492             geom.draw = f;
493             geom.colours ~= ColourID(tup.colour);
494             geom.bounds = bounds;
495  
496             return geom;
497         }
498 
499         void popFront()
500         {
501             _aes.popFront();
502         }
503 
504         @property bool empty()
505         {
506             return _aes.empty;
507         }
508 
509     private:
510         CoordType _aes;
511     }
512 
513     return GeomRange!AES(aes);
514 }
515 
516 unittest
517 {
518     auto aes = Aes!(string[], "x", string[], "y", string[], "label")(["a", "b",
519         "c", "b"], ["a", "b", "b", "a"], ["b", "b", "b", "b"]);
520 
521     auto gl = geomLabel(aes);
522     import std.range : walkLength;
523 
524     assertEqual(gl.walkLength, 4);
525 }
526 
527 // geomBox
528 /// Return the limits indicated with different alphas
529 private auto limits( RANGE )( RANGE range, double[] alphas )
530 {
531     import std.algorithm : sort, map, min, max;
532     import std.math : floor;
533     import std.conv : to;
534     auto sorted = range.sort();
535     return alphas.map!( (a) { 
536         auto id = min( sorted.length-2,
537             max(0,floor( a*(sorted.length+1) ).to!int-1 ) );
538         if (a<=0.5)
539             return sorted[id];
540         else
541             return sorted[id+1];
542     });
543 }
544 
545 unittest
546 {
547     import std.range : array, front;
548     assertEqual( [1,2,3,4,5].limits( [0.01, 0.5, 0.99] ).array, 
549             [1,3,5] );
550 
551     assertEqual( [1,2,3,4].limits( [0.41] ).front, 2 );
552     assertEqual( [1,2,3,4].limits( [0.39] ).front, 1 );
553     assertEqual( [1,2,3,4].limits( [0.61] ).front, 4 );
554     assertEqual( [1,2,3,4].limits( [0.59] ).front, 3 );
555 }
556 
557 /// Draw a boxplot. The "x" data is used. If labels are given then the data is grouped by the label
558 auto geomBox(AES)(AES aes)
559 {
560     import std.algorithm : map;
561     import std.array : array;
562     import std.range : Appender;
563 
564     Appender!(Geom[]) result;
565     auto labels = NumericLabel!(string[])( 
566         aes.map!("a.label.to!string").array );
567     auto myAes = aes.merge( Aes!(typeof(labels), "label")( labels ) );
568 
569     double delta = 0.2;
570     Tuple!(double, string)[] xTickLabels;
571 
572     foreach( grouped; myAes.group() )
573     {
574         auto lims = grouped.map!("a.x")
575             .array.limits( [0.1,0.25,0.5,0.75,0.9] ).array;
576         auto x = grouped.front.label[0];
577         xTickLabels ~= grouped.front.label;
578         result.put(
579             geomLine( [
580                 grouped.front.merge(Tuple!(double, "x", double, "y" )( 
581                     x, lims[0] )),
582                 grouped.front.merge(Tuple!(double, "x", double, "y" )( 
583                     x, lims[1] )),
584                 grouped.front.merge(Tuple!(double, "x", double, "y" )( 
585                     x+delta, lims[1] )),
586                 grouped.front.merge(Tuple!(double, "x", double, "y" )( 
587                     x+delta, lims[2] )),
588                 grouped.front.merge(Tuple!(double, "x", double, "y" )( 
589                     x-delta, lims[2] )),
590                 grouped.front.merge(Tuple!(double, "x", double, "y" )( 
591                     x-delta, lims[3] )),
592                 grouped.front.merge(Tuple!(double, "x", double, "y" )( 
593                     x, lims[3] )),
594                 grouped.front.merge(Tuple!(double, "x", double, "y" )( 
595                     x, lims[4] )),
596 
597                 grouped.front.merge(Tuple!(double, "x", double, "y" )( 
598                     x, lims[3] )),
599                 grouped.front.merge(Tuple!(double, "x", double, "y" )( 
600                     x+delta, lims[3] )),
601                 grouped.front.merge(Tuple!(double, "x", double, "y" )( 
602                     x+delta, lims[2] )),
603                 grouped.front.merge(Tuple!(double, "x", double, "y" )( 
604                     x-delta, lims[2] )),
605                 grouped.front.merge(Tuple!(double, "x", double, "y" )( 
606                     x-delta, lims[1] )),
607                 grouped.front.merge(Tuple!(double, "x", double, "y" )( 
608                     x, lims[1] ))
609              ] )
610         );
611     }
612 
613     foreach( ref g; result.data )
614     {
615         g.xTickLabels = xTickLabels;
616         g.bounds.min_x = xTickLabels.front[0] - 0.5;
617         g.bounds.max_x = xTickLabels[$-1][0] + 0.5;
618     }
619 
620     return result.data;
621 }
622 
623 ///
624 auto geomPolygon(AES)(AES aes)
625 {
626     import std.array : array;
627     import std.algorithm : map, swap;
628     import ggplotd.geometry;
629     // Turn into vertices.
630     auto vertices = aes.map!( (t) => Vertex3D( t.x, t.y, t.colour ) );
631 
632     // Find lowest, highest
633     auto triangle = vertices.array;
634     if (triangle[1].z < triangle[0].z)
635         swap( triangle[1], triangle[0] );
636     if (triangle[2].z < triangle[0].z)
637         swap( triangle[2], triangle[0] );
638     if (triangle[1].z > triangle[2].z)
639         swap( triangle[1], triangle[2] );
640 
641     if (triangle.length > 3)
642         foreach( v; triangle[3..$] )
643         {
644             if (v.z < triangle[0].z)
645                 swap( triangle[0], v );
646             else if ( v.z > triangle[2].z )
647                 swap( triangle[2], v );
648         }
649     auto gV = gradientVector( triangle[0..3] );
650 
651     immutable flags = aes.front;
652 
653     auto geom = Geom( flags );
654 
655     // Define drawFunction
656     auto f = delegate(cairo.Context context, ColourMap colourMap ) 
657     {
658         auto gradient = new cairo.LinearGradient( gV[0].x, gV[0].y, 
659             gV[1].x, gV[1].y );
660 
661         auto col0 = colourMap(ColourID(gV[0].z));
662         auto col1 = colourMap(ColourID(gV[1].z));
663         import cairo.cairo : RGBA;
664         gradient.addColorStopRGBA( 0,
665             RGBA(col0.red, col0.green, col0.blue, flags.alpha));
666         gradient.addColorStopRGBA( 1,
667             RGBA(col1.red, col1.green, col1.blue, flags.alpha));
668         context.moveTo( vertices.front.x, vertices.front.y );
669         vertices.popFront;
670         foreach( v; vertices )
671             context.lineTo( v.x, v.y );
672         context.setSource( gradient );
673         context.fill;
674         return context;
675     };
676 
677     geom.draw = f;
678 
679     geom.colours = aes.map!((t) => ColourID(t.colour)).array;
680 
681     return [geom];
682 }