/************************************************************************ *計算信噪比  ************************************************************************/#include"sim.h"#include<math.h>/********************************************************************** * *	Name:        SNRcomp *	Description:        Compares two image files using SNR *                      No conversion to 422 *	 *	Input:	       *	Returns:        *	Side effects: * *	Date: 930711	Author: Karl.Lillevold@nta.no * ***********************************************************************/void ComputeSNR(PictImage *im1, PictImage *im2, Results *res, int write){  FILE *out = NULL;  int n;  register int m;  int quad, quad_Cr, quad_Cb, diff;  PictImage *diff_image = NULL;  /* Diff. image written to diff_filename */  char *diff_filename=DEF_DIFFILENAME;  if (write) {    out = fopen(diff_filename,"ab");    diff_image = (PictImage *)malloc(sizeof(PictImage));    diff_image->lum = (unsigned char *)malloc(sizeof(char)*pels*lines);    diff_image->Cr =  (unsigned char *)malloc(sizeof(char)*pels*lines/4);    diff_image->Cb =  (unsigned char *)malloc(sizeof(char)*pels*lines/4);      }  quad = 0;  quad_Cr = quad_Cb = 0;  /* Luminance */  quad = 0;  for (n = 0; n < lines; n++)    for (m = 0; m < pels; m++) {      diff = *(im1->lum + m + n*pels) - *(im2->lum + m + n*pels);      if (write)        *(diff_image->lum + m + n*pels) = 10*diff + 128;      quad += diff * diff;    }  res->SNR_l = (float)quad/(float)(pels*lines);  if (res->SNR_l) {    res->SNR_l = (float)(255*255) / res->SNR_l;    res->SNR_l = 10 * (float)log10(res->SNR_l);  }  else res->SNR_l = (float)99.99;  /* Chrominance */  for (n = 0; n < lines/2; n++)    for (m = 0; m < pels/2; m++) {      quad_Cr += (*(im1->Cr+m + n*pels/2) - *(im2->Cr + m + n*pels/2)) *        (*(im1->Cr+m + n*pels/2) - *(im2->Cr + m + n*pels/2));      quad_Cb += (*(im1->Cb+m + n*pels/2) - *(im2->Cb + m + n*pels/2)) *        (*(im1->Cb+m + n*pels/2) - *(im2->Cb + m + n*pels/2));      if (write) {        *(diff_image->Cr + m + n*pels/2) =           (*(im1->Cr+m + n*pels/2) - *(im2->Cr + m + n*pels/2))*10+128;        *(diff_image->Cb + m + n*pels/2) =           (*(im1->Cb+m + n*pels/2) - *(im2->Cb + m + n*pels/2))*10+128;      }    }  res->SNR_Cr = (float)quad_Cr/(float)(pels*lines/4);  if (res->SNR_Cr) {    res->SNR_Cr = (float)(255*255) / res->SNR_Cr;    res->SNR_Cr = 10 * (float)log10(res->SNR_Cr);  }  else res->SNR_Cr = (float)99.99;  res->SNR_Cb = (float)quad_Cb/(float)(pels*lines/4);  if (res->SNR_Cb) {    res->SNR_Cb = (float)(255*255) / res->SNR_Cb;    res->SNR_Cb = 10 * (float)log10(res->SNR_Cb);  }  else res->SNR_Cb = (float)99.99;  if (write) {    fwrite(diff_image->lum, sizeof(char), pels*lines, out);    fwrite(diff_image->Cr,  sizeof(char), pels*lines/4, out);    fwrite(diff_image->Cb,  sizeof(char), pels*lines/4, out);    free(diff_image->lum);    free(diff_image->Cr);    free(diff_image->Cb);    free(diff_image);    fclose(out);  }  return;}