encoderquad.cpp

LDPC码的实现,包括编码器和解码器,使用了DCT.

// encoderQuad.cpp

//------------------------------------------------------------
// David Chen*, David Varodayan, Markus Flierl, Bernd Girod
// Image, Video, and Multimedia Systems Group
// Information Systems Laboratory
// Stanford University
//
// *Contact: [e-mail]
//------------------------------------------------------------

#include "encoderQuad.h"
#include <ChenImageDCT.h>

int pixelQuads_encodeFrameYUV(const short** pFrame, const char* pLadderFile, const short* pScaledQArray, uchar** pAccumSyndrome)
{
	bool bSuccess = true;

	// Encode Y component using quadrant encoder
	RECORD_SUCCESS( pixelQuads_encodeImage(pFrame[0], pLadderFile, pScaledQArray, pAccumSyndrome), bSuccess );

	// Encode U,V components using regular encoder
	RECORD_SUCCESS( pixelQuads_encodeImageQuad(pFrame[1], pLadderFile, pScaledQArray, pAccumSyndrome[4]), bSuccess );
	RECORD_SUCCESS( pixelQuads_encodeImageQuad(pFrame[2], pLadderFile, pScaledQArray, pAccumSyndrome[5]), bSuccess );

	return bSuccess ? GOOD_RETURN : BAD_RETURN;
}

int pixelQuads_encodeImage(const short* pImgX, const char* pLadderFile, const short* pScaledQArray, uchar** pAccumSyndromeQuads)
{
	bool bSuccess = true;

	// Break image into four quadrants
	short* pImgQuads[4];
	for (int nQuad = 0; nQuad < 4; nQuad++) {
		pImgQuads[nQuad] = new short[PIXELS_QUAD];
	}
	RECORD_SUCCESS( ChenImage_imageToQuadrants(WIDTH, HEIGHT, pImgX, pImgQuads), bSuccess );

	// Encode each quadrant
	for (int nQuad = 0; nQuad < 4; nQuad++) {
		RECORD_SUCCESS( pixelQuads_encodeImageQuad(pImgQuads[nQuad], pLadderFile, pScaledQArray, pAccumSyndromeQuads[nQuad]), bSuccess );
	}
	
	// Clean up
	for (int nQuad = 0; nQuad < 4; nQuad++) {
		delete [] pImgQuads[nQuad];
	}
	return GOOD_RETURN;
}

int pixelQuads_encodeImageQuad(const short* pImgX, const char* pLadderFile, const short* pScaledQArray, uchar* pAccumSyndrome)
{
	bool bSuccess = true;

	// Peform blockwise-DCT and quantization
	short* pTransQuantX = new short[PIXELS_QUAD];
	RECORD_SUCCESS( ChenImageDCT_block88DCTQuantize(WIDTH_QUAD, HEIGHT_QUAD, pImgX, pScaledQArray, pTransQuantX), bSuccess );

	// Create level-shifted coefficients
	short* pLDPCSourcePix = new short[PIXELS_QUAD];
	for (int nPix = 0; nPix < PIXELS_QUAD; nPix++) {
		pLDPCSourcePix[nPix] = (pTransQuantX[nPix] + LEVELS) % LEVELS;
	}

	// Encode coefficients
	unsigned char* pSource = new unsigned char[BITS_QUAD];
	RECORD_SUCCESS( pixelQuads_encodeUnsignedCoeffQuad(pLDPCSourcePix, pLadderFile, pSource, pAccumSyndrome), bSuccess );

	delete [] pLDPCSourcePix;
	delete [] pSource;
	return bSuccess ? GOOD_RETURN : BAD_RETURN;
}

int pixelQuads_modelChannelFrameYUV(const short** pFrameX, const short** pFrameY, const short* pScaledQArray, float** pLambda, float* pEntropy, short** pShiftsX, short** pShiftsY)
{
	bool bSuccess = true;

	// Model channel for Y component
	RECORD_SUCCESS( pixelQuads_modelChannelImage(pFrameX[0], pFrameY[0], pScaledQArray, pLambda, pEntropy, pShiftsX, pShiftsY), bSuccess );

	// Model channel for U,V components
	RECORD_SUCCESS( pixelQuads_modelChannelQuad(pFrameX[1], pFrameY[1], pScaledQArray, pLambda[4], pEntropy+4, pShiftsX[4], pShiftsY[4]), bSuccess );
	RECORD_SUCCESS( pixelQuads_modelChannelQuad(pFrameX[2], pFrameY[2], pScaledQArray, pLambda[5], pEntropy+5, pShiftsX[5], pShiftsY[5]), bSuccess );

	return bSuccess ? GOOD_RETURN : BAD_RETURN;
}

