/*
 *	This is a prototype file, not meant to be compiled directly.
 *	The including source file must define the following:
 *
 *		TYPE_ARRAY -- the base (scalar) type of all argument arrays
 *		FUNC_1D -- the name of the 1-dimensional refinement function
 *		FUNC_ND -- the name of the N-dimensional refinement function
 */
//所有用于細化處理的小波基的尺度函數類型
static void wrefine_step _PROTO((waveletfilter *wfltr, TYPE_ARRAY *a, int n));


//一維離散小波細化函數
/*
函數輸入：
TYPE_ARRAY *a：		輸入的原始數據，用于細化的小波基；
int nA：		數據數據a的大小，這個數據的大小必須是2的指數值；
int nNew：		細化后小波基的大小，同樣必須是2的指數值；
waveletfilter *wfltr：	需要進行細化的小波濾波器；
TYPE_ARRAY *aRefined：	函數的實際輸出，細化后的數據。
*/
void FUNC_1D(a, nA, nNew, wfltr, aRefined)
	TYPE_ARRAY *a;	
	int nA;				
	int nNew;		
	waveletfilter *wfltr;	
	TYPE_ARRAY *aRefined;
{
	int iA;
	double scale;

	assert(nA <= nNew);
	for (iA = 0; iA < nA; iA++)
		aRefined[iA] = a[iA];
	for (iA = 2 * nA; iA <= nNew; iA *= 2)
		wrefine_step(wfltr, aRefined, iA);

	scale = sqrt(RATIO(nNew, nA));
	for (iA = 0; iA < nNew; iA++)
		aRefined[iA] *= scale;

	return;
}

#if 0 

//完成N維小波基的細化工作
/*
函數輸入：
TYPE_ARRAY *aO：	輸入的原始數據，用于細化的小波基；
int nAO[]：		每一維上數據aO的大小，這個數據的大小必須是2的指數值；
int nDim：		小波基的維數以及數組nAO[]的大小；
bool isFwd：		小波基是否完整正變換（如果是，則取真值）
waveletfilter *wfltr：	需要進行細化的小波濾波器；
int nAR[]：		細化后小波基每一維的大小（數值必須是2的指數值）
TYPE_ARRAY *aR：	函數的實際輸出，細化后的數據。
*/
void FUNC_ND(aO, nAO, nDim, isFwd, wfltr, nAR, aR)
	TYPE_ARRAY *aO;	
	int nAO[];	
	int nDim;	
	bool isFwd;		
	waveletfilter *wfltr;	
	int nAR[];	
	TYPE_ARRAY *aR;	
{
	int k;
	int iDim, nAOTotal, nAOMax, nARMax;
	int iO1, iO2, iO3, nO, nONext
	int nOPrev = 1;
	int nRPrev = 1;
	int nR, nRNext;
	TYPE_ARRAY *aOTmp, *aRTmp;

	nAOMax = 1;
	nARMax = 1;
	nAOTotal = 1;
	for (iDim = 0; iDim < nDim; iDim++) {
		if (nAO[iDim] > nAOMax)
			nAOMax = nAO[iDim];
		if (nAR[iDim] > nARMax)
			nARMax = nAR[iDim];
		nAOTotal *= nAO[iDim];
	}

	//臨時數據空間必須足夠大，這樣能保證最大維數小波基的細化
	(void) MALLOC_LINTOK(aOTmp, nAOMax, TYPE_ARRAY);
	(void) MALLOC_LINTOK(aRTmp, nARMax, TYPE_ARRAY);

	for (iDim = 0; iDim < nDim; iDim++) {
		nO = nAO[nDim - 1 - iDim];
		nR = nAR[nDim - 1 - iDim];
		nONext = nO * nOPrev;
		nRNext = nR * nRPrev;
		if (nO >= MIN_ORDER) {
			for (iO2 = 0; iO2 < nAOTotal; iO2 += nONext) {
				for (iO1 = 0; iO1 < nOPrev; iO1++) {					
					 //將相關的行、列以及小波基其他信息進行復制
					iO3 = iO1 + iO2;
					for (k = 0; k < nO; k++) {
						aOTmp[k] = aO[iO3];
						iO3 += nOPrev;
					}
					
					//調用一維小波基的細化操作函數
					FUNC_1D(aOTmp, nO, nR, wfltr, aRTmp);

					//將細化的結果從林紓數據空間中粘貼至輸出小波基中
					i3 = i1 + i2;
				for (k = 0; k < nR; k++) {
						aR[i3] = aRTmp[k];
						i3 += nRPrev;
					}
				}
			}
		}
		nOPrev = nONext;
		nRPrev = nRNext;
	}

	FREE_LINTOK(aRTmp);
	FREE_LINTOK(aOTmp);

	return;
}
#endif

/* wrefine_step -- perform one iteration of a wavelet refinement */
//完成小波細化操作中一次迭代循環操作
static void wrefine_step(wfltr, a, n)
	waveletfilter *wfltr;
	TYPE_ARRAY *a;
	int n;
{
	int nDiv2 = n / 2;
	int i, iHtilde, j, iA;
	TYPE_ARRAY *aNew;

	(void) MALLOC_LINTOK(aNew, n, TYPE_ARRAY);
	for (i = 0; i < n; i++)
		aNew[i] = 0.0;

	/*
	 *	Wavelet refinement is identical to an inverse wavelet transform (see
	 *	"wxfrm") with the elements of Gtilde set identically to zero.
	 *
	 *	A single step of the reconstruction is summarized by:
	 *
	 *		aNew[0..n-1] = Htilde^t * a[0..n/2-1]
	 *
	 *	where Htilde^t is the transpose of Htilde.
	 */
	for (j = 0; j < nDiv2; j++) {
		for (iHtilde = 0; iHtilde < wfltr->nHtilde; iHtilde++) {
			/*
			 *	Each row of Htilde is offset by 2 from the previous one.
			 */
			iA = MOD(2 * j + iHtilde - wfltr->offHtilde, n);
			aNew[iA] += wfltr->cHtilde[iHtilde] * a[j];
		}
	}
	for (i = 0; i < n; i++)
		a[i] = aNew[i];

	FREE_LINTOK(aNew);

	return;
}