// FireworksGraphics32.cpp : This file contains the 'main' function. Program execution begins and ends there.
//
// reference for projectile calculations:
// https://scipython.com/book2/chapter-8-scipy/examples/a-projectile-with-air-resistance/

#include <iostream>
#include <fstream>

// number of firework particles to model
#define dimp 32
// dimension of the display array (assume square)
#define dim 16
// number of frames to calculate (assume about 20 frames per second)
#define numf 30

    // Drag coefficient, projectile radius(m), area(m2) and mass(kg).
double c = 0.47;
double R = 0.02;
double A = 3.141592654 * R * R;
double m = 0.004;
double gr = (1 + sqrt(5.0)) / 2.0;
double dodecmag = sqrt(3);
double icosmag = sqrt(gr * gr + 1);

// Air density(kg.m - 3), acceleration due to gravity(m.s - 2).
double rho_air = 1.28;
double g = 9.806;
// For convenience, define  this constant.
double kk = 0.5 * c * rho_air * A;

// initial velocity m/s
double v0 = 25.0;
// initial height m
double z0 = 4.0;

// variables to store calculations for each firework particle
double x[dimp];
double y[dimp];
double z[dimp];
double r;
double value;
double speed;
// use vertices of a dodecahedron plus icosahedron for the initial directions. Initially these are length sqrt(3), later they will be scaled to v0
double xdot[dimp] = { 1, 1, 1, 1, -1, -1, -1, -1, 0, 0, 0, 0, 1 / gr, 1 / gr, -1 / gr, -1 / gr, gr, gr, -gr, -gr, 0, 0, 0, 0, 1, 1, -1, -1, gr, -gr, gr, -gr };
double ydot[dimp] = { 1, 1, -1, -1, 1, 1, -1, -1, 1 / gr, 1 / gr, -1 / gr, -1 / gr, gr, -gr, gr, -gr, 0, 0, 0, 0, gr, -gr, gr, -gr, 0, 0, 0, 0, 1, 1, -1, -1 };
double zdot[dimp] = { 1, -1, 1, -1, 1, -1, 1, -1, gr, -gr, gr, -gr, 0, 0, 0, 0, 1 / gr, -1 / gr, 1 / gr, -1 / gr, 1, 1, -1, -1, gr, -gr, gr,-gr, 0, 0, 0, 0 };
// double phi[dimp];
// double theta[dimp];
// angles in radians to rotate the velocity vectors to prevent overlaps
double rotx = 0.2;
double roty = 0.15;

double vzavg = 0.0;
double vmax = 0.0;
double rmax = 0.0;
double rmin = 10000.0;
double dt = 0.05;
double dx = 2.0; // how much space does a pixel represent in meters

// space to store the image (value 1 if a pixel contains a firework particle, 0 otherwise)
int pixels[dim][dim];

int i = 0;
int j = 0;
int k = 0;
int ix = 0;
int iy = 0;
int dim2 = dim / 2;

int main()
{
    // open output file
    std::ofstream myFile("fw2_16x16x30.h");

    // Initialize positions
    for (i = 0;i < dimp;i++) {
        x[i] = -0.5;
        y[i] = 0.0;
        z[i] = z0;
    }
    // scale initial velocities
    for (i = 0;i < 20;i++) {
        xdot[i] = xdot[i] * v0 / dodecmag;
        ydot[i] = ydot[i] * v0 / dodecmag;
        zdot[i] = zdot[i] * v0 / dodecmag;
    }
    for (i = 20;i < dimp;i++) {
        xdot[i] = xdot[i] * v0 / icosmag;
        ydot[i] = ydot[i] * v0 / icosmag;
        zdot[i] = zdot[i] * v0 / icosmag;
    }
    // rotate by a small angle to prevent points from overlapping.  Rotate about x, then y.

    for (i = 0;i < dimp;i++) {
        value = ydot[i];
        ydot[i] = ydot[i] * cos(rotx) + zdot[i] * sin(rotx);
        zdot[i] = zdot[i] * cos(rotx) - value * sin(rotx);
    }
    for (i = 0;i < dimp;i++) {
        value = xdot[i];
        xdot[i] = xdot[i] * cos(roty) + zdot[i] * sin(roty);
        zdot[i] = zdot[i] * cos(roty) - value * sin(roty);
    }

    // now do the position calculations for each particle and output the results to the file.
    // numbers will be stored in binary format to save Arduino memory
    // For the initial (t=0), just output the center pixels

    myFile << "{";

    for (k = 0; k < numf; k++) {
        // initialize array to zeros
        for (i = 0; i < dim; i++) {
            for (j = 0;j < dim;j++) {
                pixels[i][j] = 0;
            }
        }
        // set array values to 1 if they correspond to a particle location
        for (i = 0; i < dimp; i++) {
            ix = int(x[i] / dx) + dim / 2;
            iy = int(z[i] / dx) + dim / 2;
            if (ix >= 0 && ix < dimp && iy >= 0 && iy < dimp) {
                pixels[ix][iy] = 1;
            }
        }
        // output the array to the file in binary format for Arduino
        myFile << "{";
        for (i = 0; i < dim; i++) {
            myFile << "0b";
            for (j = 0;j < dim;j++) {
                myFile << pixels[j][i];
            }
            if (i < dim - 1) {
                myFile << ",\n";
            }
        }
        if (k < numf - 1) {
            myFile << "},\n";
        }
        else {
            myFile << "}}\n";
        }
        // now compute particle positions and speeds at next time interval
        for (i = 0; i < dimp; i++) {
            x[i] += xdot[i] * dt;
            y[i] += ydot[i] * dt;
            z[i] += zdot[i] * dt;
            speed = sqrt(xdot[i] * xdot[i] + ydot[i] * ydot[i] + zdot[i] * zdot[i]);
            xdot[i] = xdot[i] - kk / m * speed * xdot[i] * dt;
            ydot[i] = ydot[i] - kk / m * speed * ydot[i] * dt;
            zdot[i] = zdot[i] - (kk / m * speed * zdot[i] + g) * dt;

        }

    }

    // close the file:
    myFile.close();

    for (i = 0; i < dimp; i++) {
        speed = sqrt(xdot[i] * xdot[i] + ydot[i] * ydot[i] + zdot[i] * zdot[i]);
        if (speed > vmax) vmax = speed;
        r = sqrt(x[i] * x[i] + y[i] * y[i] + z[i] * z[i]);
        if (r > rmax) rmax = r;
        if (r < rmin) rmin = r;
        vzavg += zdot[i];
    }
    vzavg = vzavg / double(dimp);

    std::cout << "Vmax = " << vmax << "\n";
    std::cout << "Vzavg = " << vzavg << "\n";
    std::cout << "Rmax = " << rmax << "\n";
    std::cout << "Rmin = " << rmin << "\n";
    std::cout << "Finished!\n";

}

// Run program: Ctrl + F5 or Debug > Start Without Debugging menu
// Debug program: F5 or Debug > Start Debugging menu

// Tips for Getting Started: 
//   1. Use the Solution Explorer window to add/manage files
//   2. Use the Team Explorer window to connect to source control
//   3. Use the Output window to see build output and other messages
//   4. Use the Error List window to view errors
//   5. Go to Project > Add New Item to create new code files, or Project > Add Existing Item to add existing code files to the project
//   6. In the future, to open this project again, go to File > Open > Project and select the .sln file