int pixelQuads_modelChannelImage(const short* pImgX, const short* pImgY, const short* pScaledQArray, float** pLambda, float* pEntropy, short** pShiftsX, short** pShiftsY)
{
	bool bSuccess = true;

	// Break image into four quadrants
	short *pImgXQuads[4], *pImgYQuads[4];
	for (int nQuad = 0; nQuad < 4; nQuad++) {
		pImgXQuads[nQuad] = new short[PIXELS_QUAD];
		pImgYQuads[nQuad] = new short[PIXELS_QUAD];
	}
	RECORD_SUCCESS( ChenImage_imageToQuadrants(WIDTH, HEIGHT, pImgX, pImgXQuads), bSuccess );
	RECORD_SUCCESS( ChenImage_imageToQuadrants(WIDTH, HEIGHT, pImgY, pImgYQuads), bSuccess );

	// Model channel for each quadrant
	for (int nQuad = 0; nQuad < 4; nQuad++) {
		RECORD_SUCCESS( pixelQuads_modelChannelQuad(pImgXQuads[nQuad], pImgYQuads[nQuad], pScaledQArray, pLambda[nQuad], pEntropy+nQuad, pShiftsX[nQuad], pShiftsY[nQuad]), bSuccess );
	}

	// Clean up
	for (int nQuad = 0; nQuad < 4; nQuad++) {
		delete [] pImgXQuads[nQuad];
		delete [] pImgYQuads[nQuad];
	}
	return bSuccess ? GOOD_RETURN : BAD_RETURN;
}

int pixelQuads_modelChannelQuad(const short* pImgX, const short* pImgY, const short* pScaledQArray, float* pLambda, float* pEntropy, short* pShiftsX, short* pShiftsY)
{
	bool bSuccess = true;

	// DCT-transform and quantize X
	short* pTransQuantX = new short[PIXELS_QUAD];
	RECORD_SUCCESS( ChenImageDCT_block88DCTQuantize(WIDTH_QUAD, HEIGHT_QUAD, pImgX, pScaledQArray, pTransQuantX), bSuccess );

	// Find best shifted block matches in pixel domain
	short* pImgYComp = new short[PIXELS_QUAD];
	short* pImgRes   = new short[PIXELS_QUAD];
	RECORD_SUCCESS( ChenImage_imageBlockMatch2D(WIDTH_QUAD, HEIGHT_QUAD, BLOCK_SIZE, MAX_SHIFT, pImgX, pImgY, pImgYComp, pImgRes, pShiftsX, pShiftsY), bSuccess );

	// DCT-transform and quantize disparity-compensated Y
	short* pTransQuantYComp = new short[PIXELS_QUAD];
	RECORD_SUCCESS( ChenImageDCT_block88DCTQuantize(WIDTH_QUAD, HEIGHT_QUAD, pImgYComp, pScaledQArray, pTransQuantYComp), bSuccess );

	// Find Laplacian constants from transform residual
	short* pTransQuantResiduals = new short[PIXELS_QUAD];
	for (int nPix = 0; nPix < PIXELS_QUAD; nPix++) {
		pTransQuantResiduals[nPix] = pTransQuantYComp[nPix] - pTransQuantX[nPix];
	}
	RECORD_SUCCESS( pixelQuads_findLaplacianConstantTransformQuad(pTransQuantResiduals, BLOCK_SIZE, pLambda, pEntropy), bSuccess );
	
	delete [] pTransQuantX;
	delete [] pImgYComp;
	delete [] pImgRes;
	delete [] pTransQuantYComp;
	delete [] pTransQuantResiduals;
	return bSuccess ? GOOD_RETURN : BAD_RETURN;
}

