#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <basic-geometry.h>

#ifdef _WIN64
#include <windows.h>
#else
#include <time.h>
#endif // _WINDOWS_

typedef struct {
    BGC_FP32_Versor versor1, versor2, result;
} structure_fp32_t;

structure_fp32_t* allocate_structures(const unsigned int amount)
{
    return calloc(amount, sizeof(structure_fp32_t));
}

structure_fp32_t* make_structures(const unsigned int amount)
{
    structure_fp32_t* list = allocate_structures(amount);

    if (list == 0) {
        return 0;
    }

    const float multiplier = 2.0f / RAND_MAX;

    for (unsigned int i = 0; i < amount; i++) {
        bgc_fp32_versor_make(
            &list[i].versor1,
            rand() * multiplier - 1.0f,
            rand() * multiplier - 1.0f,
            rand() * multiplier - 1.0f,
            rand() * multiplier - 1.0f
        );

        bgc_fp32_versor_make(
            &list[i].versor2,
            rand() * multiplier - 1.0f,
            rand() * multiplier - 1.0f,
            rand() * multiplier - 1.0f,
            rand() * multiplier - 1.0f
        );

        bgc_fp32_versor_reset(&list[i].result);
    }

    return list;
}

void print_versor_fp32(const BGC_FP32_Versor* versor)
{
    printf("Versor (s0 = %0.12f, x1 = %0.12f, x2 = %0.12f, x3 = %0.12f)\n", versor->_s0, versor->_x1, versor->_x2, versor->_x3);
}

void print_versor_fp64(const BGC_FP64_Versor* versor)
{
    printf("Versor (s0 = %0.20f, x1 = %0.20f, x2 = %0.20f, x3 = %0.20f)\n", versor->_s0, versor->_x1, versor->_x2, versor->_x3);
}

void print_vector_fp32(const BGC_FP32_Vector3* vector)
{
    printf("(%f, %f, %f) / %f\n", vector->x1,  vector->x2,  vector->x3,  bgc_fp32_vector3_get_modulus(vector));
}

void print_vector_fp64(const BGC_FP64_Vector3* vector)
{
    printf("(%lf, %lf, %lf) / %lf\n", vector->x1,  vector->x2,  vector->x3,  bgc_fp64_vector3_get_modulus(vector));
}

void list_work(const uint_fast32_t amount, structure_fp32_t* list)
{
    for (uint_fast32_t j = 0; j < 1000; j++) {
        for (uint_fast32_t i = 0; i < amount; i++) {
            bgc_fp32_versor_combine(&list[i].versor1, &list[i].versor1, &list[i].result);
        }
    }
}
/*
int main()
{
    const unsigned int amount = 1000000;
    structure_fp32_t* list = make_structures(amount);

#ifdef _WIN64
    ULONGLONG start, end;

    start = GetTickCount64();

    srand((unsigned int)(start & 0xfffffff));

    start = GetTickCount64();
#else
    struct timespec start, end;
    clock_gettime(0, &start);
    srand((unsigned int)(start.tv_nsec & 0xfffffff));

    clock_gettime(CLOCK_REALTIME, &start);
#endif // _WIN64

    list_work(amount, list);

#ifdef _WIN64
    end = GetTickCount64();

    printf("Time: %lld\n", end - start);
#else
    clock_gettime(CLOCK_REALTIME, &end);

    printf("Time: %lf\n", (end.tv_sec - start.tv_sec) * 1000.0 + (end.tv_nsec - start.tv_nsec) * 0.000001);
#endif // _WIN64

    print_versor_fp32(&list[10].versor1);
    print_versor_fp32(&list[10].versor2);
    print_versor_fp32(&list[10].result);

    free(list);

    return 0;
}
*/

