#include <stdio.h>
#include <stdlib.h>
#include <math.h>

char *file_pesi = "rete_pesi.bin";
char *file_immagini = "mnist/t10k-images.idx3-ubyte";
char *file_label = "mnist/t10k-labels.idx1-ubyte";

#include "mnist/mnist_manager.h"
// #include "cifar_10/cifar10_manager.h"
//#include "xor_manager.h"

// Siccome il char è un byte che rappresenta il valore tra 0 e 255. Per evitare confusioni definisco il tipo "byte" come in Java
typedef unsigned char byte;

double LRE = 2;
double soglia_sigmoide = 0.5;

typedef struct
{
    double *pesi;
    double bias;
    int size;
} Percettrone;

typedef struct
{
    Percettrone *percettroni;
    int size;
} Layer;

typedef struct
{
    Layer *layers;
    int size;
} ReteNeurale;

double randomico();

Percettrone inzializza_percettrone(int);
ReteNeurale inizializza_rete_neurale(int, int, int);
Layer inizializza_layer(int, int);

double sigmoide(Percettrone p, double*);
double **elabora_sigmoidi(ReteNeurale, Istanza);
void discesa_gradiente(ReteNeurale, double **, double **);
double calcola_gradiente_disceso(ReteNeurale, int, int, double **);
void correggi_pesi_percettrone_double(Percettrone *, int, double **, double);
void correggi_pesi_percettrone_byte(Percettrone *, Istanza, double, int);

int previsione(double);

void salvaReteNeurale(const char *, ReteNeurale *);
ReteNeurale *caricaReteNeurale(const char *);

/*
    ################# INIZIALIZZAZIONI ################################
*/
// Questa funzione genera un valore reale random compreso nell'intervallo [-1, 1]
double randomico()
{
    // Genero numeri nell'intervallo [-1,1]
    return ((double)(rand() % 101 * 0.01 * 2.0) -1.0);
}

Percettrone inizializza_percettrone(int n_pesi)
{
    Percettrone p;
    p.pesi = (double *)malloc(sizeof(double) * n_pesi);
    for (int i = 0; i < n_pesi; i++) {
        p.pesi[i] = randomico();
    }

    p.bias = randomico();

    p.size = n_pesi;

    return p;
}

Layer inizializza_layer(int n_percettroni, int n_pesi)
{
    Layer layer;
    layer.percettroni = (Percettrone *)malloc(sizeof(Percettrone) * n_percettroni);

    for (int i = 0; i < n_percettroni; i++)
    {
        layer.percettroni[i] = inizializza_percettrone(n_pesi);
    }

    layer.size = n_percettroni;

    return layer;
}

ReteNeurale inizializza_rete_neurale(int numero_layers, int numero_percettroni_iniziali, int numero_input)
{
    ReteNeurale r;
    r.layers = (Layer *)malloc(sizeof(Layer) * numero_layers);
    r.size = numero_layers;

    // Funzione esponenziale inversa layer 5
    for (int livello = 0; livello < numero_layers; livello++)
    {
        double esponente = (double)livello / (double)numero_layers;
        double frazione = (double)1 / (double)numero_percettroni_iniziali;

        int numero_percettroni_livello = (int)((double)numero_percettroni_iniziali * pow(frazione, esponente));
        if (livello == numero_layers - 1)
            numero_percettroni_livello = 1;

        // printf("esponente %f -> frazione: %f\n", esponente, frazione);
        printf("Layer %d -> percettroni: %d\n", livello, numero_percettroni_livello);

        if (livello == 0)
            r.layers[livello] = inizializza_layer(numero_percettroni_livello, numero_input);
        else
            r.layers[livello] = inizializza_layer(numero_percettroni_livello, r.layers[livello - 1].size);
    }

    return r;
}


/*
    ################# PREVISIONI ################################
*/


double sigmoide(Percettrone p, double *valori)
{
    double sommatoria = 0.0;
    for (int i = 0; i < p.size; i++) {
        sommatoria += (valori[i] * p.pesi[i]);
        //printf("valore [%f] peso[%f] ", valori[i], p.pesi[i]);
    }

    sommatoria += p.bias;
    double risultato = 1.0 / (1.0 + exp(-sommatoria));
    //printf(" sommatoria %f -> %f\n",sommatoria, risultato);

    return risultato;
}

int previsione(double valore)
{
    if (valore >= soglia_sigmoide)
        return 1;
    else
        return 0;
}
 
