/*
 * The Npic library and tools
 *
 * Copyright (C) 2003 Edouard Thiel <Edouard.Thiel@lif.univ-mrs.fr>
 *
 * This program is free software under the terms of the 
 * GNU Lesser General Public License (LGPL) version 2.1.
*/

/*
 * npic-sedt.c - 11/09/2008
 *
*/

#include <npic.h>

void ShowUsage ()
{
    printf (
        "npic-sedt - Compute Square Euclidean Distance Transforms.\n"
        "Usage:\n"
        "  npic-sedt  -h | -help | --help   : print help\n"
        "  npic-sedt what how in1 out1      : compute distance transforms\n"
        "\n"
        "what:  -dt  : compute SEDT\n"
        "      -idt  : idem with infinite borders\n"
        "      -rdt  : compute reverse SEDT\n"
        "\n"
        "how:  -hirata  : with Hirata's algorithm (fastest)\n"
        "      -saito   : with Saito's et al algorithm\n"
        "      -quadra  : with quadratic algorithm (slowest)\n"
        "      -verify  : compare the three algorithms; save hirata\n"
        "\n"
        "in1, out1 : image in " NPIC_KNOWN_IMAGE_EXT " format\n"
    );
}

void ArgcExit (int argc, int n)
{
    if (argc < n) {
        fprintf (stderr, "ERROR: %d argument(s) missing, "
                         "type \"npic-sedt -h\" to get help.\n", n-argc);
        exit (1);
    }
}

int main (int argc, char *argv[])
{
    char *in1, *out1, *what, *how;
    enum { E_DT, E_IDT, E_RDT} what_e;
    char *what_s[] = { "SEDT", "SEDT_inf", "SEDT_rev" };
    enum { E_HIRATA, E_SAITO, E_QUADRA, E_VERIFY } how_e;
    char *how_s[] = { "Hirata", "Saito et al", "Quadratic", "Verify" };
    Npic_image *np1;
    double t = 0;

    ArgcExit (argc, 2);
    
    if (strcmp (argv[1], "-h") == 0 || strcmp (argv[1], "-help") == 0 ||
        strcmp (argv[1], "--help") == 0) { ShowUsage (); exit (0); }
    
    ArgcExit (argc, 5);
    what = argv[1];
    how  = argv[2];
    in1  = argv[3];
    out1 = argv[4];

    if      (!strcmp (what, "-dt" )) what_e = E_DT;
    else if (!strcmp (what, "-idt")) what_e = E_IDT;
    else if (!strcmp (what, "-rdt")) what_e = E_RDT;
    else { fprintf (stderr, "ERROR: bad option \"%s\"\n", what); exit (1); }

    if      (!strcmp (how, "-hirata")) how_e = E_HIRATA;
    else if (!strcmp (how, "-saito" )) how_e = E_SAITO;
    else if (!strcmp (how, "-quadra")) how_e = E_QUADRA;
    else if (!strcmp (how, "-verify")) how_e = E_VERIFY;
    else { fprintf (stderr, "ERROR: bad option \"%s\"\n", how); exit (1); }

    printf ("Loading \"%s\"\n", in1);
    np1 = NpicReadImage (in1);
    if (np1 == NULL) exit (1);

    /* All these algorithms need L image type */
    NpicConvertImage_l (np1);

    if (how_e != E_VERIFY) {

        printf ("Computing %s %s ... ", what_s[what_e], how_s[how_e]); fflush (stdout);
        switch (what_e *10+ how_e) {
            case  E_DT *10+ E_HIRATA : t = NpicSEDT_Hirata     (np1); break;
            case E_IDT *10+ E_HIRATA : t = NpicSEDT_Hirata_inf (np1); break;
            case E_RDT *10+ E_HIRATA : t = NpicSEDT_Hirata_rev (np1); break;
            case  E_DT *10+ E_SAITO  : t = NpicSEDT_Saito      (np1); break;
            case E_IDT *10+ E_SAITO  : t = NpicSEDT_Saito_inf  (np1); break;
            case E_RDT *10+ E_SAITO  : t = NpicSEDT_Saito_rev  (np1); break;
            case  E_DT *10+ E_QUADRA : t = NpicSEDT_quadra     (np1); break;
            case E_IDT *10+ E_QUADRA : t = NpicSEDT_quadra_inf (np1); break;
            case E_RDT *10+ E_QUADRA : t = NpicSEDT_quadra_rev (np1); break;
        }
        printf ("time: %.6f s\n", t);

    } else {

        int i;
        Npic_image *np2 = NpicDupImage (np1),
                   *np3 = NpicDupImage (np1);
        NpicCopyImageI (np2, np1);
        NpicCopyImageI (np3, np1);

        printf ("Computing %s %-12s ... ", what_s[what_e], how_s[E_HIRATA]); fflush (stdout);
        switch (what_e) {
            case E_DT  : t = NpicSEDT_Hirata     (np1); break;
            case E_IDT : t = NpicSEDT_Hirata_inf (np1); break;
            case E_RDT : t = NpicSEDT_Hirata_rev (np1); break;
        }
        printf ("time: %.6f s\n", t);

        printf ("Computing %s %-12s ... ", what_s[what_e], how_s[E_SAITO]); fflush (stdout);
        switch (what_e) {
            case E_DT  : t = NpicSEDT_Saito     (np2); break;
            case E_IDT : t = NpicSEDT_Saito_inf (np2); break;
            case E_RDT : t = NpicSEDT_Saito_rev (np2); break;
        }
        printf ("time: %.6f s\n", t);

        printf ("Comparing Hirata and Saito ...\n");
        i = NpicCompare (np1, np2, 10);
        if (i != 0) printf ("### %d pixels differ\n", i);

        printf ("Computing %s %-12s ... ", what_s[what_e], how_s[E_QUADRA]); fflush (stdout);
        switch (what_e) {
            case E_DT  : t = NpicSEDT_quadra     (np3); break;
            case E_IDT : t = NpicSEDT_quadra_inf (np3); break;
            case E_RDT : t = NpicSEDT_quadra_rev (np3); break;
        }
        printf ("time: %.6f s\n", t);

        printf ("Comparing Hirata and quadratic ...\n");
        i = NpicCompare (np1, np3, 10);
        if (i != 0) printf ("### %d pixels differ\n", i);

        printf ("Comparing Saito and quadratic ...\n");
        i = NpicCompare (np2, np3, 10);
        if (i != 0) printf ("### %d pixels differ\n", i);

        NpicDestroyImage (np2);
        NpicDestroyImage (np3);
    }

    printf ("Saving \"%s\"\n", out1);
    if (NpicWriteImage (np1, out1) != NPIC_SUCCESS) exit(1);
    
    NpicDestroyImage (np1);
    exit (0);
}