#define _GNU_SOURCE
#include <stdlib.h>
#include <stdio.h>
#include <string.h>
#include <time.h>
#include <math.h>
#include <assert.h>
#include <jsfile.h>
#include <float_image.h>
#include <float_pnm_image_io.h>

#define IS_MAX 0
#define IS_MIN 1

void findMinMaxMean(float_image_t* img, double* pmin, double* pmax, double* pmean,double* prms,int c);
void plotHistogram(FILE* output_hist, FILE* output_inte, float_image_t* img, int num_steps,int c);

void findMinMaxMean(float_image_t* img, double* pmin, double* pmax, double* pmean,double* prms,int c){
	int NC,NX,NY;
	NC = img->sz[0];
	NX = img->sz[1];
	NY = img->sz[2];
	float total = NC*NX*NY;
	int x,y;
	
	double min = +INF;
	double max = -INF;
	double mean = 0;
	double rms = 0;
	//for(c = 0; c < NC; c++){
		for(x = 0; x < NX; x++){
			for(y = 0; y < NY; y++){
				float val = float_image_get_sample(img,c,x,y);
				if(min > val) min = val;
				if(max < val) max = val;
				mean+= (val/(double)total);
				rms+=(val*val)/(double)total;
			}
		}
	//}
	//mean = mean/(double)total;
	
	//rms =  sqrt(rms/(double)total);
	rms =  sqrt(rms);
	*pmean = mean;
	*pmax = max;
	*pmin = min;
	*prms = rms;

}

void plotHistogram(FILE* output_hist, FILE* output_inte, float_image_t* img, int num_steps,int c){
	double min,max,mean,rms;
	findMinMaxMean(img,&min,&max,&mean,&rms,c);
	double delta = max/(double)num_steps;
	int num_steps_real = (int) ceil(max/(double)delta);
	int* hist = (int*)malloc(sizeof(int)*num_steps_real);
	int i;
	for(i = 0; i < num_steps_real; i++){
		hist[i] = 0;
	}
	int NC,NX,NY;
	NC = img->sz[0];
	NX = img->sz[1];
	NY = img->sz[2];
	int x,y;
	double total = 0.0;
// 	for(c = 0; c < NC; c++){
		for(x = 0; x < NX; x++){
			for(y = 0; y < NY; y++){
				float val = float_image_get_sample(img,c,x,y);
				int slot = ((int)floor(val))/delta;
				if (slot < 0 ) slot = 0;
				if(slot >= num_steps_real) slot = num_steps_real - 1 ;
				hist[slot] = hist[slot]+ 1;
				total+=val;
			}
		}
// 	}
	if( output_hist != NULL){
		for(i = 0; i < num_steps_real; i++){
			double pos = (i+0.5)*delta;
			fprintf(output_hist,"%lf %d \n", pos, hist[i]);
			
		}
	}
	if( output_inte != NULL ){
		for(i = 0; i < num_steps_real; i++){
			double pos = (i+0.5)*delta;
			int val = 0;
			int j;
			for(j = num_steps_real - 1; j >= i ; j--){
				double posj = (j+0.5)*delta;
				val+= hist[j]*posj;
			}
			fprintf(output_inte,"%lf %lf \n", pos, val/total);
		}
	}
	
	free(hist);
	
}


int main(int argc, char** argv){
	if(argc < 2){
		fprintf(stderr,"Chamada do programa:\nprograma <imagem> {-verbose|-simple|-extense}\n");
		return 1;
	}
	
	char* infilename = argv[1];
	int mode = 0; //silent mode
	if(argc >=3){
	  if(strcmp(argv[2],"-verbose") == 0){
	    mode = 1;
	  }
	  if(strcmp(argv[2],"-simple") == 0){
	    mode = 2;
	  }
	  if(strcmp(argv[2],"-extense") == 0){
	    mode = 3;
	  }
	}
	FILE* imgfile = open_read(infilename,FALSE);
	float_image_t* img = float_image_read(imgfile);
	fclose(imgfile);
	int NC = img->sz[0];
	int NX = img->sz[1];
	int NY = img->sz[2];
	
	double min,max,mean,rms;
	//FILE* out_hist;
	//FILE* out_inte;
	//out_hist = open_write("hist.txt",FALSE);
	//out_inte = open_write("inte.txt",FALSE);
	//plotHistogram(out_hist,out_inte,img,1000);
	//fclose(out_hist);
	//fclose(out_inte);
	int c;
	for(c= 0; c < NC; c++){
	  findMinMaxMean(img,&min, &max,&mean,&rms,c);
	  double degmin = (180.0/M_PI)*min;
	  double degmax = (180.0/M_PI)*max;
	  double degmean = (180.0/M_PI)*mean;
	  double degrms = (180.0/M_PI)*rms;
	  if( mode == 1){
	    printf("MIN: %f ( %f deg ) MAX: %f ( %f deg ) MEAN: %f ( %f deg ) RMS: %f ( %f deg ) \n",min,degmin,max,degmax,mean,degmean,rms,degrms);
	  }else if( mode == 2){
	    printf("%d %d %d %9.6f %9.6f %9.6f %9.6f\n",NC,NX,NY,min,max,mean,rms);
	  }else if( mode == 3){
	    printf("%d %d %d %24.24f %24.24f %24.24f %24.24f\n",NC,NX,NY,min,max,mean,rms);
	  }else{
	    printf("%d %d %d %9.6e %9.6e %9.6e %9.6e\n",NC,NX,NY,min,max,mean,rms);
	  }
	}

	return 0;
}