void discesa_gradiente(ReteNeurale rete, double **sigmoidi, double **gradienti)
{
    //printf("Qui?\n");
    // For che scorre i layer dal penultimo al primo QUINI SIZE -2
    for (int indice_layer = rete.size - 2; indice_layer >= 0; indice_layer--)
    {
        //printf("layer: %d ", indice_layer);
        // printf("Mi trovo nel layer %d, ho %d percettroni\n", indice_layer, rete.layers[indice_layer].size);

        // For che scorre i percettroni del layer partendo dal primo
        // Per ogni percettrone mi devo prendere il gradiente disceso dal livello sopra e moltiplicarlo per la derivata di attivazione
        for (int indice_percettrone = 0; indice_percettrone < rete.layers[indice_layer].size; indice_percettrone++)
        {
            //printf("percettrone: %d ", indice_percettrone);
            double derivata_attivazione = sigmoidi[indice_layer][indice_percettrone] * (1.0 - sigmoidi[indice_layer][indice_percettrone]);

            // Passo anche l'indice del percettrone perchè corrisponde all'indice del peso del livello sopra
            double gradiente_disceso = calcola_gradiente_disceso(rete, indice_layer + 1, indice_percettrone, gradienti);

            gradienti[indice_layer][indice_percettrone] = gradiente_disceso * derivata_attivazione;
        }
        //printf("\n");
    }
    //printf("Qui ?\n");
}

double calcola_gradiente_disceso(ReteNeurale rete, int livello, int indice_peso, double **gradienti)
{
    // printf("Qui ci arrivo\n");
    double sommatoria = 0.0;
    // printf("Layer %d: N_percettroni: %d\n", livello, rete.layers[livello].size);
    // Calcolo la sommatoria dei gradienti dei percettroni per i pesi
    for (int indice_percettrone = 0; indice_percettrone < rete.layers[livello].size; indice_percettrone++)
    {
        sommatoria += (gradienti[livello][indice_peso] * rete.layers[livello].percettroni[indice_percettrone].pesi[indice_peso]);
    }

    return sommatoria;
}

double **elabora_sigmoidi(ReteNeurale rete, Istanza istanza)
{
    // sigmoidi è un array bidimensionale, la prima dimensione identifica il layer, la seconda il percettrone nel layer
    // sigmoidi[indice_layer][indice_percettrone]
    double **sigmoidi = (double **)malloc(sizeof(double *) * rete.size);
    double *inputs = (double *)malloc(sizeof(double *) * N_INPUTS);
    for(int i = 0; i < N_INPUTS; i++) {
        inputs[i] = (double)istanza.dati[i];
    }

    sigmoidi[0] = (double *)malloc(sizeof(double) * rete.layers[0].size);
    for(int indice_percettrone = 0; indice_percettrone < rete.layers[0].size; indice_percettrone ++) {
        sigmoidi[0][indice_percettrone] = sigmoide(rete.layers[0].percettroni[indice_percettrone], inputs);
    }

    for(int indice_layer = 1; indice_layer < rete.size; indice_layer ++) {
        sigmoidi[indice_layer] = (double *)malloc(sizeof(double) * rete.layers[indice_layer].size);
        for(int indice_percettrone = 0; indice_percettrone < rete.layers[indice_layer].size; indice_percettrone ++) {
            sigmoidi[indice_layer][indice_percettrone] = sigmoide(rete.layers[indice_layer].percettroni[indice_percettrone], sigmoidi[indice_layer-1]);
        }
    }

    return sigmoidi;
}


/*
    ################# CORREZIONI ################################
*/


void correggi_pesi_percettrone_double(Percettrone *p, int layer, double **input, double gradiente_percettrone)
{

    for (int indice_peso = 0; indice_peso < p->size; indice_peso++)
    {
        // Determino il gradiente del peso
        double gradiente_peso = gradiente_percettrone * input[layer - 1][indice_peso];

        // Modifico il peso
        p->pesi[indice_peso] += (gradiente_peso * LRE);
    }

    p->bias += (gradiente_percettrone * LRE);
}

void correggi_pesi_percettrone_byte(Percettrone *p, Istanza input, double gradiente_percettrone, int indice_percettrone)
{
    for (int indice_peso = 0; indice_peso < p->size; indice_peso++)
    {
        // Determino il gradiente del peso
        double gradiente_peso = gradiente_percettrone * (double)input.dati[indice_peso];

        // Modifico il peso Qui si impalla perchè per qualche ragione arriva size elevatissimo
        p->pesi[indice_peso] += (gradiente_peso * LRE);
    }

    p->bias += (gradiente_percettrone * LRE);
}



