/* Change Detection Defference Statistics v1.0.0 Ver.:OK (c) 12/07 - Christos Iosifides */
/*
   Name:   CD_DiffStat_8ui
   Version:  1.0.0
   Author:  Christos Iosifidis
   Date:   03/12/07 11:30
   Description: Change Detection Differennce and Statistics on 8bit unsigned int BIL datasets
   Usage: >CD_DiffStat_8ui 1st 2nd_filename_of_the_8bit_unsigned_int_BIL_dataset Signed16bitInt_BIL_dataset rows columns bands Statistics_Report_file

   This program is free software: you can redistribute it and/or modify
   it under the terms of the GNU General Public License as published by
   the Free Software Foundation, either version 3 of the License, or
   (at your option) any later version.

   This program is distributed in the hope that it will be useful,
   but WITHOUT ANY WARRANTY; without even the implied warranty of
   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
   GNU General Public License for more details.

   You should have received a copy of the GNU General Public License
   along with this program.  If not, see <http://www.gnu.org/licenses/>.
 */

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

void error(const char *);

int main(int argc, char *argv[]) {
    FILE *fpin1, *fpin2, *fpout;
    int rows, cols, bands, X1, X2, *stat_a1, s;
    unsigned char *in_buffer1, *in_buffer2;
    signed short int *out_buffer;
    register int r, c, i, j, k;

    if (argc!=7)
        error("Usage parameters: 1st 2nd_filename_of_the_8bit_unsigned_int_BIL_dataset Signed16bitInt_BIL_dataset rows columns bands");
    rows=atoi(argv[4]);
    cols=atoi(argv[5]);
    bands=atoi(argv[6]);
    if (rows<3)
        error("rows must be larger than 3");
    if (cols<3)
        error("columns must be larger than 3");
    if (bands<1)
        error("bands must be more than 1");

    /* Allocate statistics array */
    if ((stat_a1=(int *)malloc(sizeof(int)*256*256*bands))==NULL)
        error("Not enough memory avaliable for statistics array");
    /* Initialize statistics array - calloc can be used instead */
    for (i=0;i<256;i++)
        for (j=0;j<256;j++)
            for (k=0;k<bands;k++) {
                stat_a1[ (i*256+j)*bands+k ]=0;
            }

    /* Allocate 1st input buffer */
    if ((in_buffer1=(unsigned char *)malloc(cols))==NULL)
        error("Not enough memory avaliable for 1st input buffer");
    /* Allocate 2nd input buffer */
    if ((in_buffer2=(unsigned char *)malloc(cols))==NULL)
        error("Not enough memory avaliable for 2nd input buffer");
    /* Open 1st input image */
    if ((fpin1=fopen(argv[1],"rb"))==NULL)
        error("Bad 1st input file name");
    /* Open 2nd input image */
    if ((fpin2=fopen(argv[2],"rb"))==NULL)
        error("Bad 2nd input file name");

    /* Allocate output file buffer */
    if ((out_buffer=(signed short int *)malloc(cols*sizeof(signed short int)))==NULL)
        error("Not enough memory avaliable for output buffer");
    /* Open output image */
    if ((fpout=fopen(argv[3],"wb"))==NULL)
        error("Bad output file open");

    /* For every row in [ 0, rows ) */
    for (r=0;r<rows;r++)
        /* For every band in [ 0, bands ) */
        for (k=0;k<bands;k++) {
            /* Read line in buffers */
            if (fread(in_buffer1, 1, cols, fpin1)!=cols)
                error("Bad read from 1st input file");
            if (fread(in_buffer2, 1, cols, fpin2)!=cols)
                error("Bad read from 2nd input file");
            /* For every column in [ 0, cols ) */
            for (c=0;c<cols;c++) {
                X1=(int)in_buffer1[c];
                X2=(int)in_buffer2[c];
                stat_a1[(X1*256+X2)*bands+k]++;
                out_buffer[c]=(signed short int)( X1-X2 );
            }
            /* Write results */
            if (fwrite(out_buffer, sizeof(signed short int), cols, fpout)!=cols)
                error("Bad write to output file");
        }

    /* Close files and ... */
    fclose(fpin1);
    fclose(fpin2);
    fclose(fpout);
    free(in_buffer1);
    free(in_buffer2);
    free(out_buffer);

    /* Statistics Report */
    printf("\nCD_DiffStat_8UI reports:\n");
    printf("First input dataset [A]:\t\t%s\n",argv[1]);
    printf("Second input dataset [B]:\t\t%s\n",argv[2]);
    printf("First-Second difference dataset:\t%s\n",argv[3]);
    printf("Rows:\t\t%8d\n",rows);
    printf("Columns:\t%8d\n",cols);
    printf("Bands:\t\t%8d\n",bands);
    printf("\nStatistics:\n");
    /* For every band in [ 0, bands ) */
    for (k=0;k<bands;k++) {
        printf("\nBand:\t%8d\n",k+1);
        printf("_B:__A:_");
        for (i=0;i<256;i++)
            printf("%8d",i+1);
        for (i=0;i<256;i++) {
            printf("\n%7d:",i+1);
            s=0;
            for (j=0;j<256;j++) {
                printf("%8d",stat_a1[(i*256+j)*bands+k]);
                s+=stat_a1[(i*256+j)*bands+k];
            }
            printf("%8d",s);
        }
        for (i=0;i<256;i++) {
            printf("\n________");
            s=0;
            for (j=0;j<256;j++)
                s+=stat_a1[(j*256+i)*bands+k];
            printf("%8d",s);
        }
        s=0;
        for (i=0;i<256;i++)
            s+=stat_a1[(i*256+i)*bands+k];
        printf("%8d",s);
        s=0;
        for (i=0;i<256;i++)
            for (j=0;j<256;j++)
                s+=stat_a1[(i*256+j)*bands+k];
        printf("%8d\n",s);
    }
    free(stat_a1);

    return 0;
}

void error(const char *s) {
    printf("\nCD_DiffStat_8UI reports: Error: <%s>.\n",s);
    exit(1);
}