y4munsharp.c

Motion JPEG编解码器源代码

			interlaced ? uvheight / 2 : uvheight);
		blur_line(ctable_uv, cmatrix_uv, cmatrix_uv_len,
			cur_col,
			dest_col,
			interlaced ? uvheight / 2 : uvheight);
		put_column(dest_col, &o_yuv[2][col],
			interlaced ? 2 * uvwidth : uvwidth,
			interlaced ? uvheight / 2 : uvheight);
		if	(interlaced)
			{
			get_column(&o_yuv[2][col + uvwidth], cur_col,
				2 * uvwidth,
				uvheight / 2);
			blur_line(ctable_uv, cmatrix_uv, cmatrix_uv_len,
				cur_col,
				dest_col,
				interlaced ? uvheight / 2 : uvheight);
			put_column(dest_col, &o_yuv[2][col + uvwidth],
				2 * uvwidth,
				uvheight / 2);
			}
		}
merging:
	mjpeg_log(LOG_DEBUG, "Merging luma frame %d", frameno);
	for	(row = 0, i_ptr = i_yuv[0], o_ptr = o_yuv[0]; row < yheight; row++)
		{
		for	(i = 0; i < ywidth; i++, i_ptr++, o_ptr++)
			{
			diff = *i_ptr - *o_ptr;
			if	(abs(2 * diff) < y_threshold)
				diff = 0;
			value = *i_ptr + (y_amount * diff);
/*
 * For video the limits are 16 and 235 for the luma rather than 0 and 255!
*/
			if	(value < lowy)
				value = lowy;
			else if	(value > highy)
				value = highy;
			*o_ptr = value;
			}
		}

	if	(uv_radius == -1.0)
		goto done;

	mjpeg_log(LOG_DEBUG, "Merging chroma frame %d", frameno);
	for	(row = 0, i_ptr = i_yuv[1], o_ptr = o_yuv[1]; row < uvheight; row++)
		{
		for	(i = 0; i < uvwidth; i++, i_ptr++, o_ptr++)
			{
			diff = *i_ptr - *o_ptr;
			if	(abs(2 * diff) < uv_threshold)
				diff = 0;
			value = *i_ptr + (uv_amount * diff);
/*
 * For video the limits are 16 and 240 for the chroma rather than 0 and 255!
*/
			if	(value < lowuv)
				value = lowuv;
			else if	(value > highuv)
				value = highuv;
			*o_ptr = value;
			}
		}
	for	(row = 0, i_ptr = i_yuv[2], o_ptr = o_yuv[2]; row < uvheight; row++)
		{
		for	(i = 0; i < uvwidth; i++, i_ptr++, o_ptr++)
			{
			diff = *i_ptr - *o_ptr;
			if	(abs(2 * diff) < uv_threshold)
				diff = 0;
			value = *i_ptr + (uv_amount * diff);
/*
 * For video the limits are 16 and 240 for the chroma rather than 0 and 255!
*/
			if	(value < 16)
				value = 16;
			else if	(value > highuv)
				value = highuv;
			*o_ptr = value;
			}
		}
done:
	return;
	}

void
get_column(u_char *in, u_char *out, int stride, int numrows)
	{
	int	i;

	for	(i = 0; i < numrows; i++)
		{
		*out++ = *in;
		in += stride;
		}
	}

void
put_column(u_char *in, u_char *out, int stride, int numrows)
	{
	int	i;

	for	(i = 0; i < numrows; i++)
		{
		*out = *in++;
		out += stride;
		}
	}

/* 
 * The blur_line(), gen_convolve_matrix() and gen_lookup_table() functions 
 * were lifted almost intact from the GIMP unsharp plugin.  malloc was used 
 * instead of g_new() and the style was cleaned up a little but the logic 
 * was left untouched.
*/

