#define _GNU_SOURCE
#include <math.h>
#include <geraYuvGrad.h>
#include <estimaDistancia_b.h>

#include <sample_conv.h>
#include <float_image_gradient_2.h>
#include <float_image_gradient.h>
#include <operacoes_b.h>
#include <string.h>
#include <ia.h>
#include <reducaoRGB.h>

typedef float_image_t * float_image_tp;

//PROTOTIPOS

void BuscaMusis(float_image_tp imgs[],int res_max,int cumul,QDIST_t qualDist,int usa_ia,int usa_md_sd,float base);
void leituraImagens(char *nomeA,char *nomeB,char *nomeC,float_image_tp imgs[]);
void liberaListaImagens(float_image_tp imgs[],int res_max);
void criaHiMdSd(float_image_t *Md,float_image_t *Sd);
void geraPiramide(float_image_t *img,float_image_tp p1[],float_image_tp p2[],float_image_tp p3[],float_image_tp p4[],float_image_tp p5[],float_image_tp p6[],int res_max,int usa_ia,int usa_md_sd);


//IMPLEMENTACOES

int main(int argc, char** argv)
{
  char *nomeA = argv[1];
  char *nomeB = argv[2];
  char *nomeC = argv[3];
  int tipo=atoi(argv[4]);
  int res_max = atoi(argv[5]);
  int usa_ia=atoi(argv[6]);
  int usa_md_sd=atoi(argv[7]);
  int musis=atoi(argv[8]);
  QDIST_t qualDist = atoi(argv[9]);
  int cumul = atoi(argv[10]);
  float base = atof(argv[11]);

  float_image_tp imgs[3];
  leituraImagens(nomeA,nomeB,nomeC,imgs);

  int i;
  if (tipo==1)
  {
	for(i=0;i<3;i++)
		imgs[i]=criaImagem_Gradiente(imgs[i]);
  }
  else if (tipo==2)
  {
	for(i=0;i<3;i++)
		imgs[i]=criaImagemYUV(imgs[i]);
  }

  if (musis)
  {
      if (cumul)
      {
	if (qualDist==QDIST_EURGB)
        {
          if (tipo==1)
            fprintf(stderr,"$k$ & $\\fdist(A,B)$ & $\\fdist(A,C)$ & $\\lambda_k$ \\\\ \\hline \n");
          else
	    fprintf(stderr,"$k$ & $\\adist(A,B)$ & $\\adist(A,C) & $\\lambda_k$ \\\\ \\hline \n");         
        }
        else if (qualDist==QDIST_EUGRADYUV)
          fprintf(stderr,"$k$ & $\\fdist(\\h{A},\\h{B})$ & $\\fdist(\\h{A},\\h{C})$ & $\\lambda_k$ \\\\ \\hline \n");
      }
      else
      {
	 fprintf(stderr,"$k$ & $\\edist^{(k)}(A,B)$ & $\\edist^{(k)}$ \\\\ \n");
      }
      BuscaMusis(imgs,res_max,cumul,qualDist,usa_ia,usa_md_sd,base);
  }
  else
  {
        fprintf(stderr,"$k$ & $\\edist(A,B)$ & $\\edist(A,C)$ \\\\ \n");
	int ndistExP;
	float ab=calcula_dist_Monoescala(imgs[0],imgs[1],0,qualDist,&ndistExP);
	float ac=calcula_dist_Monoescala(imgs[0],imgs[2],0,qualDist,&ndistExP);
	fprintf(stderr,"0 & %8.6f & %8.6f \\\\ \n",ab,ac);
  }

  liberaListaImagens(imgs,3);	

  return 0;
}