int pixelQuads_encodeUnsignedCoeffQuad(const short* pImgX, const char* pLadderFile, unsigned char* pSource, unsigned char* pAccumSyndrome) 
{
	bool bSuccess = true;	

	// Converts pixels of image into bit sequence
	RECORD_SUCCESS( pixelQuads_unsignedCoefficientsToBitstreamQuad(pImgX, pSource), bSuccess );
	double* pSourceD = new double[BITS_QUAD];
	for (int nBit = 0; nBit < BITS_QUAD; nBit++) pSourceD[nBit] = pSource[nBit];

	// Encode using LDPC
	double* pAccumSyndromeD = new double[BITS_QUAD];
    RECORD_SUCCESS( pixelQuads_encodeBitsLDPCQuad(pSourceD, pLadderFile, pAccumSyndromeD), bSuccess );
	for (int nBit = 0; nBit < BITS_QUAD; nBit++) pAccumSyndrome[nBit] = (uchar)pAccumSyndromeD[nBit];

	delete [] pSourceD;
	delete [] pAccumSyndromeD;
	return bSuccess ? GOOD_RETURN : BAD_RETURN;
}

int pixelQuads_encodeBitsLDPCQuad(const double* pSource, const char* pLadderFile, double* pParity)
{
    
    int nCodeLength, nSubcodeLength, nCodeEdges;
    int nNumCodes, nTotalNumInc, nNumInc;
    
    FILE* pFile = fopen(pLadderFile, "r");
	if (pFile == NULL) return BAD_RETURN;
 
    fscanf(pFile, "%d", &nNumCodes);
    fscanf(pFile, "%d", &nCodeLength);
    fscanf(pFile, "%d", &nCodeEdges);
    fscanf(pFile, "%d", &nTotalNumInc);
        
    int* pIR = new int[nCodeEdges];
	int* pJC = new int[nCodeLength+1];
    int* pTxSeq = new int[nTotalNumInc];

	for (int k = 0; k < nCodeLength+1; k++) {
        fscanf(pFile, "%d", pJC+k);
	}
    for (int k = 0; k < nNumCodes; k++)
    {
        fscanf(pFile, "%d", &nNumInc);
		for (int l = 0; l < nNumInc; l++) {
            fscanf(pFile, "%d", pTxSeq+l);
		}
		for (int l = 0; l < nCodeEdges; l++) {
            fscanf(pFile, "%d", pIR+l);
		}
    }
    nSubcodeLength = (nCodeLength/nTotalNumInc)*nNumInc;
    
    memset(pParity, 0, nSubcodeLength*sizeof(double));
    
    // source * H'
	for (int k = 0; k < nCodeLength; k++) {
		for (int l = pJC[k]; l < pJC[k+1]; l++) {
            pParity[pIR[l]] += pSource[k];
		}
	}
    
    // accumulate
	for (int k = 1; k < nSubcodeLength; k++) {
        pParity[k] += pParity[k-1];
	}
    
    // mod 2
	for (int k = 0; k < nSubcodeLength; k++) {
        pParity[k] = (double) ((int) pParity[k] % 2);
	}

	delete [] pIR;
	delete [] pJC;
    delete [] pTxSeq;
	return GOOD_RETURN;
}

int pixelQuads_unsignedCoefficientsToBitstreamQuad(const short* pImg, uchar* pBitstream) {
	for (int nRow = 0; nRow < HEIGHT_QUAD; nRow++) {
		for (int nCol = 0; nCol < WIDTH_QUAD; nCol++) {
			for (int nBit = 0; nBit < BITPLANES; nBit++) {
				pBitstream[(nRow*WIDTH_QUAD + nCol)*BITPLANES + nBit] = (uchar)((pImg[nRow*WIDTH_QUAD + nCol] >> (BITPLANES-1-nBit)) % 2);
			}
		}
	}
	return GOOD_RETURN;
}

int pixelQuads_bitstreamToSignedCoefficientsQuad(const uchar* pBitstream, short* pImg) {
	short nPixel, nPart;
	for (int nRow = 0; nRow < HEIGHT_QUAD; nRow++) {
		for (int nCol = 0; nCol < WIDTH_QUAD; nCol++) {
			nPixel = 0;
			// Decode magnitude
			for (int nBit = 0; nBit < BITPLANES; nBit++) {
				nPart = (short) pBitstream[(nRow*WIDTH_QUAD + nCol)*BITPLANES + nBit];
				nPart = nPart << (BITPLANES-1-nBit);
				nPixel += nPart;
			}
			// Decode sign (assume AC coefficients <= 128)
			if (nRow % BLOCK_SIZE != 0 || nCol % BLOCK_SIZE != 0) {
				if (nPixel > 128) nPixel -= 256;
			}
			pImg[nRow*WIDTH_QUAD + nCol] = nPixel;
		}
	}
	return GOOD_RETURN;
}