/*
    ################# IMPORT EXPORT ################################
*/
void salvaReteNeurale(const char *filename, ReteNeurale *rete)
{
    FILE *file = fopen(filename, "wb");
    if (!file)
    {
        perror("Errore nell'apertura del file");
        exit(EXIT_FAILURE);
    }

    // Scrivi il numero di layer
    fwrite(&rete->size, sizeof(int), 1, file);

    // Scrivi ogni layer
    for (int i = 0; i < rete->size; i++)
    {
        Layer *layer = &rete->layers[i];
        fwrite(&layer->size, sizeof(int), 1, file);

        // Scrivi ogni percettrone nel layer
        for (int j = 0; j < layer->size; j++)
        {
            Percettrone *perc = &layer->percettroni[j];
            fwrite(&perc->size, sizeof(int), 1, file);
            fwrite(perc->pesi, sizeof(double), perc->size, file);
            fwrite(&perc->bias, sizeof(double), 1, file);
        }
    }

    fclose(file);
}

ReteNeurale *caricaReteNeurale(const char *filename)
{
    FILE *file = fopen(filename, "rb");

    if (!file)
    {
        perror("Errore nell'apertura del file");
        return NULL;
    }

    ReteNeurale *rete = malloc(sizeof(ReteNeurale));
    if (!rete)
    {
        perror("Errore nell'allocazione della memoria");
        return NULL;
    }

    // Leggi il numero di layer
    fread(&rete->size, sizeof(int), 1, file);
    rete->layers = malloc(rete->size * sizeof(Layer));
    if (!rete->layers)
    {
        perror("Errore nell'allocazione della memoria");
        return NULL;
    }

    // Leggi ogni layer
    for (int i = 0; i < rete->size; i++)
    {
        Layer *layer = &rete->layers[i];
        fread(&layer->size, sizeof(int), 1, file);
        layer->percettroni = malloc(layer->size * sizeof(Percettrone));
        if (!layer->percettroni)
        {
            perror("Errore nell'allocazione della memoria");
            return NULL;
        }

        // Leggi ogni percettrone nel layer
        for (int j = 0; j < layer->size; j++)
        {
            Percettrone *perc = &layer->percettroni[j];
            fread(&perc->size, sizeof(int), 1, file);
            perc->pesi = malloc(perc->size * sizeof(double));
            if (!perc->pesi)
            {
                perror("Errore nell'allocazione della memoria");
                return NULL;
            }
            fread(perc->pesi, sizeof(double), perc->size, file);
            fread(&perc->bias, sizeof(double), 1, file);
        }
    }

    fclose(file);
    return rete;
}









/* double sigmoide_byte(Percettrone, byte *, int);
double sigmoide_double(Percettrone, double *, int);
double *funzioni_attivazione_layer_byte(Layer, byte *);
double *funzioni_attivazione_layer_double(Layer, double *); */
/* void correggi_layer_interni(ReteNeurale *, double **, double **);
void correggi_layer_input(Layer *, double **, double **, byte *, int); */
/* // Questa funzione calcola tutte le funzioni di attivazione dei percettroni del layer che prende dei byte come inputs
double *funzioni_attivazione_layer_byte(Layer layer, byte *inputs)
{

    double *funzioni = (double *)malloc(sizeof(double) * layer.size);

    for (int i = 0; i < layer.size; i++)
    {
        funzioni[i] = sigmoide_byte(layer.percettroni[i], inputs, layer.percettroni[i].size);
        // printf("\tsigmoide layer input %f\n", funzioni[i]);
    }

    return funzioni;
}

// Questa funzione calcola tutte le funzioni di attivazione dei percettroni del layer che prende dei double come inputs (le sigmoidi del livello precedente)
double *funzioni_attivazione_layer_double(Layer layer, double *inputs)
{

    double *funzioni = (double *)malloc(sizeof(double) * layer.size);

    for (int i = 0; i < layer.size; i++)
    {
        funzioni[i] = sigmoide_double(layer.percettroni[i], inputs, layer.percettroni[i].size);
        // printf("\tsigmoide layer %d: %f\n", i, funzioni[i]);
    }

    return funzioni;
}
 */