void BuscaMusis(float_image_tp* imgs,int res_max,int cumul,QDIST_t qualDist,int usa_ia,int usa_md_sd,float base)
{
  float_image_tp aMdLo[res_max],aMdHi[res_max],aSdLo[res_max],aSdHi[res_max],aLo[res_max],aHi[res_max];
  float_image_tp bMdLo[res_max],bMdHi[res_max],bSdLo[res_max],bSdHi[res_max],bLo[res_max],bHi[res_max];
  float_image_tp cMdLo[res_max],cMdHi[res_max],cSdLo[res_max],cSdHi[res_max],cLo[res_max],cHi[res_max];

  geraPiramide(imgs[0],aMdLo,aSdLo,aMdHi,aSdHi,aLo,aHi,res_max,usa_ia,usa_md_sd);
  geraPiramide(imgs[1],bMdLo,bSdLo,bMdHi,bSdHi,bLo,bHi,res_max,usa_ia,usa_md_sd);
  geraPiramide(imgs[2],cMdLo,cSdLo,cMdHi,cSdHi,cLo,cHi,res_max,usa_ia,usa_md_sd);

  int res,i;
  float distab[3],distac[3];
  float distAcc=0,distBcc=0;
  Interval distEstIA = (Interval){ 0.0, 1.0 };
  Interval distEstIB = (Interval){ 0.0, 1.0 };

    double lambda[res_max+1];
    for (i=0;i<=res_max;i++)
    {
      lambda[i]=calculaLambda(cumul,i,res_max,base);
      //fprintf(stderr, "Lambda[%d] = %10.7f\n", i, lambda[i]);
    }        

  int ndistExP=0,ndistIAP=0;

  for(res=res_max;res>=0;res--)
  {

     float_image_t *ALo=aLo[res];
     float_image_t *AHi=aHi[res];
     float_image_t *AMdLo=aMdLo[res];
     float_image_t *AMdHi=aMdHi[res];
     float_image_t *ASdLo=aSdLo[res];
     float_image_t *ASdHi=aSdHi[res];

     float_image_t *BLo=bLo[res];
     float_image_t *BHi=bHi[res];
     float_image_t *BMdLo=bMdLo[res];
     float_image_t *BMdHi=bMdHi[res];
     float_image_t *BSdLo=bSdLo[res];
     float_image_t *BSdHi=bSdHi[res];

     float_image_t *CLo=cLo[res];
     float_image_t *CHi=cHi[res];
     float_image_t *CMdLo=cMdLo[res];
     float_image_t *CMdHi=cMdHi[res];
     float_image_t *CSdLo=cSdLo[res];
     float_image_t *CSdHi=cSdHi[res];

     if (cumul)
     {
       estima_dist_Multiescala
        ( ALo, AHi, AMdLo, AMdHi, ASdLo, ASdHi,
          BLo, BHi, BMdLo, BMdHi, BSdLo, BSdHi,
          res, res_max, 0, usa_ia, usa_md_sd, 
          qualDist, &distAcc, lambda, distab, &distEstIA,
          &ndistExP, &ndistIAP
        );
       
       estima_dist_Multiescala
        ( ALo, AHi, AMdLo, AMdHi, ASdLo, ASdHi,
          CLo, CHi, CMdLo, CMdHi, CSdLo, CSdHi,
          res, res_max, 0, usa_ia, usa_md_sd, 
          qualDist, &distBcc, lambda, distac, &distEstIB,
          &ndistExP, &ndistIAP
        );   
     }
     else
     {		
	estima_dist_Monoescala(ALo,AHi,AMdLo,AMdHi,ASdLo,ASdHi,BLo,BHi,BMdLo,BMdHi,BSdLo,BSdHi,res,res_max,0,usa_ia,usa_md_sd,qualDist,distab,&ndistExP,&ndistIAP);

	estima_dist_Monoescala(ALo,AHi,AMdLo,AMdHi,ASdLo,ASdHi,CLo,CHi,CMdLo,CMdHi,CSdLo,CSdHi,res,res_max,0,usa_ia,usa_md_sd,qualDist,distac,&ndistExP,&ndistIAP);
     }
     fprintf(stderr,"%d & %8.6f & %8.6f & %8.6f \\\\ \\hline \n",res,distab[0],distac[0],lambda[res]);
  }
  
  liberaListaImagens(aLo,res_max);
  liberaListaImagens(bLo,res_max);
  liberaListaImagens(cLo,res_max);
  liberaListaImagens(aHi,res_max);
  liberaListaImagens(bHi,res_max);
  liberaListaImagens(cHi,res_max);
  liberaListaImagens(aMdLo,res_max);
  liberaListaImagens(bMdLo,res_max);
  liberaListaImagens(cMdLo,res_max);
  liberaListaImagens(aMdHi,res_max);
  liberaListaImagens(bMdHi,res_max);
  liberaListaImagens(cMdHi,res_max);
  liberaListaImagens(aSdLo,res_max);
  liberaListaImagens(bSdLo,res_max);
  liberaListaImagens(cSdLo,res_max);
  liberaListaImagens(aSdHi,res_max);
  liberaListaImagens(bSdHi,res_max);
  liberaListaImagens(cSdHi,res_max);
}