/*
int main() {
    BGC_FP32_Complex complex, exponent, result;

    bgc_fp32_complex_make(0, 1, &complex);

    bgc_fp32_complex_make(4, 0, &exponent);

    bgc_fp32_complex_get_exponation(&complex, exponent.real, exponent.imaginary, &result);

    printf("(%f, %f) ^ (%f, %f) = (%f, %f)\n", complex.real, complex.imaginary, exponent.real, exponent.imaginary, result.real, result.imaginary);

    return 0;
}
*/
/*
int main() {
    BGC_FP32_Versor start = { 1.0f, 0.0f, 0.0f, 0.0f };
    BGC_FP32_Versor end = { 0.0f, 1.0f, 0.0f, 0.0f };
    BGC_FP32_Versor result;
    bgc_fp32_versor_spherical_interpolation(&start, &end, 0.5f, &result);
    printf("Result: %0.12f, %0.12f, %0.12f, %0.12f\n", result.s0, result.x1, result.x2, result.x3);
    return 0;
}
*/

void test_basis_difference_fp32()
{

    BGC_FP32_Vector3 initial_primary, initial_auxiliary;
    BGC_FP32_Vector3 final_primary, final_auxiliary;
    BGC_FP32_Versor turn;

    // No turn
    bgc_fp32_vector3_make(&initial_primary,   1.0f, 0.0f, 0.0f);
    bgc_fp32_vector3_make(&initial_auxiliary, 0.0f, 1.0f, 0.0f);

    bgc_fp32_vector3_make(&final_primary,   1.0f, 0.0f, 0.0f);
    bgc_fp32_vector3_make(&final_auxiliary, 0.0f, 1.0f, 0.0f);

    bgc_fp32_versor_make_basis_difference(&turn, &initial_primary, &initial_auxiliary, &final_primary, &final_auxiliary);

    printf("\nNo turn:\n");
    print_versor_fp32(&turn);

    // Turn around (1, 1, 0) axis on 180 degrees
    bgc_fp32_vector3_make(&initial_primary,   1.0f, 0.0f, 0.0f);
    bgc_fp32_vector3_make(&initial_auxiliary, 0.0f, 1.0f, 0.0f);

    bgc_fp32_vector3_make(&final_primary,   0.0f, 1.0f, 0.0f);
    bgc_fp32_vector3_make(&final_auxiliary, 1.0f, 0.0f, 0.0f);

    bgc_fp32_versor_make_basis_difference(&turn, &initial_primary, &initial_auxiliary, &final_primary, &final_auxiliary);

    printf("\nTurn around (1, 1, 0) axis on 180 degrees:\n");
    print_versor_fp32(&turn);

    // 180 degree turn
    bgc_fp32_vector3_make(&initial_primary,   1.0f, 0.0f, 0.0f);
    bgc_fp32_vector3_make(&initial_auxiliary, 0.0f, 1.0f, 0.0f);

    bgc_fp32_vector3_make(&final_primary,  -1.0f, 0.0f, 0.0f);
    bgc_fp32_vector3_make(&final_auxiliary, 0.0f, 1.0f, 0.0f);

    bgc_fp32_versor_make_basis_difference(&turn, &initial_primary, &initial_auxiliary, &final_primary, &final_auxiliary);

    printf("\n180 degree turn around (0, 1, 0):\n");
    print_versor_fp32(&turn);

    // 90 degree turn around x3 axis
    bgc_fp32_vector3_make(&initial_primary,   2.0f, 0.0f, 0.0f);
    bgc_fp32_vector3_make(&initial_auxiliary, 0.0f, 3.1f, 0.0f);

    bgc_fp32_vector3_make(&final_primary,  0.0f, 10.0f, 0.0f);
    bgc_fp32_vector3_make(&final_auxiliary,-1.0f, 0.0f, 0.0f);

    bgc_fp32_versor_make_basis_difference(&turn, &initial_primary, &initial_auxiliary, &final_primary, &final_auxiliary);

    printf("\n90 degree turn around (0, 0, 1):\n");
    print_versor_fp32(&turn);

    // Unorthogonal pairs turn at 90 degrees around x3 axis
    bgc_fp32_vector3_make(&initial_primary,    2.0f, 0.0f, 0.0f);
    bgc_fp32_vector3_make(&initial_auxiliary, -2.0f, 3.1f, 0.0f);

    bgc_fp32_vector3_make(&final_primary, 0.0f, 10.0f, 0.0f);
    bgc_fp32_vector3_make(&final_auxiliary, -1.0f, 5.0f, 0.0f);

    bgc_fp32_versor_make_basis_difference(&turn, &initial_primary, &initial_auxiliary, &final_primary, &final_auxiliary);

    printf("\nUnorthogonal pairs turn at 90 degrees around (0, 0, 1):\n");
    print_versor_fp32(&turn);

    // Zero vectors
    bgc_fp32_vector3_make(&initial_primary,   0.0f, 0.0f, 0.0f);
    bgc_fp32_vector3_make(&initial_auxiliary, 0.0f, 1.0f, 0.0f);
    bgc_fp32_vector3_make(&final_primary,     1.0f, 0.0f, 0.0f);
    bgc_fp32_vector3_make(&final_auxiliary,   0.0f, 1.0f, 0.0f);

    int code;

    code = bgc_fp32_versor_make_basis_difference(&turn, &initial_primary, &initial_auxiliary, &final_primary, &final_auxiliary);

    if (code >= 0) {
        printf("\nZero vectors: this cannot be!\n");
        print_versor_fp32(&turn);
    }
    else {
        printf("\nZero vector validation works fine\n");
    }

    // Parallel vectors
    bgc_fp32_vector3_make(&initial_primary,   1.0f, 0.0f, 0.0f);
    bgc_fp32_vector3_make(&initial_auxiliary, 2.0f, 0.0f, 0.0f);
    bgc_fp32_vector3_make(&final_primary,     1.0f, 0.0f, 0.0f);
    bgc_fp32_vector3_make(&final_auxiliary,   0.0f, 1.0f, 0.0f);

    code = bgc_fp32_versor_make_basis_difference(&turn, &initial_primary, &initial_auxiliary, &final_primary, &final_auxiliary);

    if (code >= 0) {
        printf("\nParallel vectors: this cannot be!\n");
        print_versor_fp32(&turn);
    }
    else {
        printf("\nParallelism validation works fine\n");
    }

    // Small angle turn (about 1 degree):
    bgc_fp32_vector3_make(&initial_primary,   1.0f, 0.0f, 0.0f);
    bgc_fp32_vector3_make(&initial_auxiliary, 0.0f, 1.0f, 0.0f);

    bgc_fp32_vector3_make(&final_primary,    0.999848f, 0.017452f, 0.0f);
    bgc_fp32_vector3_make(&final_auxiliary, -0.017452f, 0.999848f, 0.0f);

    bgc_fp32_versor_make_basis_difference(&turn , &initial_primary, &initial_auxiliary, &final_primary, &final_auxiliary);

    printf("\nSmall angle turn (about 1 degree):\n");
    print_versor_fp32(&turn);

    // About 179 degrees turn
    bgc_fp32_vector3_make(&initial_primary,   1.0f, 0.0f, 0.0f);
    bgc_fp32_vector3_make(&initial_auxiliary, 0.0f, 1.0f, 0.0f);

    bgc_fp32_vector3_make(&final_primary,  -0.999848f, -0.017452f, 0.0f);
    bgc_fp32_vector3_make(&final_auxiliary, 0.017452f, -0.999848f, 0.0f);

    bgc_fp32_versor_make_basis_difference(&turn, &initial_primary, &initial_auxiliary, &final_primary, &final_auxiliary);

    printf("\nAbout 179 degrees turn:\n");
    print_versor_fp32(&turn);

    // 120 degrees around (-1, -1, 1)
    bgc_fp32_vector3_make(&initial_primary,   1.0f, 0.0f, 0.0f);
    bgc_fp32_vector3_make(&initial_auxiliary, 0.0f, 1.0f, 0.0f);

    bgc_fp32_vector3_make(&final_primary,   0.0f, 1.0f, 0.0f);
    bgc_fp32_vector3_make(&final_auxiliary, 0.0f, 0.0f, -1.0f);

    bgc_fp32_versor_make_basis_difference(&turn, &initial_primary, &initial_auxiliary, &final_primary, &final_auxiliary);

    printf("\n120 degees turn:\n");
    print_versor_fp32(&turn);


    // About 1 degree turn difference between initial_primary and initial_auxiliary directions
    bgc_fp32_vector3_make(&initial_primary,   1.0f, 0.0f, 0.0f);
    bgc_fp32_vector3_make(&initial_auxiliary, 0.999848f, 0.017452f, 0.0f);
    bgc_fp32_vector3_make(&final_primary, 0.0f, 1.0f, 0.0f);
    bgc_fp32_vector3_make(&final_auxiliary, -1.0f, 0.0f, 0.0f);

    bgc_fp32_versor_make_basis_difference(&turn, &initial_primary, &initial_auxiliary, &final_primary, &final_auxiliary);

    printf("\nAbout 1 degree turn difference between initial_primary and initial_auxiliary directions:\n");
    print_versor_fp32(&turn);

    // About 0.01 degree turn difference between initial_primary and initial_auxiliary directions
    bgc_fp32_vector3_make(&initial_primary,   1.0f, 0.0f, 0.0f);
    bgc_fp32_vector3_make(&initial_auxiliary, 1.0f, 0.000001f, 0.0f);
    bgc_fp32_vector3_make(&final_primary,     0.0f, -1.0f, 0.0f);
    bgc_fp32_vector3_make(&final_auxiliary,   1.0f, 0.0f, 0.0f);

    bgc_fp32_versor_make_basis_difference(&turn, &initial_primary, &initial_auxiliary, &final_primary, &final_auxiliary);

    printf("\nAbout 0.01 degree turn difference between initial_primary and initial_auxiliary directions:\n");
    print_versor_fp32(&turn);
}