/* this function is written as if it is blurring a column at a time,
 * even though it can operate on rows, too.  There is no difference
 * in the processing of the lines, at least to the blur_line function.
*/
static void
blur_line (double *ctable, double *cmatrix, int cmatrix_length,
	   u_char  *cur_col, u_char  *dest_col, int y)
  {
  double scale, sum, *cmatrix_p, *ctable_p;
  int i=0, j=0, row, cmatrix_middle = cmatrix_length/2;
  u_char  *cur_col_p, *cur_col_p1, *dest_col_p;

  /* this first block is the same as the non-optimized version --
   * it is only used for very small pictures, so speed isn't a
   * big concern.
   */
  if (cmatrix_length > y)
     {
     for (row = 0; row < y ; row++)
	 {
	 scale=0;
	  /* find the scale factor */
	 for  (j = 0; j < y ; j++)
	      {
	      /* if the index is in bounds, add it to the scale counter */
	      if  ((j + cmatrix_length/2 - row >= 0) &&
		   (j + cmatrix_length/2 - row < cmatrix_length))
		  scale += cmatrix[j + cmatrix_length/2 - row];
	      }
	  for  (i = 0; i< 1; i++)
	       {
	       sum = 0;
	       for  (j = 0; j < y; j++)
                    {
		    if  ((j >= row - cmatrix_length/2) &&
		         (j <= row + cmatrix_length/2))
		    	sum += cur_col[j + i] * cmatrix[j];
		    }
	        dest_col[row + i] = (u_char) ROUND (sum / scale);
	        }
	   }
      }
  else
      {
      /* for the edge condition, we only use available info and scale to one */
      for (row = 0; row < cmatrix_middle; row++)
	  {
	  /* find scale factor */
	  scale=0;
	  for  (j = cmatrix_middle - row; j<cmatrix_length; j++)
	       scale += cmatrix[j];
	  for  (i = 0; i < 1; i++)
	       {
	       sum = 0;
	       for  (j = cmatrix_middle - row; j<cmatrix_length; j++)
		    sum += cur_col[(row + j-cmatrix_middle) + i] * cmatrix[j];
	        dest_col[row + i] = (u_char) ROUND (sum / scale);
	       }
	   }
      /* go through each pixel in each col */
      dest_col_p = dest_col + row;
      for  (; row < y-cmatrix_middle; row++)
	   {
	   cur_col_p = (row - cmatrix_middle) + cur_col;
	   for  (i = 0; i < 1; i++)
	        {
	        sum = 0;
	        cmatrix_p = cmatrix;
	        cur_col_p1 = cur_col_p;
	        ctable_p = ctable;
	        for  (j = cmatrix_length; j>0; j--)
		     {
		     sum += *(ctable_p + *cur_col_p1);
		     cur_col_p1 += 1;
		     ctable_p += 256;
		     }
	        cur_col_p++;
	        *(dest_col_p++) = ROUND (sum);
	        }
	    }
	
      /* for the edge condition, we only use available info, and scale to one */
       for  (; row < y; row++)
	    {
	    /* find scale factor */
	    scale=0;
	    for  (j = 0; j< y-row + cmatrix_middle; j++)
	         scale += cmatrix[j];
	    for  (i = 0; i < 1; i++)
	         {
	         sum = 0;
	         for  (j = 0; j<y-row + cmatrix_middle; j++)
		      sum += cur_col[(row + j-cmatrix_middle) + i] * cmatrix[j];
	         dest_col[row + i] = (u_char) ROUND (sum / scale);
	         }
	    }
      }
  }

