My attempt at Roger Bagulas L-system

/*
read all for further reference:

https://plus.google.com/101799841244447089430/posts/QboKwsJRJYb
*/

void rb_dual_spiral(
    ct::plot2d& plot, // our visual plotter
    unsigned int irecur,
    unsigned int irecur_max,
    ct_complex p0,
    ct_complex p1
) {
    //F+F+FF+FF+FFFFF+FF+FF+F+F

    if (irecur >= irecur_max) return;

    ct_complex dif = p1 - p0;
    ct_float unit = (std::abs(dif) / 5) / std::abs(dif);


    dif = p1 - p0;
    ct_complex ml0 = {
        p0.real() + dif.real() * unit * 1,
        p0.imag() + dif.imag() * unit * 1
    };

    ct_complex ml1 = {
        p0.real() + dif.real() * unit * 2,
        p0.imag() + dif.imag() * unit * 2
    };

    ct_complex ml2 = {
        p0.real() + dif.real() * unit * 3,
        p0.imag() + dif.imag() * unit * 3
    };

    ct_complex ml3 = {
        p0.real() + dif.real() * unit * 4,
        p0.imag() + dif.imag() * unit * 4
    };


    dif = p1 - p0;
    ct_complex lp0 = {
        p0.real() + -dif.imag() * unit * 2,
        p0.imag() + dif.real() * unit * 2
    };


    dif = lp0 - p0;
    ct_complex lp1 = {
        lp0.real() + dif.imag() * 1,
        lp0.imag() + -dif.real() * 1
    };


    dif = lp1 - lp0;
    ct_complex lp2 = {
        lp1.real() + dif.imag() * .5,
        lp1.imag() + -dif.real() * .5
    };


    dif = lp2 - lp1;
    ct_complex lp3 = {
        lp2.real() + dif.imag() * 1,
        lp2.imag() + -dif.real() * 1
    };


    ////
    dif = p1 - p0;
    ct_complex rp0 = {
        p1.real() + -dif.imag() * unit * 2,
        p1.imag() + dif.real() * unit * 2
    };


    dif = rp0 - p1;
    ct_complex rp1 = {
        rp0.real() + -dif.imag() * 1,
        rp0.imag() + dif.real() * 1
    };


    dif = rp1 - rp0;
    ct_complex rp2 = {
        rp1.real() + -dif.imag() * .5,
        rp1.imag() + dif.real() * .5
    };


    dif = rp2 - rp1;
    ct_complex rp3 = {
        rp2.real() + -dif.imag() * 1,
        rp2.imag() + dif.real() * 1
    };

    plot.line(p0, ml0, CT_RGBF(0, 0, 0));
    plot.line(ml0, ml1, CT_RGBF(0, 0, 0));
    plot.line(ml1, ml2, CT_RGBF(0, 0, 0));
    plot.line(ml2, ml3, CT_RGBF(0, 0, 0));
    plot.line(ml3, p1, CT_RGBF(0, 0, 0));

    plot.line(p0, lp0, CT_RGBF(1, 0, 0));
    plot.line(lp0, lp1, CT_RGBF(1, 0, 0));
    plot.line(lp1, lp2, CT_RGBF(1, 0, 0));
    plot.line(lp2, lp3, CT_RGBF(1, 0, 0));

    plot.line(p1, rp0, CT_RGBF(1, 0, 1));
    plot.line(rp0, rp1, CT_RGBF(1, 0, 1));
    plot.line(rp1, rp2, CT_RGBF(1, 0, 1));
    plot.line(rp2, rp3, CT_RGBF(1, 0, 1));

    rb_dual_spiral(plot, irecur + 1, irecur_max, p0, ml0);
    rb_dual_spiral(plot, irecur + 1, irecur_max, ml0, ml1);
    rb_dual_spiral(plot, irecur + 1, irecur_max, ml1, ml2);
    rb_dual_spiral(plot, irecur + 1, irecur_max, ml2, ml3);
    rb_dual_spiral(plot, irecur + 1, irecur_max, ml3, p1);

    rb_dual_spiral(plot, irecur + 1, irecur_max, lp0, p0);
    rb_dual_spiral(plot, irecur + 1, irecur_max, lp1, lp0);
    rb_dual_spiral(plot, irecur + 1, irecur_max, lp2, lp1);
    rb_dual_spiral(plot, irecur + 1, irecur_max, lp3, lp2);

    rb_dual_spiral(plot, irecur + 1, irecur_max, p1, rp0);
    rb_dual_spiral(plot, irecur + 1, irecur_max, rp0, rp1);
    rb_dual_spiral(plot, irecur + 1, irecur_max, rp1, rp2);
    rb_dual_spiral(plot, irecur + 1, irecur_max, rp2, rp3);

}