/* // Questa funzione prende la matrice dei gradienti e la matrice delle sigmoidi per correggere tutti i layer tranne quello di ingresso
void correggi_layer_interni(ReteNeurale *rete, double **gradienti, double **sigmoidi)
{

    for (int indice_layer = rete->size - 1; indice_layer > 0; indice_layer--)
    {
        for (int indice_percettrone = 0; indice_percettrone < rete->layers[indice_layer].size; indice_percettrone++)
        { // Numero percettroni

            for (int indice_peso = 0; indice_peso < rete->layers[indice_layer].percettroni[indice_percettrone].size; indice_peso++)
            { // Numero pesi
                gradienti[indice_layer][indice_percettrone] = gradienti[rete->size - 1][0] * (sigmoidi[indice_layer][indice_percettrone] * (1 - sigmoidi[indice_layer][indice_percettrone]));
                rete->layers[indice_layer].percettroni[indice_percettrone].pesi[indice_peso] += (gradienti[indice_layer][indice_percettrone] * LRE * sigmoidi[indice_layer - 1][indice_percettrone]);
                // rete->layers[indice_layer].percettroni[indice_percettrone].pesi[indice_peso] += (gradienti[rete->size-1][0] * LRE * sigmoidi[indice_layer-1][indice_percettrone]);
            }
            rete->layers[indice_layer].percettroni[indice_percettrone].bias += (gradienti[indice_layer][indice_percettrone] * LRE);
            // printf("bias: %f\n", rete->layers[indice_layer].percettroni[indice_percettrone].bias);
        }
    }
}
// Questa funzione prende tutti i parametri della precedente + gli input passati dal dataset per correggere il layer di ingresso
void correggi_layer_input(Layer *layer, double **gradienti, double **sigmoidi, byte *inputs, int n_layers)
{
    // L'indice del layer d'ingresso che prende byte per input
    int indice_layer = 0;
    for (int indice_percettrone = 0; indice_percettrone < layer->size; indice_percettrone++)
    { // Numero percettroni
        for (int indice_peso = 0; indice_peso < layer->percettroni->size; indice_peso++)
        { // Numero pesi

            gradienti[indice_layer][indice_percettrone] = gradienti[n_layers - 1][0] * (sigmoidi[indice_layer][indice_percettrone] * (1 - sigmoidi[indice_layer][indice_percettrone]));
            layer->percettroni[indice_percettrone].pesi[indice_peso] += (gradienti[indice_layer][indice_percettrone] * LRE * inputs[indice_peso]);
            // layer->percettroni[indice_percettrone].pesi[indice_peso] += (gradienti[n_layers-1][0] * LRE * inputs[indice_peso]);
        }
        layer->percettroni[indice_percettrone].bias += (gradienti[n_layers - 1][0] * LRE);
    }
}
 */
/* 
// Questa funzione viene usata per il primo livello perchè ha un vettore di byte (unsigned char) in input
double sigmoide_byte(Percettrone p, byte *valori, int n_input)
{

    double sommatoria = 0.0;
    // printf("valori: [%d][%d]", valori[0], valori[1]);
    // printf("pesi: [%f][%f]", p.pesi[0], p.pesi[1]);

    for (int i = 0; i < n_input; i++)
    {
        sommatoria += ((double)valori[i] * p.pesi[i]);
    }
    // printf("sommatoria= %f\n", sommatoria);
    double funzione = sommatoria + p.bias;
    double potenza_e = exp(-funzione);
    // printf("potenza_e= %f\n", potenza_e);
    // formula sigmoide
    double risultato = 1.0 / (1.0 + potenza_e);
    // printf("risultato= %f\n", risultato);
    return risultato;
}

// Questa funzione viene usata per gli altri livelli dove gli input sono double, ossia i valori della sigmoide dei livelli precedenti
double sigmoide_double(Percettrone p, double *valori, int n_input)
{
    double sommatoria = 0.0;
    for (int i = 0; i < n_input; i++)
    {
        sommatoria += (valori[i] * p.pesi[i]);
    }

    double funzione = sommatoria + p.bias;
    // printf("sommatoria= %f\n", sommatoria);
    double potenza_e = exp(-funzione);
    // printf("potenza_e= %f\n", potenza_e);
    // formula sigmoide
    double risultato = 1.0 / (1.0 + potenza_e);
    // printf("risultato= %f\n", risultato);

    return risultato;
}
 */