package sample.svg.contributed;
import static util.SVGUtil.*;

public class Nagae2 {
    private static final double  CENTER_X = 148.5;
    private static final double  CENTER_Y = 105.0;
    private static final int  SCALE = 4;

    private static final int  COLOR_DELTA = 17;
    private static final int  COLOR_MAX = 0xFF;
    private static final int  COLOR_MIN = 0x00;

    private static final double  P = 10.0;
    private static final double  R = 28.0;
    private static final double  B = 2.6666667;

    private static double dfXdT(double t, double x, double y, double z) {
        return (P)*((y) - (x)) + (t);
    }

    private static double dfYdT(double t, double x, double y, double z) {
        return (R)*(x) - (x)*(z) - (y);
    }

    private static double dfZdT(double t, double x, double y, double z) {
        return (x)*(y) - (B)*(z) - (t);
    } 

    private static double mySin = Math.sin(-PI/4.0);
    private static double myCos = Math.cos(-PI/4.0);

    private static int cg = COLOR_MAX;
    private static int cb = COLOR_MIN;
    private static int mode = 0;

    // 座標変換と出力
    private static void output(
                               double xold, double yold, double zold,
                               double xnew, double ynew, double znew
                               ) {
        double xo = xold*myCos - yold*mySin;
        double yo = xold*mySin + yold*myCos;
        double xn = xnew*myCos - ynew*mySin;
        double yn = xnew*mySin + ynew*myCos;

        xo *= SCALE;
        yo *= SCALE;
        xn *= SCALE;
        yn *= SCALE;

        xo += CENTER_X;
        yo += CENTER_Y;
        xn += CENTER_X;
        yn += CENTER_Y;

        switch (mode) {
        case 0: cb += COLOR_DELTA; if(cb == COLOR_MAX) mode++; break;
        case 1: cg -= COLOR_DELTA; if(cg == COLOR_MIN) mode++; break;
        case 2: cb -= COLOR_DELTA; if(cb == COLOR_MIN) mode++; break;
        case 3: cg += COLOR_DELTA; if(cg == COLOR_MAX) mode=0; break;
        }

        stroke(rgb255(0, cg, cb));
        line(xo, yo, xn, yn);
    }

    // ルンゲクッタ法
    private static void rungekutta(
                                   double t0, double tx, int cnt,
                                   double x0, double y0, double z0
                                   ) {
        int i;
        double dt = (tx - t0) / cnt;

        double       told = t0;
        double xnew, xold = x0;
        double ynew, yold = y0;
        double znew, zold = z0;

        for (i = 1; i <= cnt; i++) {
            double kx1 = dt*dfXdT(told, xold, yold, zold);
            double ky1 = dt*dfYdT(told, xold, yold, zold);
            double kz1 = dt*dfZdT(told, xold, yold, zold);

            double kx2 = dt*dfXdT(told + dt/2, xold + kx1/2, yold + ky1/2, zold + kz1/2);
            double ky2 = dt*dfYdT(told + dt/2, xold + kx1/2, yold + ky1/2, zold + kz1/2);
            double kz2 = dt*dfZdT(told + dt/2, xold + kx1/2, yold + ky1/2, zold + kz1/2);

            double kx3 = dt*dfXdT(told + dt/2, xold + kx2/2, yold + ky2/2, zold + kz2/2);
            double ky3 = dt*dfYdT(told + dt/2, xold + kx2/2, yold + ky2/2, zold + kz2/2);
            double kz3 = dt*dfZdT(told + dt/2, xold + kx2/2, yold + ky2/2, zold + kz2/2);

            double kx4 = dt*dfXdT(told + dt/2, xold + kx3, yold + ky3, zold + kz3);
            double ky4 = dt*dfYdT(told + dt/2, xold + kx3, yold + ky3, zold + kz3);
            double kz4 = dt*dfZdT(told + dt/2, xold + kx3, yold + ky3, zold + kz3);

            xnew = xold + (kx1 + kx2*2 + kx3*2 + kx4)/6;
            ynew = yold + (ky1 + ky2*2 + ky3*2 + ky4)/6;
            znew = zold + (kz1 + kz2*2 + kz3*2 + kz4)/6;

            output(xold, yold, zold, xnew, ynew, znew);

            xold = xnew;
            yold = ynew;
            zold = znew;

            told = t0 + dt*i;
        }
    }

    public static void main(String[] args) {
        start();
        rulers();
        
        strokeWeight(0.5);

        translate(25, 25);
        scale(0.9, 0.9);

        rungekutta(0, 150, 30000, 1, 1, 1);

        finish();
        return;
    }
}