/*
 * generates a 1-D convolution matrix to be used for each pass of 
 * a two-pass gaussian blur.  Returns the length of the matrix.
*/
static int
gen_convolve_matrix(double radius, double **cmatrix_p)
	{
	int matrix_length, matrix_midpoint, i, j;
	double *cmatrix, std_dev, sum, base_x;
	
  /* we want to generate a matrix that goes out a certain radius
   * from the center, so we have to go out ceil(rad-0.5) pixels,
   * inlcuding the center pixel.  Of course, that's only in one direction,
   * so we have to go the same amount in the other direction, but not count
   * the center pixel again.  So we double the previous result and subtract
   * one.
   * The radius parameter that is passed to this function is used as
   * the standard deviation, and the radius of effect is the
   * standard deviation * 2.  It's a little confusing.
   */
	radius = fabs(radius) + 1.0;
	
	std_dev = radius;
	radius = std_dev * 2;

	/* go out 'radius' in each direction */
	matrix_length = 2 * ceil(radius-0.5) + 1;
	if (matrix_length <= 0) matrix_length = 1;
	matrix_midpoint = matrix_length/2 + 1;
	*cmatrix_p = (double *)malloc(sizeof (double) * matrix_length);
	cmatrix = *cmatrix_p;

  /*  Now we fill the matrix by doing a numeric integration approximation
   * from -2*std_dev to 2*std_dev, sampling 50 points per pixel.
   * We do the bottom half, mirror it to the top half, then compute the
   * center point.  Otherwise asymmetric quantization errors will occur.
   *  The formula to integrate is e^-(x^2/2s^2).
   */

  /* first we do the top (right) half of matrix */
	for	(i = matrix_length/2 + 1; i < matrix_length; i++)
		{
		base_x = i - floor(matrix_length/2) - 0.5;
		sum = 0;
		for	(j = 1; j <= 50; j++)
			{
			if	(base_x+0.02*j <= radius)
				sum += exp (-(base_x+0.02*j)*(base_x+0.02*j) / 
					(2*std_dev*std_dev));
			}
		cmatrix[i] = sum/50;
		}

  /* mirror the thing to the bottom half */
	for	(i=0; i<=matrix_length/2; i++)
		cmatrix[i] = cmatrix[matrix_length-1-i];
	
  /* find center val -- calculate an odd number of quanta to make it symmetric,
   * even if the center point is weighted slightly higher than others. */
	sum = 0;
	for	(j=0; j<=50; j++)
  		sum += exp (-(0.5+0.02*j)*(0.5+0.02*j) / (2*std_dev*std_dev));
	cmatrix[matrix_length/2] = sum/51;
	
  /* normalize the distribution by scaling the total sum to one */
	sum=0;
	for	(i = 0; i < matrix_length; i++)
		sum += cmatrix[i];
	for	(i=0; i<matrix_length; i++)
		cmatrix[i] = cmatrix[i] / sum;
	return(matrix_length);
	}

/* generates a lookup table for every possible product of 0-255 and
   each value in the convolution matrix.  The returned array is
   indexed first by matrix position, then by input multiplicand (?) value.
*/
static double *
gen_lookup_table(double *cmatrix, int cmatrix_length)
	{
	int i, j;
	double* lookup_table = (double *)malloc(sizeof (double) * cmatrix_length * 256);
	double* lookup_table_p = lookup_table, *cmatrix_p = cmatrix;

	for	(i=0; i<cmatrix_length; i++)
		{
		for	(j=0; j<256; j++)
	  		*(lookup_table_p++) = *cmatrix_p * (double)j;
		cmatrix_p++;
		}
	return(lookup_table);
	}

void usage(void)
	{
	fprintf(stderr, "%s: usage: [-v 0|1|2] [-N] [-L radius,amount,threshold] [-C radius,amount,threshold]\n", __progname);
	fprintf(stderr, "%s:\tradius and amount are floating point numbers\n",
		__progname);
	fprintf(stderr, "%s:\tthreshold is integer.\n", __progname);
	fprintf(stderr, "%s:\tdefault for -L is 3.0,0.3,4\n", __progname);
	fprintf(stderr, "%s:\tchroma not filtered UNLESS -C used, no default\n",
		__progname);
	fprintf(stderr, "%s:-v verbose 0=quiet 1=normal 2=debug (default: 1)\n", __progname);
	fprintf(stderr, "%s:-N disables the Y' 16-235 clip/core\n", __progname);
	fprintf(stderr, "%s:   disables the CbCr 16-240 clip/core\n", __progname);
	fprintf(stderr, "%s:   (allow the entire 0 to 255 range).\n", __progname);
	exit(1);
	}