int pixelQuads_findLaplacianConstantTransformQuad(const short* pTransResidual, const int nBlockSize, float* pLapParam, float* pEntropy)
{
	bool bSuccess = true;

	// Initialize
	int* pLevelsVec = new int[2*LEVELS-1];
	for (int nLevel = 0; nLevel < 2*LEVELS-1; nLevel++) {
		pLevelsVec[nLevel] = nLevel - LEVELS + 1;
	}
	int* pDiffHist = new int[2*LEVELS-1];
	float* pDiffPMF = new float[2*LEVELS-1];
	float* pQVec = new float[2*LEVELS-1];

	// Calculate PMF of residuals for each subband
	int nBlockCols = WIDTH_QUAD / nBlockSize;
	int nBlockRows = HEIGHT_QUAD / nBlockSize;
	int nBlocks = nBlockCols * nBlockRows;
	int nCoeffsPerBand = nBlocks;
	int nBands = nBlockSize * nBlockSize;
	short* pSubband = new short[nBlocks];
	float* pSubbandEntropy = new float[nBands];
	for (int nBandX = 0; nBandX < nBlockSize; nBandX++) 
	{
		for (int nBandY = 0; nBandY < nBlockSize; nBandY++) 
		{
			// Extract samples from subband
			for (int nBlockRow = 0; nBlockRow < nBlockRows; nBlockRow++) {
				int nRowReal = nBlockRow*nBlockSize + nBandY;
				for (int nBlockCol = 0; nBlockCol < nBlockCols; nBlockCol++) {
					int nColReal = nBlockCol*nBlockSize + nBandX;
					pSubband[nBlockRow*nBlockCols + nBlockCol] = pTransResidual[nRowReal*WIDTH_QUAD + nColReal];
				}
			}

			// Calculate histogram and PMF
			RECORD_SUCCESS( ChenImage_diffImagePMF(nBlockCols, nBlockRows, pSubband, pDiffPMF), bSuccess );

			// Find best Laplacian parameter to find residual statistics
			float fQSum, fEntropySum, fMinEntropySum, fBestLaplacian;
			fMinEntropySum = LARGE;
			for (float fLaplacianTry = (float)MIN_LAPLACIAN; fLaplacianTry <= (float)MAX_LAPLACIAN; fLaplacianTry += (float)0.01) {
				// Form Q vector
				fQSum = 0.0;
				for (int nLevel = 0; nLevel < 2*LEVELS-1; nLevel++) {
					pQVec[nLevel] = exp( -fLaplacianTry * abs(pLevelsVec[nLevel]) );
					fQSum += pQVec[nLevel];
				}
				for (int nLevel = 0; nLevel < 2*LEVELS-1; nLevel++) {
					pQVec[nLevel] /= fQSum;
				}

				// Calculate conditional entropy
				fEntropySum = 0.0;
				for (int nLevel = 0; nLevel < 2*LEVELS-1; nLevel++) {
					if (pDiffPMF[nLevel] > 0) {
						fEntropySum -= (float)( pDiffPMF[nLevel] * log(pQVec[nLevel])/log(2.0) / 8.0 );
					}
				}

				// Compare against best entropy
				if (fEntropySum < fMinEntropySum) {
					fMinEntropySum = fEntropySum;
					fBestLaplacian = fLaplacianTry;
				}
			}
			pLapParam[nBandY*nBlockSize + nBandX] = fBestLaplacian;
			pSubbandEntropy[nBandY*nBlockSize + nBandX] = fMinEntropySum;

		} // end nBandY
	} // end nBandX
	
	// Calculate entropy
	float fRateSum = 0;
	for (int nCoeff = 0; nCoeff < nBands; nCoeff++) {
		fRateSum += pSubbandEntropy[nCoeff];
	}
	*pEntropy = fRateSum / nBands;

	// Clean up
	delete [] pLevelsVec;
	delete [] pDiffHist;
	delete [] pQVec;
	delete [] pSubbandEntropy;
	return bSuccess ? GOOD_RETURN : BAD_RETURN;
}