void test_basis_difference_fp64()
{

    BGC_FP64_Vector3 initial_primary, initial_auxiliary;
    BGC_FP64_Vector3 final_primary, final_auxiliary;
    BGC_FP64_Versor turn;

    // No turn
    bgc_fp64_vector3_make(&initial_primary,   1.0, 0.0, 0.0);
    bgc_fp64_vector3_make(&initial_auxiliary, 0.0, 1.0, 0.0);
    bgc_fp64_vector3_make(&final_primary,     1.0, 0.0, 0.0);
    bgc_fp64_vector3_make(&final_auxiliary,   0.0, 1.0, 0.0);

    bgc_fp64_versor_make_basis_difference(&turn, &initial_primary, &initial_auxiliary, &final_primary, &final_auxiliary);

    printf("\nNo turn:\n");
    print_versor_fp64(&turn);

    // Turn around (1, 1, 0) axis on 180 degrees
    bgc_fp64_vector3_make(&initial_primary,   1.0, 0.0, 0.0);
    bgc_fp64_vector3_make(&initial_auxiliary, 0.0, 1.0, 0.0);
    bgc_fp64_vector3_make(&final_primary,     0.0, 1.0, 0.0);
    bgc_fp64_vector3_make(&final_auxiliary,   1.0, 0.0, 0.0);

    bgc_fp64_versor_make_basis_difference(&turn, &initial_primary, &initial_auxiliary, &final_primary, &final_auxiliary);

    printf("\nTurn around (1, 1, 0) axis on 180 degrees:\n");
    print_versor_fp64(&turn);

    // 180 degree turn
    bgc_fp64_vector3_make(&initial_primary,   1.0, 0.0, 0.0);
    bgc_fp64_vector3_make(&initial_auxiliary, 0.0, 1.0, 0.0);

    bgc_fp64_vector3_make(&initial_auxiliary, -1.0, 0.0, 0.0);
    bgc_fp64_vector3_make(&final_auxiliary,    0.0, 1.0, 0.0);

    bgc_fp64_versor_make_basis_difference(&turn, &initial_primary, &initial_auxiliary, &final_primary, &final_auxiliary);

    printf("\n180 degree turn around (0, 1, 0):\n");
    print_versor_fp64(&turn);

    // 90 degree turn around x3 axis
    bgc_fp64_vector3_make(&initial_primary,   2.0, 0.0, 0.0);
    bgc_fp64_vector3_make(&initial_auxiliary, 0.0, 3.1, 0.0);

    bgc_fp64_vector3_make(&final_primary,   0.0, 10.0, 0.0);
    bgc_fp64_vector3_make(&final_auxiliary, -1.0, 0.0, 0.0);

    bgc_fp64_versor_make_basis_difference(&turn, &initial_primary, &initial_auxiliary, &final_primary, &final_auxiliary);

    printf("\n90 degree turn around (0, 0, 1):\n");
    print_versor_fp64(&turn);

    // Unorthogonal pairs turn at 90 degrees around x3 axis
    bgc_fp64_vector3_make(&initial_primary,    2.0, 0.0, 0.0);
    bgc_fp64_vector3_make(&initial_auxiliary, -2.0, 3.1, 0.0);

    bgc_fp64_vector3_make(&final_primary,   0.0, 10.0, 0.0);
    bgc_fp64_vector3_make(&final_auxiliary, -1.0, 5.0, 0.0);

    bgc_fp64_versor_make_basis_difference(&turn, &initial_primary, &initial_auxiliary, &final_primary, &final_auxiliary);

    printf("\nUnorthogonal pairs turn at 90 degrees around (0, 0, 1):\n");
    print_versor_fp64(&turn);

    // Zero vectors
    bgc_fp64_vector3_make(&initial_primary,   0.0, 0.0, 0.0);
    bgc_fp64_vector3_make(&initial_auxiliary, 0.0, 1.0, 0.0);
    bgc_fp64_vector3_make(&final_primary,     1.0, 0.0, 0.0);
    bgc_fp64_vector3_make(&final_auxiliary,   0.0, 1.0, 0.0);

    int code;

    code = bgc_fp64_versor_make_basis_difference(&turn, &initial_primary, &initial_auxiliary, &final_primary, &final_auxiliary);

    if (code >= 0) {
        printf("\nZero vectors: this cannot be!\n");
        print_versor_fp64(&turn);
    }
    else {
        printf("\nZero vector validation works fine\n");
    }

    // Parallel vectors
    bgc_fp64_vector3_make(&initial_primary,   1.0, 0.0, 0.0);
    bgc_fp64_vector3_make(&initial_auxiliary, 2.0, 0.0, 0.0);
    bgc_fp64_vector3_make(&final_primary,     1.0, 0.0, 0.0);
    bgc_fp64_vector3_make(&final_auxiliary,   0.0, 1.0, 0.0);

    code = bgc_fp64_versor_make_basis_difference(&turn, &initial_primary, &initial_auxiliary, &final_primary, &final_auxiliary);

    if (code >= 0) {
        printf("\nParallel vectors: this cannot be!\n");
        print_versor_fp64(&turn);
    }
    else {
        printf("\nParallelism validation works fine\n");
    }

    // Small angle turn (about 1 degree):
    bgc_fp64_vector3_make(&initial_primary,   1.0, 0.0, 0.0);
    bgc_fp64_vector3_make(&initial_auxiliary, 0.0, 1.0, 0.0);

    bgc_fp64_vector3_make(&final_primary,    0.999848, 0.017452, 0.0);
    bgc_fp64_vector3_make(&final_auxiliary, -0.017452, 0.999848, 0.0);

    bgc_fp64_versor_make_basis_difference(&turn, &initial_primary, &initial_auxiliary, &final_primary, &final_auxiliary);

    printf("\nSmall angle turn (about 1 degree):\n");
    print_versor_fp64(&turn);

    // About 179 degrees turn
    bgc_fp64_vector3_make(&initial_primary,   1.0, 0.0, 0.0);
    bgc_fp64_vector3_make(&initial_auxiliary, 0.0, 1.0, 0.0);

    bgc_fp64_vector3_make(&final_primary,  -0.999848, -0.017452, 0.0);
    bgc_fp64_vector3_make(&final_auxiliary, 0.017452, -0.999848, 0.0);

    bgc_fp64_versor_make_basis_difference(&turn, &initial_primary, &initial_auxiliary, &final_primary, &final_auxiliary);

    printf("\nAbout 179 degrees turn:\n");
    print_versor_fp64(&turn);

    // 120 degrees around (-1, -1, 1)
    bgc_fp64_vector3_make(&initial_primary,   1.0, 0.0, 0.0);
    bgc_fp64_vector3_make(&initial_auxiliary, 0.0, 1.0, 0.0);

    bgc_fp64_vector3_make(&final_primary,   0.0, 1.0, 0.0);
    bgc_fp64_vector3_make(&final_auxiliary, 0.0, 0.0, -1.0);

    bgc_fp64_versor_make_basis_difference(&turn, &initial_primary, &initial_auxiliary, &final_primary, &final_auxiliary);

    printf("\n120 degees turn:\n");
    print_versor_fp64(&turn);


    // About 1 degree turn difference between initial_primary and initial_auxiliary directions
    bgc_fp64_vector3_make(&initial_primary,   1.0, 0.0, 0.0);
    bgc_fp64_vector3_make(&initial_auxiliary, 0.999848, 0.017452, 0.0);
    bgc_fp64_vector3_make(&final_primary,     0.0, 1.0, 0.0);
    bgc_fp64_vector3_make(&final_auxiliary,  -1.0, 0.0, 0.0);

    bgc_fp64_versor_make_basis_difference(&turn, &initial_primary, &initial_auxiliary, &final_primary, &final_auxiliary);

    printf("\nAbout 1 degree turn difference between initial_primary and initial_auxiliary directions:\n");
    print_versor_fp64(&turn);

    // About 0.001 degree turn difference between initial_primary and initial_auxiliary directions
    bgc_fp64_vector3_make(&initial_primary,   1.0, 0.0, 0.0);
    bgc_fp64_vector3_make(&initial_auxiliary, 1.0, 0.000001, 0.0);
    bgc_fp64_vector3_make(&final_primary,     0.0, -1.0, 0.0);
    bgc_fp64_vector3_make(&final_auxiliary,   1.0, 0.0, 0.0);

    bgc_fp64_versor_make_basis_difference(&turn, &initial_primary, &initial_auxiliary, &final_primary, &final_auxiliary);

    printf("\nAbout 0.01 degree turn difference between initial_primary and initial_auxiliary directions:\n");
    print_versor_fp64(&turn);
}

#include "affine3.h"

int main()
{
    //BGC_FP32_Versor start = { 1.0f, 0.0f, 0.0f, 0.0f };
    //BGC_FP32_Versor end = { 0.0f, 1.0f, 0.0f, 0.0f };
    
    BGC_FP32_Versor start = { 1.0f, 0.0f, 0.0f, 0.0f };
    BGC_FP32_Versor end = { 0.9999f, 0.01414f, 0.0f, 0.0f };
    BGC_FP32_Slerp slerp;
    BGC_FP32_Versor result;
    bgc_fp32_slerp_make_full(&slerp, &start, &end);
    bgc_fp32_slerp_get_phase_versor(&result, &slerp, 0.5f);

    //print_versor_fp32(&result);
    
    test_basis_difference_fp64();
    
    //printf("Affine3 performance test: %f\n", test_bgc_affine3_performance(10000000, 10));

    //printf("sizeof(BGC_FP32_Affine3)   = %zu\n", sizeof(BGC_FP32_Affine3));
    //printf("offsetof(shift)         = %zu\n", offsetof(BGC_FP32_Affine3, shift));
    //printf("sizeof(BGC_FP32_Matrix3x3) = %zu\n", sizeof(BGC_FP32_Matrix3x3));

    return 0;
}