void leituraImagens(char *nomeA,char *nomeB,char *nomeC,float_image_tp imgs[])
{
	float_image_t *A = getImageOpenRGB(nomeA,2);
	float_image_t *B = getImageOpenRGB(nomeB,2);
	float_image_t *C = getImageOpenRGB(nomeC,2);
	imgs[0]=A;
	imgs[1]=B;
	imgs[2]=C;
}

void liberaListaImagens(float_image_tp imgs[],int res_max)
{
	int i;
	float_image_t *img=NULL;
	for (i=0;i<res_max;i++)
	{
		img=imgs[i];
		float_image_free(img);
		img=NULL;
	}
}

void criaHiMdSd(float_image_t *Md,float_image_t *Sd)
{
	int i,j,c;
	float valormd=0,valorsd=0;
	for (i=0;i<Md->sz[1];i++)
		for (j=0;j<Md->sz[2];j++)
			for (c=0;c<Md->sz[0];c++)
			{
				valormd=float_image_get_sample(Md, c, i, j);
				valorsd=float_image_get_sample(Sd, c, i, j);
				valormd+=1/255;
				valorsd+=1/255;
				float_image_set_sample(Md,c,i,j,valormd);
				float_image_set_sample(Sd,c,i,j,valorsd);
			}
}



void geraPiramide(float_image_t *img,float_image_tp p1[],float_image_tp p2[],float_image_tp p3[],float_image_tp p4[],float_image_tp p5[],float_image_tp p6[],int res_max,int usa_ia,int usa_md_sd)
{
  int res=0;

 float_image_t *Md=float_image_copy(img);
 float_image_t *Vr=float_image_copy(img);
 float_image_fill(Vr, 0.0);
 float_image_t *Lo=float_image_copy(img);
 float_image_t *Hi=float_image_copy(img);
 
  while (1)
    {
	p1[res]=float_image_copy(Md);
	float_image_t *Sd = calculaDesvioPadrao(Vr);
	p2[res]=float_image_copy(Sd);
	float_image_t *MdHi=float_image_copy(Md);
	float_image_t *SdHi=float_image_copy(Sd);	   
	criaHiMdSd(MdHi,SdHi);
	p3[res]=MdHi;
	p4[res]=SdHi;
	float_image_free(Sd);
	
	p5[res]=float_image_copy(Lo);
	p6[res]=float_image_copy(Hi);     
	res++;
        if (res>res_max) break;        
        calculaMediaVarianciaRGB(Md,Vr,2);
	calculaIntervaloLoRGB(Lo,2,0.0,0);
        calculaIntervaloHiRGB(Hi,2,0.0,0);
    }
  float_image_free(Md);
  float_image_free(Vr);
  float_image_free(Lo);
  float_image_free(Hi);
}