shithub: aacenc

Download patch

ref: 53eb89ba87c0d1ae8d61a4cf12af83427f5ce4b4
parent: aaf7c645525630b7e298c163b84b42cc89b8e2d8
author: menno <menno>
date: Wed Jan 5 16:41:37 EST 2000

Added Long Term Prediction (LTP), new decoder is needed for this

--- a/Makefile
+++ b/Makefile
@@ -1,4 +1,4 @@
-SOURCE=aac_back_pred.c aac_qc.c aac_se_enc.c bitstream.c enc_tf.c encoder.c imdct.c is.c mc_enc.c ms.c psych.c pulse.c tns.c transfo.c fastfft.c
+SOURCE=aac_qc.c aac_se_enc.c bitstream.c enc_tf.c encoder.c imdct.c is.c mc_enc.c ms.c psych.c pulse.c tns.c transfo.c fastfft.c nok_ltp_enc.c nok_pitch.c
 
 
 OBJ = $(SOURCE:.c=.o)
--- a/aac_qc.h
+++ b/aac_qc.h
@@ -7,6 +7,7 @@
 #include "tf_main.h"
 #include "tns.h"
 #include "all.h"
+#include "nok_ltp_common.h"
 
 #ifdef __cplusplus
 extern "C" {
@@ -56,6 +57,7 @@
   int reset_group_number;               /* Prediction reset group number */
   TNS_INFO* tnsInfo;                    /* Ptr to tns data */
   AACPulseInfo pulseInfo;
+  NOK_LT_PRED_STATUS *ltpInfo;          /* Prt to LTP data */
 } AACQuantInfo;
 
 
--- a/aac_se_enc.c
+++ b/aac_se_enc.c
@@ -281,54 +281,25 @@
       BsPutBit(fixed_stream,max_sfb,LEN_MAX_SFBL);
     }
     bit_count += LEN_MAX_SFBL;
-    bit_count += WritePredictorData(quantInfo,fixed_stream,writeFlag);
+    bit_count += WriteLTP_PredictorData(quantInfo,fixed_stream,writeFlag);
   }
 
   return bit_count;
 }
 
-
 /*****************************************************************************/
-/* WritePredictorData(...), write predictor data.                            */
+/* WriteLTP_PredictorData(...), write LTP predictor data.                    */
 /*****************************************************************************/
-int WritePredictorData(AACQuantInfo* quantInfo,    /* AACQuantInfo structure */
-		       BsBitStream* fixed_stream,  /* Pointer to bitstream */
-		       int writeFlag)              /* 1 means write, 0 means count only */  
+int WriteLTP_PredictorData(AACQuantInfo* quantInfo,    /* AACQuantInfo structure */
+                           BsBitStream* fixed_stream,  /* Pointer to bitstream */
+                           int writeFlag)              /* 1 means write, 0 means count only */  
 {
 	int bit_count = 0;
 
-	/* Write global predictor data present */
-	short predictorDataPresent = quantInfo->pred_global_flag;
-	int numBands = min(max_pred_sfb,quantInfo->nr_of_sfb);
+	bit_count += nok_ltp_encode (fixed_stream, quantInfo->block_type, quantInfo->nr_of_sfb, 
+		quantInfo->ltpInfo, writeFlag);
 
-	if (writeFlag) {
-		BsPutBit(fixed_stream,predictorDataPresent,LEN_PRED_PRES);  /* predictor_data_present */
-		if (predictorDataPresent) {
-			int b;
-
-			/* Code segment added by JB */
-			if (quantInfo->reset_group_number == -1)
-				BsPutBit(fixed_stream,0,LEN_PRED_RST); /* No prediction reset */
-			else
-			{
-				BsPutBit(fixed_stream,1,LEN_PRED_RST);
-				BsPutBit(fixed_stream,(unsigned long)quantInfo->reset_group_number,
-					LEN_PRED_RSTGRP);
-			}
-			/* End of code segment */
-
-			for (b=0;b<numBands;b++) {
-				BsPutBit(fixed_stream,quantInfo->pred_sfb_flag[b],LEN_PRED_ENAB);
-			}
-		}
-	}
-	bit_count = LEN_PRED_PRES;
-	bit_count += (predictorDataPresent) ?
-		(LEN_PRED_RST + 
-		((quantInfo->reset_group_number)!=-1)*LEN_PRED_RSTGRP + 
-		numBands*LEN_PRED_ENAB) : 0;
-
-	return bit_count;
+	return (bit_count);
 }
 
 /*****************************************************************************/
--- a/aac_se_enc.h
+++ b/aac_se_enc.h
@@ -46,6 +46,7 @@
 #include "bitstream.h"
 #include "all.h"
 #include "aac_qc.h"
+#include "nok_ltp_enc.h"
 
 extern int max_pred_sfb;
 
@@ -98,6 +99,12 @@
 		 BsBitStream* fixedStream, /* Pointer to bitstream */
 		 int writeFlag);           /* 1 means write, 0 means count only */
 
+/*****************************************************************************/
+/* WriteLTP_PredictorData(...), write LTP predictor data.                    */
+/*****************************************************************************/
+int WriteLTP_PredictorData(AACQuantInfo* quantInfo,  /* Pointer to AACQuantInfo structure */
+                           BsBitStream* fixedStream, /* Pointer to bitstream */
+                           int writeFlag);           /* 1 means write, 0 means count only */
 
 /*****************************************************************************/
 /* WritePredictorData(...), write predictor data.                            */
--- a/all.h
+++ b/all.h
@@ -60,5 +60,11 @@
     MS_Info ms_info;    /* MS information */
 } Ch_Info;
 
+#ifdef MPEG4V1
+#define _WINDOW_TYPE_EXT
+typedef byte WINDOW_TYPE_EXT;
+#include "nok_lt_prediction.h"
+#endif /* MPEG4V1 */
+
 #endif	/* _all_h_ */
 
--- a/enc_tf.c
+++ b/enc_tf.c
@@ -9,7 +9,6 @@
 #include "block.h"
 #include "tf_main.h"
 #include "psych.h"
-#include "aac_back_pred.h"
 #include "mc_enc.h"
 #include "ms.h"
 #include "is.h"
@@ -16,6 +15,8 @@
 #include "aac_qc.h"
 #include "all.h"
 #include "aac_se_enc.h"
+#include "nok_ltp_enc.h"
+
 #define SQRT2 sqrt(2)
 
 /* AAC tables */
@@ -38,6 +39,7 @@
 double *overlap_buffer[MAX_TIME_CHANNELS];
 double *DTimeSigBuf[MAX_TIME_CHANNELS];
 double *DTimeSigLookAheadBuf[MAX_TIME_CHANNELS+2];
+double *nok_tmp_DTimeSigBuf[MAX_TIME_CHANNELS]; /* temporary fix to the buffer size problem. */
 
 /* static variables used by the T/F mapping */
 enum QC_MOD_SELECT qc_select = AAC_QC;                   /* later f(encPara) */
@@ -49,6 +51,7 @@
 /* Additional variables for AAC */
 int aacAllowScalefacs = 1;              /* Allow AAC scalefactors to be nonconstant */
 TNS_INFO tnsInfo[MAX_TIME_CHANNELS];
+NOK_LT_PRED_STATUS nok_lt_status[MAX_TIME_CHANNELS];
 
 AACQuantInfo quantInfo[MAX_TIME_CHANNELS];               /* Info structure for AAC quantization and coding */
 
@@ -148,6 +151,9 @@
 		overlap_buffer[chanNum] = (double*)malloc(sizeof(double)*block_size_samples);
 		memset(overlap_buffer[chanNum],0,(block_size_samples)*sizeof(double));
 		block_type[chanNum] = ONLY_LONG_WINDOW;
+		nok_lt_status[chanNum].delay =  (int*)malloc(MAX_SHORT_WINDOWS*sizeof(int));
+		nok_tmp_DTimeSigBuf[chanNum]  = (double*)malloc(2*block_size_samples*sizeof(double));
+		memset(nok_tmp_DTimeSigBuf[chanNum],0,(2*block_size_samples)*sizeof(double));
 	}
 	for (chanNum=0;chanNum<MAX_TIME_CHANNELS+2;chanNum++) {
 		DTimeSigLookAheadBuf[chanNum]   = (double*)malloc((block_size_samples)*sizeof(double));
@@ -154,8 +160,6 @@
 		memset(DTimeSigLookAheadBuf[chanNum],0,(block_size_samples)*sizeof(double));
 	}
 
-	PredInit();
-
 	/* initialize psychoacoustic module */
 	EncTf_psycho_acoustic_init();
 
@@ -168,6 +172,12 @@
 		TnsInit(sampling_rate,profile,&tnsInfo[chanNum]);
 		quantInfo[chanNum].tnsInfo = &tnsInfo[chanNum];         /* Set pointer to TNS data */
 	}
+
+	/* Init LTP predictor */
+	for (chanNum=0;chanNum<MAX_TIME_CHANNELS;chanNum++) {
+		nok_init_lt_pred (&nok_lt_status[chanNum]);
+		quantInfo[chanNum].ltpInfo = &nok_lt_status[chanNum];  /* Set pointer to LTP data */
+	}
 }
 
 /*****************************************************************************************
@@ -222,12 +232,7 @@
 	CH_PSYCH_OUTPUT_LONG chpo_long[MAX_TIME_CHANNELS+2];
 	CH_PSYCH_OUTPUT_SHORT chpo_short[MAX_TIME_CHANNELS+2][MAX_SHORT_WINDOWS];
 
-//	memset(chpo_long,0,sizeof(CH_PSYCH_OUTPUT_LONG)*(MAX_TIME_CHANNELS+2));
-//	memset(chpo_short,0,sizeof(CH_PSYCH_OUTPUT_SHORT)*(MAX_TIME_CHANNELS+2)*MAX_SHORT_WINDOWS);
-//	memset(p_ratio_long,0,sizeof(double)*(MAX_TIME_CHANNELS)*MAX_SCFAC_BANDS);
-//	memset(p_ratio_short,0,sizeof(double)*(MAX_TIME_CHANNELS)*MAX_SCFAC_BANDS);
-
-	{ /* convert float input to double, which is the internal format */
+	{
 		/* store input data in look ahead buffer which may be necessary for the window switching decision */
 		int i;
 		int chanNum;
@@ -234,6 +239,11 @@
 		
 		for (chanNum=0;chanNum<max_ch;chanNum++) {
 			for( i=0; i<block_size_samples; i++ ) {
+				/* temporary fix: a linear buffer for LTP containing the whole time frame */
+				nok_tmp_DTimeSigBuf[chanNum][i] = DTimeSigBuf[chanNum][i];
+				nok_tmp_DTimeSigBuf[chanNum][block_size_samples + i] = DTimeSigLookAheadBuf[chanNum][i];
+			}
+			for( i=0; i<block_size_samples; i++ ) {
 				/* last frame input data are encoded now */
 				DTimeSigBuf[chanNum][i] = DTimeSigLookAheadBuf[chanNum][i];
 				DTimeSigLookAheadBuf[chanNum][i] = as->inputBuffer[chanNum][i];
@@ -329,8 +339,7 @@
 //	block_type[1] = ONLY_LONG_WINDOW;
 //	block_type[0] = ONLY_SHORT_WINDOW;
 //	block_type[1] = ONLY_SHORT_WINDOW;
-//	if (as->use_MS)
-//		block_type[1] = block_type[0];
+	block_type[1] = block_type[0];
 
 	{
 		int chanNum;
@@ -364,7 +373,7 @@
 				quantInfo[chanNum].window_group_length[7] = 0;
 #endif
 				break;
-				
+
 			default:
 				no_sub_win   = 1;
 				sub_win_size = block_size_samples;
@@ -446,24 +455,8 @@
 		}
 	}
 
-//	if (as->use_MS) {
-		MSPreprocess(p_ratio_long, p_ratio_short, chpo_long, chpo_short,
-			channelInfo, block_type, quantInfo, as->use_MS, max_ch);
-//	} else {
-//		int chanNum;
-//		for (chanNum=0;chanNum<max_ch;chanNum++) {
-//
-//			/* Save p_ratio from psychoacoustic model for next frame.  */
-//			/* Psycho model is using a look-ahead window for block switching */
-//			if (as->use_MS) {
-//				memcpy( (char*)p_ratio_long[chanNum], (char*)chpo_long[chanNum+2].p_ratio, (NSFB_LONG)*sizeof(double) );
-//				memcpy( (char*)p_ratio_short[chanNum],(char*)chpo_short[chanNum+2][0].p_ratio,(MAX_SHORT_WINDOWS*NSFB_SHORT)*sizeof(double) );
-//			} else {
-//				memcpy( (char*)p_ratio_long[chanNum], (char*)chpo_long[chanNum].p_ratio, (NSFB_LONG)*sizeof(double) );
-//				memcpy( (char*)p_ratio_short[chanNum],(char*)chpo_short[chanNum][0].p_ratio,(MAX_SHORT_WINDOWS*NSFB_SHORT)*sizeof(double) );
-//			}
-//		}
-//	}
+	MSPreprocess(p_ratio_long, p_ratio_short, chpo_long, chpo_short,
+		channelInfo, block_type, quantInfo, as->use_MS, max_ch);
 
 	MSEnergy(spectral_line_vector, energy, chpo_long, chpo_short, sfb_width_table,
 		channelInfo, block_type, quantInfo, as->use_MS, max_ch);
@@ -532,27 +525,72 @@
 				max_ch);
 		}
 
-		/***********************************************************************/
-		/* If prediction is used, compute predictor info and residual spectrum */
-		/***********************************************************************/
-		for (chanNum=0;chanNum<max_ch;chanNum++) {
-//			if (qc_select == AAC_PRED) {
-			if (0) {
-				int numPredBands;
-				max_pred_sfb = 40;
-				numPredBands = min(max_pred_sfb,nr_of_sfb[chanNum]);
-				PredCalcPrediction( spectral_line_vector[chanNum],
-					reconstructed_spectrum[chanNum],
-					(int)block_type[chanNum],
-					numPredBands,
-					sfb_width_table[chanNum],
-					&(quantInfo[chanNum].pred_global_flag),
-					quantInfo[chanNum].pred_sfb_flag,
-					&(quantInfo[chanNum].reset_group_number),
-					chanNum);
-			} else {
-				quantInfo[chanNum].pred_global_flag = 0;
-			}
+		/*******************************************************************************/
+		/* If LTP prediction is used, compute LTP predictor info and residual spectrum */
+		/*******************************************************************************/
+		for(chanNum=0;chanNum<max_ch;chanNum++) 
+		{
+			if(block_type[chanNum] != ONLY_SHORT_WINDOW) 
+			{
+				if(channelInfo[chanNum].cpe)
+				{
+					if(channelInfo[chanNum].ch_is_left) 
+					{
+						int i;
+						int leftChan=chanNum;
+						int rightChan=channelInfo[chanNum].paired_ch;
+                        
+						nok_ltp_enc(spectral_line_vector[leftChan], 
+							nok_tmp_DTimeSigBuf[leftChan], 
+							block_type[leftChan], 
+							WS_FHG,
+							block_size_samples,
+							block_size_samples/2,
+							block_size_samples/short_win_in_long, 
+							&sfb_offset_table[leftChan][0], 
+							nr_of_sfb[leftChan],
+							&nok_lt_status[leftChan]);
+
+						nok_lt_status[rightChan].global_pred_flag = 
+							nok_lt_status[leftChan].global_pred_flag;
+						for(i = 0; i < NOK_MAX_BLOCK_LEN_LONG; i++)
+							nok_lt_status[rightChan].pred_mdct[i] = 
+							nok_lt_status[leftChan].pred_mdct[i];
+						for(i = 0; i < MAX_SCFAC_BANDS; i++)
+							nok_lt_status[rightChan].sfb_prediction_used[i] = 
+							nok_lt_status[leftChan].sfb_prediction_used[i];
+						nok_lt_status[rightChan].weight = nok_lt_status[leftChan].weight;
+						nok_lt_status[rightChan].delay[0] = nok_lt_status[leftChan].delay[0];
+
+						if(block_type[leftChan] != block_type[rightChan])
+							nok_ltp_enc(spectral_line_vector[rightChan],
+							nok_tmp_DTimeSigBuf[rightChan], 
+							block_type[rightChan], 
+							WS_FHG,
+							block_size_samples,
+							block_size_samples/2,
+							block_size_samples/short_win_in_long, 
+							&sfb_offset_table[rightChan][0], 
+							nr_of_sfb[rightChan],
+							&nok_lt_status[rightChan]);
+
+					} /* if(channelInfo[chanNum].ch_is_left) */
+				} /* if(channelInfo[chanNum].cpe) */
+				else
+					nok_ltp_enc(spectral_line_vector[chanNum], 
+					nok_tmp_DTimeSigBuf[chanNum], 
+					block_type[chanNum], 
+					WS_FHG,
+					block_size_samples,
+					block_size_samples/2,
+					block_size_samples/short_win_in_long, 
+					&sfb_offset_table[chanNum][0], 
+					nr_of_sfb[chanNum],
+					&nok_lt_status[chanNum]);
+
+			} /* if(channelInfo[chanNum].present... */
+			else
+				quantInfo[chanNum].ltpInfo->global_pred_flag = 0;
 		} /* for(chanNum... */
 
 		/******************************************/
@@ -607,13 +645,49 @@
 				return error;
 		}
 
-		/* If short window, reconstruction not needed for prediction */
+		/**********************************************************/
+		/* Reconstruct MS Stereo bands for prediction            */
+		/**********************************************************/
+		if (as->use_MS != -1) {
+			MSReconstruct(reconstructed_spectrum,
+				channelInfo,
+				sfb_offset_table,
+				block_type,
+				quantInfo,
+				max_ch);
+		}
+
+		/**********************************************************/
+		/* Reconstruct Intensity Stereo bands for prediction     */
+		/**********************************************************/
+		if (as->use_IS) {
+			ISReconstruct(reconstructed_spectrum,
+				channelInfo,
+				sfb_offset_table,
+				block_type,
+				quantInfo,
+				max_ch);
+		}
+
+		/**********************************************************/
+		/* Update LTP history buffer                              */
+		/**********************************************************/
 		for (chanNum=0;chanNum<max_ch;chanNum++) {
-			if ((block_type[chanNum]==ONLY_SHORT_WINDOW)) {
-				int sind;
-				for (sind=0;sind<1024;sind++) {
-					reconstructed_spectrum[chanNum][sind]=0.0;
-				}
+			nok_ltp_reconstruct(reconstructed_spectrum[chanNum], 
+				block_type[chanNum], 
+				WS_FHG, block_size_samples,
+				block_size_samples/2,
+				block_size_samples/short_win_in_long, 
+				&sfb_offset_table[chanNum][0], 
+				nr_of_sfb[chanNum],
+				&nok_lt_status[chanNum]);
+		}
+
+		/* If short window, reconstruction not needed for prediction */
+		if ((block_type[chanNum]==ONLY_SHORT_WINDOW)) {
+			int sind;
+			for (sind=0;sind<1024;sind++) {
+				reconstructed_spectrum[chanNum][sind]=0.0;
 			}
 		}
 
--- a/faac.dsp
+++ b/faac.dsp
@@ -51,7 +51,7 @@
 LINK32=link.exe
 # ADD BASE LINK32 kernel32.lib user32.lib gdi32.lib winspool.lib comdlg32.lib advapi32.lib shell32.lib ole32.lib oleaut32.lib uuid.lib odbc32.lib odbccp32.lib kernel32.lib user32.lib gdi32.lib winspool.lib comdlg32.lib advapi32.lib shell32.lib ole32.lib oleaut32.lib uuid.lib odbc32.lib odbccp32.lib /nologo /subsystem:console /machine:I386
 # ADD LINK32 libsndfile.lib kernel32.lib user32.lib gdi32.lib winspool.lib comdlg32.lib advapi32.lib shell32.lib ole32.lib oleaut32.lib uuid.lib odbc32.lib odbccp32.lib /nologo /subsystem:console /machine:I386
-# SUBTRACT LINK32 /profile
+# SUBTRACT LINK32 /profile /nodefaultlib
 
 !ELSEIF  "$(CFG)" == "faac - Win32 Debug"
 
@@ -88,10 +88,6 @@
 # PROP Default_Filter "cpp;c;cxx;rc;def;r;odl;idl;hpj;bat"
 # Begin Source File
 
-SOURCE=.\aac_back_pred.c
-# End Source File
-# Begin Source File
-
 SOURCE=.\aac_qc.c
 # End Source File
 # Begin Source File
@@ -129,6 +125,14 @@
 # Begin Source File
 
 SOURCE=.\ms.c
+# End Source File
+# Begin Source File
+
+SOURCE=.\nok_ltp_enc.c
+# End Source File
+# Begin Source File
+
+SOURCE=.\nok_pitch.c
 # End Source File
 # Begin Source File
 
--- a/faac_dll.dsp
+++ b/faac_dll.dsp
@@ -92,10 +92,6 @@
 # PROP Default_Filter "cpp;c;cxx;rc;def;r;odl;idl;hpj;bat"
 # Begin Source File
 
-SOURCE=.\aac_back_pred.c
-# End Source File
-# Begin Source File
-
 SOURCE=.\aac_qc.c
 # End Source File
 # Begin Source File
@@ -133,6 +129,14 @@
 # Begin Source File
 
 SOURCE=.\ms.c
+# End Source File
+# Begin Source File
+
+SOURCE=.\nok_ltp_enc.c
+# End Source File
+# Begin Source File
+
+SOURCE=.\nok_pitch.c
 # End Source File
 # Begin Source File
 
--- a/imdct.c
+++ b/imdct.c
@@ -111,90 +111,300 @@
   Mdct_in      overlap_select     /*  YT 970615 for son_PP */
 )
 {
-	double         transf_buf[ 2*MAX_SHIFT_LEN_LONG ];
-	double         windowed_buf[ 2*MAX_SHIFT_LEN_LONG ];
-	double         *p_o_buf;
-	int            k;
+  double         transf_buf[ 2*MAX_SHIFT_LEN_LONG ];
+  double         windowed_buf[ 2*MAX_SHIFT_LEN_LONG ];
+  double         *p_o_buf;
+  int            k;
 
-	double           window_long[MAX_SHIFT_LEN_LONG]; 
-	double           window_long_prev[MAX_SHIFT_LEN_LONG]; 
-	double           window_short[MAX_SHIFT_LEN_LONG]; 
-	double           window_short_prev[MAX_SHIFT_LEN_LONG]; 
-	double           *window_short_prev_ptr;
+  double           window_long[MAX_SHIFT_LEN_LONG]; 
+  double           window_long_prev[MAX_SHIFT_LEN_LONG]; 
+  double           window_med[MAX_SHIFT_LEN_LONG]; 
+  double           window_med_prev[MAX_SHIFT_LEN_LONG]; 
+  double           window_short[MAX_SHIFT_LEN_LONG]; 
+  double           window_short_prev[MAX_SHIFT_LEN_LONG]; 
+  double           *window_short_prev_ptr;
 
-	int nflat_ls    = (nlong-nshort)/ 2; 
-	int transfak_ls =  nlong/nshort; 
-	
-	static Window_shape wfun_select_prev=WS_FHG;
-	static int firstTime=1;
-	window_short_prev_ptr = window_short_prev ; 
 
-	calc_window( window_long,      nlong, wfun_select ); 
-	calc_window( window_long_prev, nlong, wfun_select_prev ); 
-	calc_window( window_short,      nshort, wfun_select ); 
-	calc_window( window_short_prev, nshort, wfun_select_prev ); 
-	
-	/* create / shift old values */
-	/* We use p_overlap here as buffer holding the last frame time signal*/
-	/* YT 970615 for son_pp */
-	if (firstTime){
-		firstTime=0;
-		vcopy( &zero, transf_buf, 0, 1, nlong );
-	} else
-		vcopy( p_overlap, transf_buf, 1, 1, nlong );
+  int nflat_ls    = (nlong-nshort)/ 2; 
+  int transfak_ls =  nlong/nshort; 
+  int nflat_lm    = (nlong-nmed)  / 2; 
+  int transfak_lm =  nlong/nmed; 
+  int nflat_ms    = (nmed-nshort) / 2; 
+  int transfak_ms =  nmed/nshort; 
+ 
+  static Window_shape wfun_select_prev=WS_FHG;
+  static int firstTime=1;
+  window_short_prev_ptr = window_short_prev ; 
 
-	/* Append new data */
-	vcopy( p_in_data, transf_buf+nlong, 1, 1, nlong );
-	vcopy( p_in_data, p_overlap,        1, 1, nlong );
 
-	/* Set ptr to transf-Buffer */
-	p_o_buf = transf_buf;
-	
-	
-	/* Separate action for each Block Type */
-	switch( block_type ) {
-	case ONLY_LONG_WINDOW :
-		vmult1( p_o_buf, window_long_prev, windowed_buf, nlong );
-		vmult2( p_o_buf+nlong, window_long+nlong-1, windowed_buf+nlong, nlong );
-		mdct( windowed_buf, p_out_mdct, 2*nlong );    
-		break;
-		
-	case LONG_SHORT_WINDOW :
-		vmult1( p_o_buf, window_long_prev, windowed_buf, nlong );
-		vcopy( p_o_buf+nlong, windowed_buf+nlong, 1, 1, nflat_ls );
-		vmult2( p_o_buf+nlong+nflat_ls, window_short+nshort-1, windowed_buf+nlong+nflat_ls, nshort );
-		vcopy( &zero, windowed_buf+2*nlong-1, 0, -1, nflat_ls );
-		mdct( windowed_buf, p_out_mdct, 2*nlong );
-		break;
+// if( (nlong%nshort) || (nlong > MAX_SHIFT_LEN_LONG) || (nshort > MAX_SHIFT_LEN_LONG/2) ) { 
+//    CommonExit( 1, "mdct_analysis: Problem with window length" ); } 
+//  if( (nlong%nmed) || (nmed%nshort) || (nmed > MAX_SHIFT_LEN_LONG/2) ) { 
+//    CommonExit( 1, "mdct_analysis: Problem with window length" );  } 
+//  if( transfak_lm%2 ) { 
+//    CommonExit( 1, "mdct_analysis: Problem with window length" );  } 
+  
+  /* --- PATCH: Use WS_FHG window shape if we start with SHORT windows --- */
+  if  (block_type==LONG_SHORT_WINDOW || block_type==ONLY_SHORT_WINDOW ){
+    wfun_select=WS_FHG; } 
+  /* --- PATCH  End  --- */
+  calc_window( window_long,      nlong, wfun_select ); 
+  calc_window( window_long_prev, nlong, wfun_select_prev ); 
+  calc_window( window_med, nmed, wfun_select ); 
+  calc_window( window_med_prev, nmed, wfun_select_prev ); 
+  calc_window( window_short,      nshort, wfun_select ); 
+  calc_window( window_short_prev, nshort, wfun_select_prev ); 
+  
+  /* create / shift old values */
+  /* We use p_overlap here as buffer holding the last frame time signal*/
+/* YT 970615 for son_pp */
+if(overlap_select != MNON_OVERLAPPED){
+  if (firstTime){
+    firstTime=0;
+    vcopy( &zero, transf_buf, 0, 1, nlong );
+  }
+  else
+    vcopy( p_overlap, transf_buf, 1, 1, nlong );
 
-	case SHORT_LONG_WINDOW :
-		vcopy( &zero, windowed_buf, 0, 1, nflat_ls );
-		vmult1( p_o_buf+nflat_ls, window_short_prev_ptr, windowed_buf+nflat_ls, nshort );
-		vcopy( p_o_buf+nflat_ls+nshort, windowed_buf+nflat_ls+nshort, 1, 1, nflat_ls );
-		vmult2( p_o_buf+nlong, window_long+nlong-1, windowed_buf+nlong, nlong );
-		mdct( windowed_buf, p_out_mdct, 2*nlong );
-		break;
+  /* Append new data */
+  vcopy( p_in_data, transf_buf+nlong, 1, 1, nlong );
+  vcopy( p_in_data, p_overlap,        1, 1, nlong );
+}
+else{
+	vcopy( p_in_data, transf_buf, 1, 1, 2*nlong );
+}
 
-	case ONLY_SHORT_WINDOW :
-		p_o_buf += nflat_ls;
-		for (k=transfak_ls-1; k-->=0; ) {
-			vmult1( p_o_buf, window_short_prev_ptr, windowed_buf, nshort );
-			vmult2( p_o_buf+nshort, window_short+nshort-1, windowed_buf+nshort, nshort );
-			mdct( windowed_buf, p_out_mdct, 2*nshort );
+  /* Set ptr to transf-Buffer */
+  p_o_buf = transf_buf;
+  
+  
+  /* Separate action for each Block Type */
+  switch( block_type ) {
+   case ONLY_LONG_WINDOW :
+     vmult1( p_o_buf, window_long_prev, windowed_buf, nlong );
+     vmult2( p_o_buf+nlong, window_long+nlong-1, windowed_buf+nlong, nlong );
+     mdct( windowed_buf, p_out_mdct, 2*nlong );    
+     break;
+   
+   case LONG_SHORT_WINDOW :
+    vmult1( p_o_buf, window_long_prev, windowed_buf, nlong );
+    vcopy( p_o_buf+nlong, windowed_buf+nlong, 1, 1, nflat_ls );
+    vmult2( p_o_buf+nlong+nflat_ls, window_short+nshort-1, windowed_buf+nlong+nflat_ls, nshort );
+    vcopy( &zero, windowed_buf+2*nlong-1, 0, -1, nflat_ls );
+    mdct( windowed_buf, p_out_mdct, 2*nlong );
+    break;
 
-			p_out_mdct += nshort;
-			/* YT 970615 for sonPP*/
-			p_o_buf += nshort;
-			window_short_prev_ptr = window_short; 
-		}
-		break;
+   case SHORT_LONG_WINDOW :
+    vcopy( &zero, windowed_buf, 0, 1, nflat_ls );
+    vmult1( p_o_buf+nflat_ls, window_short_prev_ptr, windowed_buf+nflat_ls, nshort );
+    vcopy( p_o_buf+nflat_ls+nshort, windowed_buf+nflat_ls+nshort, 1, 1, nflat_ls );
+    vmult2( p_o_buf+nlong, window_long+nlong-1, windowed_buf+nlong, nlong );
+    mdct( windowed_buf, p_out_mdct, 2*nlong );
+    break;
+
+   case ONLY_SHORT_WINDOW :
+	if(overlap_select != MNON_OVERLAPPED){ /* YT 970615 for sonPP */
+    	p_o_buf += nflat_ls;
 	}
+    for (k=transfak_ls-1; k-->=0; ) {
+      vmult1( p_o_buf,        window_short_prev_ptr, windowed_buf, nshort );
+      vmult2( p_o_buf+nshort, window_short+nshort-1, windowed_buf+nshort, nshort );
+      mdct( windowed_buf, p_out_mdct, 2*nshort );
 
-	/* Set output data 
-	vcopy(transf_buf, p_out_mdct,1, 1, nlong); */
-	
-	/* --- Save old window shape --- */
-	wfun_select_prev = wfun_select;
+      p_out_mdct += nshort;
+	/* YT 970615 for sonPP*/
+	  if(overlap_select != MNON_OVERLAPPED) p_o_buf += nshort;
+      else 	p_o_buf += 2*nshort;
+      window_short_prev_ptr = window_short; 
+    }
+    break;
+   default :
+    break;
+//      CommonExit( 1, "mdct_synthesis: Unknown window type" ); 
+  }
+
+  /* Set output data 
+  vcopy(transf_buf, p_out_mdct,1, 1, nlong); */
+  
+  /* --- Save old window shape --- */
+  wfun_select_prev = wfun_select;
 }
+
+void imdct(double in_data[], double out_data[], int len)
+{
+  vcopy(in_data, out_data, 1, 1, len/2);
+  IMDCT(out_data, len);
+}
+
+void freq2buffer(
+  double           p_in_data[], 
+  double           p_out_data[],
+  double           p_overlap[],
+  enum WINDOW_TYPE block_type,
+  int              nlong,            /* shift length for long windows   */
+  int              nmed,             /* shift length for medium windows */
+  int              nshort,           /* shift length for short windows  */
+  Window_shape     wfun_select,      /* offers the possibility to select different window functions */
+  Mdct_in	   overlap_select    /* select imdct output *TK*	*/
+				     /* switch (overlap_select) {	*/
+				     /* case OVERLAPPED:		*/
+				     /*   p_out_data[]			*/
+				     /*   = overlapped and added signal */
+				     /*		(bufferlength: nlong)	*/
+				     /* case MNON_OVERLAPPED:		*/
+				     /*   p_out_data[]			*/
+				     /*   = non overlapped signal	*/
+				     /*		(bufferlength: 2*nlong)	*/
+  )
+{
+  double           *o_buf, transf_buf[ 2*MAX_SHIFT_LEN_LONG ];
+  double           overlap_buf[ 2*MAX_SHIFT_LEN_LONG ]; 
+ 
+  double           window_long[MAX_SHIFT_LEN_LONG]; 
+  double           window_long_prev[MAX_SHIFT_LEN_LONG]; 
+  double           window_med[MAX_SHIFT_LEN_LONG]; 
+  double           window_med_prev[MAX_SHIFT_LEN_LONG]; 
+  double           window_short[MAX_SHIFT_LEN_LONG]; 
+  double           window_short_prev[MAX_SHIFT_LEN_LONG]; 
+  double           *window_short_prev_ptr;
+ 
+  double  *fp; 
+  int     k; 
+ 
+  int nflat_ls    = (nlong-nshort)/ 2; 
+  int transfak_ls =  nlong/nshort; 
+  int nflat_lm    = (nlong-nmed)  / 2; 
+  int transfak_lm =  nlong/nmed; 
+  int nflat_ms    = (nmed-nshort) / 2; 
+  int transfak_ms =  nmed/nshort;  
+ 
+  static Window_shape wfun_select_prev=WS_FHG;
+  window_short_prev_ptr=window_short_prev ; 
+
+
+  if( (nlong%nshort) || (nlong > MAX_SHIFT_LEN_LONG) || (nshort > MAX_SHIFT_LEN_LONG/2) ) { 
+//    CommonExit( 1, "mdct_synthesis: Problem with window length" ); 
+  } 
+  if( (nlong%nmed) || (nmed%nshort) || (nmed > MAX_SHIFT_LEN_LONG/2) ) { 
+//    CommonExit( 1, "mdct_synthesis: Problem with window length" ); 
+  } 
+  if( transfak_lm%2 ) { 
+//    CommonExit( 1, "mdct_synthesis: Problem with window length" ); 
+  } 
+
+
+
+  /* --- PATCH: Use WS_FHG window shape if we start 
+                with SHORT windows            --- */
+  if  (block_type==LONG_SHORT_WINDOW ||
+       block_type==ONLY_SHORT_WINDOW ){
+    wfun_select=WS_FHG;
+  } 
+  /* --- PATCH  End  --- */
+  
+  calc_window( window_long,      nlong, wfun_select ); 
+  calc_window( window_long_prev, nlong, wfun_select_prev ); 
+  calc_window( window_med, nmed, wfun_select ); 
+  calc_window( window_med_prev, nmed, wfun_select_prev );
+  calc_window( window_short,      nshort, wfun_select ); 
+  calc_window( window_short_prev_ptr, nshort, wfun_select_prev ); 
+  
+ 
+  /* Assemble overlap buffer */ 
+  vcopy( p_overlap, overlap_buf, 1, 1, nlong ); 
+  o_buf = overlap_buf; 
+ 
+
+  /* Separate action for each Block Type */ 
+   switch( block_type ) { 
+   case ONLY_LONG_WINDOW : 
+    imdct( p_in_data, transf_buf, 2*nlong ); 
+    vmult1( transf_buf, window_long_prev, transf_buf, nlong ); 
+    if (overlap_select != MNON_OVERLAPPED) {
+      vadd( transf_buf, o_buf, o_buf, 1, 1, 1, nlong ); 
+      vmult2( transf_buf+nlong, window_long+nlong-1, o_buf+nlong, nlong ); 
+    }
+    else { /* overlap_select == NON_OVERLAPPED */
+      vmult2( transf_buf+nlong, window_long+nlong-1, transf_buf+nlong, nlong );
+    }
+
+    break; 
+ 
+  case LONG_SHORT_WINDOW : 
+    imdct( p_in_data, transf_buf, 2*nlong ); 
+    vmult1( transf_buf, window_long_prev, transf_buf, nlong ); 
+    if (overlap_select != MNON_OVERLAPPED) {
+      vadd( transf_buf, o_buf, o_buf, 1, 1, 1, nlong ); 
+      vcopy( transf_buf+nlong, o_buf+nlong, 1, 1, nflat_ls ); 
+      vmult2( transf_buf+nlong+nflat_ls, window_short+nshort-1, o_buf+nlong+nflat_ls, nshort ); 
+      vcopy( &zero, o_buf+2*nlong-1, 0, -1, nflat_ls ); 
+    }
+    else { /* overlap_select == NON_OVERLAPPED */
+      vmult2( transf_buf+nlong+nflat_ls, window_short+nshort-1, transf_buf+nlong+nflat_ls, nshort ); 
+      vcopy( &zero, transf_buf+2*nlong-1, 0, -1, nflat_ls ); 
+    }
+    break; 
+    
+  case SHORT_LONG_WINDOW : 
+    imdct( p_in_data, transf_buf, 2*nlong ); 
+    vmult1( transf_buf+nflat_ls, window_short_prev_ptr, transf_buf+nflat_ls, nshort ); 
+    if (overlap_select != MNON_OVERLAPPED) {
+      vadd( transf_buf+nflat_ls, o_buf+nflat_ls, o_buf+nflat_ls, 1, 1, 1, nshort ); 
+      vcopy( transf_buf+nflat_ls+nshort, o_buf+nflat_ls+nshort, 1, 1, nflat_ls ); 
+      vmult2( transf_buf+nlong, window_long+nlong-1, o_buf+nlong, nlong ); 
+    }
+    else { /* overlap_select == NON_OVERLAPPED */
+      vcopy( &zero, transf_buf, 0, 1, nflat_ls);
+      vmult2( transf_buf+nlong, window_long+nlong-1, transf_buf+nlong, nlong);
+    }
+    break; 
+ 
+  case ONLY_SHORT_WINDOW : 
+    if (overlap_select != MNON_OVERLAPPED) {
+      fp = o_buf + nflat_ls; 
+    }
+    else { /* overlap_select == NON_OVERLAPPED */
+      fp = transf_buf;
+    }
+    for( k = transfak_ls-1; k-->= 0; ) { 
+      if (overlap_select != MNON_OVERLAPPED) {
+        imdct( p_in_data, transf_buf, 2*nshort ); 
+        vmult1( transf_buf, window_short_prev_ptr, transf_buf, nshort ); 
+        vadd( transf_buf, fp, fp, 1, 1, 1, nshort ); 
+        vmult2( transf_buf+nshort, window_short+nshort-1, fp+nshort, nshort ); 
+        p_in_data += nshort; 
+        fp        += nshort; 
+        window_short_prev_ptr = window_short; 
+      }
+      else { /* overlap_select == NON_OVERLAPPED */
+        imdct( p_in_data, fp, 2*nshort );
+        vmult1( fp, window_short_prev_ptr, fp, nshort ); 
+        vmult2( fp+nshort, window_short+nshort-1, fp+nshort, nshort ); 
+        p_in_data += nshort; 
+        fp        += 2*nshort;
+        window_short_prev_ptr = window_short; 
+      }
+    } 
+    vcopy( &zero, o_buf+2*nlong-1, 0, -1, nflat_ls ); 
+    break;     
+   default :
+	   break;
+    //CommonExit( 1, "mdct_synthesis: Unknown window type" ); 
+  } 
+ 
+  if (overlap_select != MNON_OVERLAPPED) {
+    vcopy( o_buf, p_out_data, 1, 1, nlong ); 
+  }
+  else { /* overlap_select == NON_OVERLAPPED */
+    vcopy( transf_buf, p_out_data, 1, 1, 2*nlong ); 
+  }
+  
+  /* save unused output data */ 
+  vcopy( o_buf+nlong, p_overlap, 1, 1, nlong ); 
+
+  /* --- Save old window shape --- */
+  wfun_select_prev = wfun_select;
+  
+} 
 
 /***********************************************************************************************/ 
--- a/is.c
+++ b/is.c
@@ -141,3 +141,62 @@
 	}  /* for (chanNum... */
 } 
 
+
+void ISReconstruct(double *spectrum[MAX_TIME_CHANNELS],   /* array of pointers to spectral data */
+				   Ch_Info *channelInfo,                  /* Pointer to Ch_Info */
+				   int sfb_offset_table[][MAX_SCFAC_BANDS+1],
+				   enum WINDOW_TYPE block_type[MAX_TIME_CHANNELS], /* Block type */
+				   AACQuantInfo* quantInfo,
+				   int numberOfChannels)                 /* Number of channels */
+{
+	int chanNum;
+	int sfbNum;
+	int lineNum,line_offset;
+	int startWindow,stopWindow,w;
+
+	/* Look for channel_pair_elements */
+	for (chanNum=0;chanNum<numberOfChannels;chanNum++) {
+		if (channelInfo[chanNum].present) {
+			if ((channelInfo[chanNum].cpe)&&(channelInfo[chanNum].ch_is_left)) {
+				int leftChan=chanNum;
+				int rightChan=channelInfo[chanNum].paired_ch;
+
+				/* If intensity stereo is present in right channel, reconstruct spectrum */
+				IS_Info *isInfo;
+				isInfo = &(channelInfo[rightChan].is_info);
+				if (isInfo->is_present) {
+					int numGroups = quantInfo[rightChan].num_window_groups;
+					int maxSfb = quantInfo[rightChan].max_sfb;
+					int g,b;
+
+					/* Currently, enable Intensity Stereo above band IS_MIN_BAND */
+					startWindow=0;
+					for (g=0;g<numGroups;g++) {
+
+						/* Enable IS bands */
+						int numWindows = quantInfo[leftChan].window_group_length[g];
+						stopWindow = startWindow + numWindows;
+						for (sfbNum=0;sfbNum<maxSfb;sfbNum++) {
+							b = g*maxSfb+sfbNum;
+							if (isInfo->is_used[b]) {
+								double scale = pow( 0.5, 0.25 * ((double)isInfo->fac[b])); 
+								scale = scale * (-2.0 * isInfo->sign[b] + 1.0);
+								for (w=startWindow;w<stopWindow;w++) {
+									line_offset = w*BLOCK_LEN_SHORT;
+									for (lineNum=sfb_offset_table[rightChan][sfbNum];
+									lineNum<sfb_offset_table[rightChan][sfbNum+1];
+									lineNum++) {
+										spectrum[rightChan][line_offset+lineNum] = 
+											spectrum[leftChan][line_offset+lineNum] * scale;
+									}
+								}
+							}  /* if (isInfo... */
+						}  /* for (sfbNum... */
+						startWindow = stopWindow;
+					}
+				}  /* if (block_type... */
+			}  /* if ((channelInfo... */
+		}  /* if (channelInfo... */
+	}  /* for (chanNum... */
+} 
+
--- a/is.h
+++ b/is.h
@@ -15,7 +15,7 @@
 #include "all.h"
 #include "tf_main.h"
 #include "aac_qc.h"
-#include "math.h"
+#include <math.h>
 
 #define IS_MIN_BAND_L 28
 #define IS_MIN_BAND_S 7
--- /dev/null
+++ b/nok_ltp_common.h
@@ -1,0 +1,92 @@
+/**************************************************************************
+
+This software module was originally developed by
+
+Mikko Suonio (Nokia)
+
+in the course of development of the MPEG-2 NBC/MPEG-4 Audio standard
+ISO/IEC 13818-7, 14496-1,2 and 3. This software module is an
+implementation of a part of one or more MPEG-2 NBC/MPEG-4 Audio tools
+as specified by the MPEG-2 NBC/MPEG-4 Audio standard. ISO/IEC gives
+users of the MPEG-2 NBC/MPEG-4 Audio standards free license to this
+software module or modifications thereof for use in hardware or
+software products claiming conformance to the MPEG-2 NBC/ MPEG-4 Audio
+standards. Those intending to use this software module in hardware or
+software products are advised that this use may infringe existing
+patents. The original developer of this software module and his/her
+company, the subsequent editors and their companies, and ISO/IEC have
+no liability for use of this software module or modifications thereof
+in an implementation. Copyright is not released for non MPEG-2
+NBC/MPEG-4 Audio conforming products. The original developer retains
+full right to use the code for his/her own purpose, assign or donate
+the code to a third party and to inhibit third party from using the
+code for non MPEG-2 NBC/MPEG-4 Audio conforming products. This
+copyright notice must be included in all copies or derivative works.
+
+Copyright (c) 1997.
+
+***************************************************************************/
+
+#ifndef _NOK_LTP_COMMON_H
+#define _NOK_LTP_COMMON_H
+
+/*
+  Macro:	MAX_SHORT_WINDOWS
+  Purpose:	Number of short windows in one long window.
+  Explanation:	-  */
+#ifndef MAX_SHORT_WINDOWS
+#define MAX_SHORT_WINDOWS NSHORT
+#endif
+
+/*
+  Macro:	MAX_SCFAC_BANDS
+  Purpose:	Maximum number of scalefactor bands in one frame.
+  Explanation:	-  */
+#ifndef MAX_SCFAC_BANDS
+#define MAX_SCFAC_BANDS MAXBANDS
+#endif
+
+/*
+  Macro:	BLOCK_LEN_LONG
+  Purpose:	Length of one long window
+  Explanation:	-  */
+#ifndef BLOCK_LEN_LONG
+#define BLOCK_LEN_LONG LN2
+#endif
+
+
+/*
+  Macro:	NOK_MAX_BLOCK_LEN_LONG
+  Purpose:	Informs the routine of the maximum block size used.
+  Explanation:	This is needed since the TwinVQ long window
+  		is different from the AAC long window.  */
+#define	NOK_MAX_BLOCK_LEN_LONG BLOCK_LEN_LONG //(2 * BLOCK_LEN_LONG) 
+
+/*
+  Macro:	NOK_LT_BLEN
+  Purpose:	Length of the history buffer.
+  Explanation:	Has to hold two long windows of time domain data.  */
+#ifndef	NOK_LT_BLEN
+#define NOK_LT_BLEN (4 * NOK_MAX_BLOCK_LEN_LONG)
+#endif
+
+/*
+  Type:		NOK_LT_PRED_STATUS
+  Purpose:	Type of the struct holding the LTP encoding parameters.
+  Explanation:	-  */
+typedef struct
+  {
+    short buffer[NOK_LT_BLEN];
+    double pred_mdct[2 * NOK_MAX_BLOCK_LEN_LONG];
+    int weight_idx;
+    double weight;
+    int sbk_prediction_used[MAX_SHORT_WINDOWS];
+    int sfb_prediction_used[MAX_SCFAC_BANDS];
+    int *delay;
+    int global_pred_flag;
+    int side_info;
+  }
+NOK_LT_PRED_STATUS;
+
+
+#endif /* _NOK_LTP_COMMON_H */
--- /dev/null
+++ b/nok_ltp_common_internal.h
@@ -1,0 +1,115 @@
+/**************************************************************************
+
+This software module was originally developed by
+
+Mikko Suonio (Nokia)
+
+in the course of development of the MPEG-2 NBC/MPEG-4 Audio standard
+ISO/IEC 13818-7, 14496-1,2 and 3. This software module is an
+implementation of a part of one or more MPEG-2 NBC/MPEG-4 Audio tools
+as specified by the MPEG-2 NBC/MPEG-4 Audio standard. ISO/IEC gives
+users of the MPEG-2 NBC/MPEG-4 Audio standards free license to this
+software module or modifications thereof for use in hardware or
+software products claiming conformance to the MPEG-2 NBC/ MPEG-4 Audio
+standards. Those intending to use this software module in hardware or
+software products are advised that this use may infringe existing
+patents. The original developer of this software module and his/her
+company, the subsequent editors and their companies, and ISO/IEC have
+no liability for use of this software module or modifications thereof
+in an implementation. Copyright is not released for non MPEG-2
+NBC/MPEG-4 Audio conforming products. The original developer retains
+full right to use the code for his/her own purpose, assign or donate
+the code to a third party and to inhibit third party from using the
+code for non MPEG-2 NBC/MPEG-4 Audio conforming products. This
+copyright notice must be included in all copies or derivative works.
+
+Copyright (c) 1997.
+
+***************************************************************************/
+
+#ifndef _NOK_LTP_COMMON_INTERNAL_H
+#define _NOK_LTP_COMMON_INTERNAL_H
+
+
+/*
+  Purpose:      Number of LTP coefficients. */
+#define LPC 1
+
+/*
+  Purpose:      Maximum LTP lag.  */
+#define DELAY 2048
+
+/*
+  Purpose:	Length of the bitstream element ltp_data_present.  */
+#define	LEN_LTP_DATA_PRESENT 1
+
+/*
+  Purpose:	Length of the bitstream element ltp_lag.  */
+#define	LEN_LTP_LAG 11
+
+/*
+  Purpose:	Length of the bitstream element ltp_coef.  */
+#define	LEN_LTP_COEF 3
+
+/*
+  Purpose:	Length of the bitstream element ltp_short_used.  */
+#define	LEN_LTP_SHORT_USED 1
+
+/*
+  Purpose:	Length of the bitstream element ltp_short_lag_present.  */
+#define	LEN_LTP_SHORT_LAG_PRESENT 1
+
+/*
+  Purpose:	Length of the bitstream element ltp_short_lag.  */
+#define	LEN_LTP_SHORT_LAG 5
+
+/*
+  Purpose:	Offset of the lags written in the bitstream.  */
+#define	NOK_LTP_LAG_OFFSET 16
+
+/*
+  Purpose:	Length of the bitstream element ltp_long_used.  */
+#define	LEN_LTP_LONG_USED 1
+
+/*
+  Purpose:	Upper limit for the number of scalefactor bands
+   		which can use lt prediction with long windows.
+  Explanation:	Bands 0..NOK_MAX_LT_PRED_SFB-1 can use lt prediction.  */
+#define	NOK_MAX_LT_PRED_LONG_SFB 40
+
+/*
+  Purpose:	Upper limit for the number of scalefactor bands
+   		which can use lt prediction with short windows.
+  Explanation:	Bands 0..NOK_MAX_LT_PRED_SFB-1 can use lt prediction.  */
+#define	NOK_MAX_LT_PRED_SHORT_SFB 13
+
+/*
+   Purpose:      Buffer offset to maintain block alignment.
+   Explanation:  This is only used for a short window sequence.  */
+#define SHORT_SQ_OFFSET (BLOCK_LEN_LONG-(BLOCK_LEN_SHORT*4+BLOCK_LEN_SHORT/2))
+
+/*
+  Purpose:	Number of codes for LTP weight. */
+#define CODESIZE 8
+
+/*
+   Purpose:      Float type for external data
+   Explanation:  - */
+typedef double float_ext;
+
+/*
+  Purpose:	Codebook for LTP weight coefficients.  */
+static double codebook[CODESIZE] =
+{
+  0.570829,
+  0.696616,
+  0.813004,
+  0.911304,
+  0.984900,
+  1.067894,
+  1.194601,
+  1.369533
+};
+
+
+#endif /* _NOK_LTP_COMMON_INTERNAL_H */
--- /dev/null
+++ b/nok_ltp_enc.c
@@ -1,0 +1,502 @@
+/**************************************************************************
+
+This software module was originally developed by
+Nokia in the course of development of the MPEG-2 AAC/MPEG-4 
+Audio standard ISO/IEC13818-7, 14496-1, 2 and 3.
+This software module is an implementation of a part
+of one or more MPEG-2 AAC/MPEG-4 Audio tools as specified by the
+MPEG-2 aac/MPEG-4 Audio standard. ISO/IEC  gives users of the
+MPEG-2aac/MPEG-4 Audio standards free license to this software module
+or modifications thereof for use in hardware or software products
+claiming conformance to the MPEG-2 aac/MPEG-4 Audio  standards. Those
+intending to use this software module in hardware or software products
+are advised that this use may infringe existing patents. The original
+developer of this software module, the subsequent
+editors and their companies, and ISO/IEC have no liability for use of
+this software module or modifications thereof in an
+implementation. Copyright is not released for non MPEG-2 aac/MPEG-4
+Audio conforming products. The original developer retains full right to
+use the code for the developer's own purpose, assign or donate the code to a
+third party and to inhibit third party from using the code for non
+MPEG-2 aac/MPEG-4 Audio conforming products. This copyright notice
+must be included in all copies or derivative works.
+Copyright (c)1997.  
+
+***************************************************************************/
+/**************************************************************************
+  nok_ltp_enc.c  -  Performs Long Term Prediction for the MPEG-4 
+                    T/F Encoder.
+  Author(s): Mikko Suonio, Juha Ojanpera
+  	     Nokia Research Center, Speech and Audio Systems, 1997.
+  *************************************************************************/
+
+
+/**************************************************************************
+  Version Control Information			Method: CVS
+  Identifiers:
+  $Revision: 1.1 $
+  $Date: 2000/01/05 21:41:37 $ (check in)
+  $Author: menno $
+  *************************************************************************/
+
+
+/**************************************************************************
+  External Objects Needed
+  *************************************************************************/
+/*
+  Standard library declarations.  */
+#include <stdio.h>
+
+/*
+  Interface to related modules.  */
+#include "tf_main.h"
+#include "interface.h"
+#include "bitstream.h"
+#include "nok_ltp_common.h"
+#include "nok_pitch.h"
+
+/**************************************************************************
+  External Objects Provided
+  *************************************************************************/
+#include "nok_ltp_enc.h"
+
+/**************************************************************************
+  Internal Objects
+  *************************************************************************/
+#include "nok_ltp_common_internal.h"
+
+static short double_to_int (double sig_in);
+
+
+/**************************************************************************
+  Title:	nok_init_lt_pred
+
+  Purpose:	Initialize the history buffer for long term prediction
+
+  Usage:        nok_init_lt_pred (lt_status)
+
+  Input:	lt_status
+  			- buffer: history buffer
+                        - pred_mdct: prediction transformed to frequency 
+                          domain
+                        - weight_idx : 
+                          3 bit number indicating the LTP coefficient in 
+                          the codebook 
+                        - sbk_prediction_used: 
+                          1 bit for each subblock indicating wheather
+                          LTP is used in that subblock 
+                        - sfb_prediction_used:
+                          1 bit for each scalefactor band (sfb) where LTP 
+                          can be used indicating whether LTP is switched 
+                          on (1) /off (0) in that sfb.
+                        - delay: LTP lag
+                        - side_info: LTP side information
+
+  Output:	lt_status
+  			- buffer: filled with 0
+                        - pred_mdct: filled with 0
+                        - weight_idx : filled with 0
+                        - sbk_prediction_used: filled with 0
+                        - sfb_prediction_used: filled with 0
+                        - delay: filled with 0
+                        - side_info: filled with 1
+
+  References:	-
+
+  Explanation:	-
+
+  Author(s):	Mikko Suonio
+  *************************************************************************/
+void
+nok_init_lt_pred (NOK_LT_PRED_STATUS *lt_status)
+{
+	int i;
+
+	for (i = 0; i < NOK_LT_BLEN; i++)
+		lt_status->buffer[i] = 0;
+
+	for (i = 0; i < NOK_MAX_BLOCK_LEN_LONG; i++)
+		lt_status->pred_mdct[i] = 0;
+
+	lt_status->weight_idx = 0;
+	for(i=0; i<NSHORT; i++)
+		lt_status->sbk_prediction_used[i] = lt_status->delay[i] = 0;
+
+	for(i=0; i<MAX_SCFAC_BANDS; i++)
+		lt_status->sfb_prediction_used[i] = 0;
+
+	lt_status->side_info = LEN_LTP_DATA_PRESENT;
+}
+
+
+/**************************************************************************
+  Title:	nok_ltp_enc
+
+  Purpose:      Performs long term prediction.
+
+  Usage:	nok_ltp_enc(p_spectrum, p_time_signal, win_type, win_shape, 
+                            sfb_offset, num_of_sfb, lt_status, buffer_update)
+
+  Input:        p_spectrum    - spectral coefficients
+                p_time_signal - time domain input samples
+                win_type      - window sequence (frame, block) type
+                win_shape     - shape of the mdct window
+                sfb_offset    - scalefactor band boundaries
+                num_of_sfb    - number of scalefactor bands in each block
+
+  Output:	p_spectrum    - residual spectrum
+                lt_status     - buffer: history buffer
+                              - pred_mdct:prediction transformed to frequency domain
+                                for subsequent use
+                              - weight_idx : 
+                                3 bit number indicating the LTP coefficient in 
+                                the codebook 
+                              - sbk_prediction_used: 
+                                1 bit for each subblock indicating wheather
+                                LTP is used in that subblock 
+                              - sfb_prediction_used:
+                                1 bit for each scalefactor band (sfb) where LTP 
+                                can be used indicating whether LTP is switched 
+                                on (1) /off (0) in that sfb.
+                              - delay: LTP lag
+                              - side_info: LTP side information
+                        
+  References:	1.) estimate_delay in nok_pitch.c
+                2.) pitch in nok_pitch.c
+                3.) buffer2freq
+                4.) snr_pred in nok_pitch.c
+                5.) freq2buffer
+                6.) double_to_int
+
+  Explanation:  -
+
+  Author(s):	Juha Ojanpera
+  *************************************************************************/
+
+int
+nok_ltp_enc(double *p_spectrum, double *p_time_signal, enum WINDOW_TYPE win_type, 
+            Window_shape win_shape, int block_size_long, int block_size_medium,
+	    int block_size_short, int *sfb_offset, int num_of_sfb,
+            NOK_LT_PRED_STATUS *lt_status)
+{
+    int i;
+    int last_band;
+	double num_bit[MAX_SHORT_WINDOWS];
+    double predicted_samples[2 * NOK_MAX_BLOCK_LEN_LONG];
+
+    lt_status->global_pred_flag = 0;
+    lt_status->side_info = 1;
+
+    switch(win_type)
+    {
+	case ONLY_LONG_WINDOW:
+	case LONG_SHORT_WINDOW:
+	case SHORT_LONG_WINDOW:
+		last_band = (num_of_sfb < NOK_MAX_LT_PRED_LONG_SFB) ? num_of_sfb : NOK_MAX_LT_PRED_LONG_SFB;
+
+		lt_status->delay[0] = estimate_delay (p_time_signal, lt_status->buffer, 2 * block_size_long);
+
+//		fprintf(stderr, "(LTP) lag : %i ", lt_status->delay[0]);
+
+		pitch (p_time_signal, predicted_samples, lt_status->buffer, 
+			&lt_status->weight_idx, lt_status->delay[0], 2 * block_size_long);
+
+		/* Transform prediction to frequency domain and save it for subsequent use. */
+		buffer2freq (predicted_samples, lt_status->pred_mdct, NULL, win_type,
+			win_shape, block_size_long, block_size_medium,
+			block_size_short, MNON_OVERLAPPED);
+
+		lt_status->side_info = LEN_LTP_DATA_PRESENT + last_band + LEN_LTP_LAG + LEN_LTP_COEF;
+    
+		num_bit[0] = snr_pred (p_spectrum, lt_status->pred_mdct,
+			lt_status->sfb_prediction_used, sfb_offset, 
+			win_type, lt_status->side_info, last_band);
+
+//		fprintf(stderr, " bit gain : %f\n", num_bit[0]);
+
+		lt_status->global_pred_flag = (num_bit[0] == 0.0) ? 0 : 1;
+
+		if(lt_status->global_pred_flag)
+			for (i = 0; i < sfb_offset[last_band]; i++)
+				p_spectrum[i] -= lt_status->pred_mdct[i];
+		else
+			lt_status->side_info = 1;
+        break;
+
+	case ONLY_SHORT_WINDOW:
+		break;
+    }
+    
+    return (lt_status->global_pred_flag);
+}
+
+
+/**************************************************************************
+  Title:	nok_ltp_reconstruct
+
+  Purpose:      Updates LTP history buffer.
+
+  Usage:	nok_ltp_reconstruct(p_spectrum, win_type,  win_shape, 
+                                    block_size_long, block_size_medium,
+                                    block_size_short, sfb_offset, 
+                                    num_of_sfb, lt_status)
+
+  Input:        p_spectrum    - reconstructed spectrum
+                win_type      - window sequence (frame, block) type
+                win_shape     - shape of the mdct window
+                sfb_offset    - scalefactor band boundaries
+                num_of_sfb    - number of scalefactor bands in each block
+
+  Output:	p_spectrum    - reconstructed spectrum
+                lt_status     - buffer: history buffer
+                           
+  References:	1.) buffer2freq
+                2.) freq2buffer
+                3.) double_to_int
+
+  Explanation:  -
+
+  Author(s):	Juha Ojanpera
+  *************************************************************************/
+
+void
+nok_ltp_reconstruct(double *p_spectrum, enum WINDOW_TYPE win_type, 
+                    Window_shape win_shape,
+		    int block_size_long, int block_size_medium,
+                    int block_size_short, int *sfb_offset, int num_of_sfb,
+                    NOK_LT_PRED_STATUS *lt_status)
+{
+	int i, j, last_band;
+	double predicted_samples[2 * NOK_MAX_BLOCK_LEN_LONG];
+	double overlap_buffer[2 * NOK_MAX_BLOCK_LEN_LONG];
+
+    
+	switch(win_type)
+	{
+	case ONLY_LONG_WINDOW:
+	case LONG_SHORT_WINDOW:
+	case SHORT_LONG_WINDOW:
+		last_band = (num_of_sfb < NOK_MAX_LT_PRED_LONG_SFB) ? num_of_sfb : NOK_MAX_LT_PRED_LONG_SFB;
+		
+		if(lt_status->global_pred_flag)
+			for (i = 0; i < sfb_offset[last_band]; i++)
+				p_spectrum[i] += lt_status->pred_mdct[i];
+
+		/* Finally update the time domain history buffer. */
+		freq2buffer (p_spectrum, predicted_samples, overlap_buffer, win_type,
+			block_size_long, block_size_medium, block_size_short,
+			win_shape, MNON_OVERLAPPED);
+		
+		for (i = 0; i < NOK_LT_BLEN - block_size_long; i++)
+			lt_status->buffer[i] = lt_status->buffer[i + block_size_long];
+
+		j = NOK_LT_BLEN - 2 * block_size_long;
+		for (i = 0; i < block_size_long; i++)
+		{
+			lt_status->buffer[i + j] =
+				double_to_int (predicted_samples[i] + lt_status->buffer[i + j]);
+			lt_status->buffer[NOK_LT_BLEN - block_size_long + i] =
+				double_to_int (predicted_samples[i + block_size_long]);
+		}
+		break;
+	
+    case ONLY_SHORT_WINDOW:
+#if 0
+		for (i = 0; i < NOK_LT_BLEN - block_size_long; i++)
+			lt_status->buffer[i] = lt_status->buffer[i + block_size_long];
+		
+		for (i = NOK_LT_BLEN - block_size_long; i < NOK_LT_BLEN; i++)
+			lt_status->buffer[i] = 0;
+
+		for (i = 0; i < block_size_long; i++)
+			overlap_buffer[i] = 0;
+
+		/* Finally update the time domain history buffer. */
+		freq2buffer (p_spectrum, predicted_samples, overlap_buffer, win_type, block_size_long, 
+			block_size_medium, block_size_short, win_shape, MNON_OVERLAPPED);
+
+		for(sw = 0; sw < MAX_SHORT_WINDOWS; sw++)
+		{
+			i = NOK_LT_BLEN - 2 * block_size_long + SHORT_SQ_OFFSET + sw * block_size_short;
+			for (j = 0; j < 2 * block_size_short; j++)
+				lt_status->buffer[i + j] = double_to_int (predicted_samples[sw * block_size_short * 2 + j] + 
+				lt_status->buffer[i + j]);
+		}
+#endif
+		break;
+	}
+
+	return;
+}                      
+
+
+/**************************************************************************
+  Title:	nok_ltp_encode
+
+  Purpose:      Writes LTP parameters to the bit stream.
+
+  Usage:	nok_ltp_encode (bs, win_type, num_of_sfb, lt_status)
+
+  Input:        bs         - bit stream
+                win_type   - window sequence (frame, block) type
+                num_of_sfb - number of scalefactor bands
+                lt_status  - side_info:
+			     1, if prediction not used in this frame
+			     >1 otherwise
+                           - weight_idx : 
+                             3 bit number indicating the LTP coefficient in 
+                             the codebook 
+                           - sfb_prediction_used:
+                             1 bit for each scalefactor band (sfb) where LTP 
+                             can be used indicating whether LTP is switched 
+                             on (1) /off (0) in that sfb.
+                           - delay: LTP lag
+
+  Output:	-
+
+  References:	1.) BsPutBit
+
+  Explanation:  -
+
+  Author(s):	Juha Ojanpera
+  *************************************************************************/
+
+int
+nok_ltp_encode (BsBitStream *bs, enum WINDOW_TYPE win_type, int num_of_sfb, 
+                NOK_LT_PRED_STATUS *lt_status, int write_flag)
+{
+	int i, last_band;
+	int first_subblock;
+	int prev_subblock;
+	int bit_count = 0;
+
+
+	bit_count += 1;
+	
+	if (lt_status->side_info > 1)
+	{
+		if(write_flag)
+			BsPutBit (bs, 1, 1);    	/* LTP used */
+
+		switch(win_type)
+		{
+		case ONLY_LONG_WINDOW:
+		case LONG_SHORT_WINDOW:
+		case SHORT_LONG_WINDOW:
+			bit_count += LEN_LTP_LAG;
+			bit_count += LEN_LTP_COEF;
+			if(write_flag)
+			{
+				BsPutBit (bs, lt_status->delay[0], LEN_LTP_LAG);
+				BsPutBit (bs, lt_status->weight_idx,  LEN_LTP_COEF);
+			}
+
+			last_band = (num_of_sfb < NOK_MAX_LT_PRED_LONG_SFB) ? num_of_sfb : NOK_MAX_LT_PRED_LONG_SFB;
+			bit_count += last_band;
+			if(write_flag)
+			{
+				for (i = 0; i < last_band; i++)
+					BsPutBit (bs, lt_status->sfb_prediction_used[i], LEN_LTP_LONG_USED);
+			}
+			break;
+			
+		case ONLY_SHORT_WINDOW:
+			for(i=0; i < MAX_SHORT_WINDOWS; i++)
+			{
+				if(lt_status->sbk_prediction_used[i])
+				{
+					first_subblock = i;
+					break;
+				}
+			}
+			bit_count += LEN_LTP_LAG;
+			bit_count += LEN_LTP_COEF;
+			
+			if(write_flag)
+			{
+				BsPutBit (bs, lt_status->delay[first_subblock], LEN_LTP_LAG);
+				BsPutBit (bs, lt_status->weight_idx,  LEN_LTP_COEF);
+			}
+
+			prev_subblock = first_subblock;
+			for(i = 0; i < MAX_SHORT_WINDOWS; i++)
+			{
+				bit_count += LEN_LTP_SHORT_USED;
+				if(write_flag)
+					BsPutBit (bs, lt_status->sbk_prediction_used[i], LEN_LTP_SHORT_USED);
+
+				if(lt_status->sbk_prediction_used[i])
+				{
+					if(i > first_subblock)
+					{
+						int diff;
+						
+						diff = lt_status->delay[prev_subblock] - lt_status->delay[i];
+						if(diff)
+						{
+							bit_count += 1;
+							bit_count += LEN_LTP_SHORT_LAG;
+							if(write_flag)
+							{
+								BsPutBit (bs, 1, 1);
+								BsPutBit (bs, diff + NOK_LTP_LAG_OFFSET, LEN_LTP_SHORT_LAG);
+							}
+						}
+						else
+						{
+							bit_count += 1;
+							if(write_flag)
+								BsPutBit (bs, 0, 1);
+						}
+					}
+				}
+			}  
+			break;
+
+		default:
+			//        CommonExit(1, "nok_ltp_encode : unsupported window sequence %i", win_type);
+			break;
+		}
+	}
+	else
+		if(write_flag)
+			BsPutBit (bs, 0, 1);    	/* LTP not used */
+
+	return (bit_count);
+}
+
+
+/**************************************************************************
+  Title:	double_to_int
+
+  Purpose:      Converts floating point format to integer (16-bit).
+
+  Usage:	y = double_to_int(sig_in)
+
+  Input:	sig_in  - floating point number
+
+  Output:	y  - integer number
+
+  References:	-
+
+  Explanation:  -
+
+  Author(s):	Juha Ojanpera
+  *************************************************************************/
+
+static short
+double_to_int (double sig_in)
+{
+	short sig_out;
+
+	if (sig_in > 32767)
+		sig_out = 32767;
+	else if (sig_in < -32768)
+		sig_out = -32768;
+	else if (sig_in > 0.0)
+		sig_out = (short) (sig_in + 0.5);
+	else if (sig_in <= 0.0)
+		sig_out = (short) (sig_in - 0.5);
+
+	return (sig_out);
+}
--- /dev/null
+++ b/nok_ltp_enc.h
@@ -1,0 +1,50 @@
+/**************************************************************************
+
+This software module was originally developed by
+Nokia in the course of development of the MPEG-2 AAC/MPEG-4 
+Audio standard ISO/IEC13818-7, 14496-1, 2 and 3.
+This software module is an implementation of a part
+of one or more MPEG-2 AAC/MPEG-4 Audio tools as specified by the
+MPEG-2 aac/MPEG-4 Audio standard. ISO/IEC  gives users of the
+MPEG-2aac/MPEG-4 Audio standards free license to this software module
+or modifications thereof for use in hardware or software products
+claiming conformance to the MPEG-2 aac/MPEG-4 Audio  standards. Those
+intending to use this software module in hardware or software products
+are advised that this use may infringe existing patents. The original
+developer of this software module, the subsequent
+editors and their companies, and ISO/IEC have no liability for use of
+this software module or modifications thereof in an
+implementation. Copyright is not released for non MPEG-2 aac/MPEG-4
+Audio conforming products. The original developer retains full right to
+use the code for the developer's own purpose, assign or donate the code to a
+third party and to inhibit third party from using the code for non
+MPEG-2 aac/MPEG-4 Audio conforming products. This copyright notice
+must be included in all copies or derivative works.
+Copyright (c)1997.  
+
+***************************************************************************/
+
+#ifndef _NOK_LTP_ENC_H
+#define _NOK_LTP_ENC_H
+
+#include "nok_ltp_common.h"
+
+extern void nok_init_lt_pred (NOK_LT_PRED_STATUS * lt_status);
+
+extern int nok_ltp_enc(double *p_spectrum, double *p_time_signal,
+		       enum WINDOW_TYPE win_type, Window_shape win_shape,
+		       int block_size_long, int block_size_medium,
+		       int block_size_short, int *sfb_offset, int num_of_sfb,
+		       NOK_LT_PRED_STATUS *lt_status);
+
+extern void nok_ltp_reconstruct(double *p_spectrum, enum WINDOW_TYPE win_type, 
+                                Window_shape win_shape, int block_size_long, 
+                                int block_size_medium, int block_size_short, 
+                                int *sfb_offset, int num_of_sfb,
+                                NOK_LT_PRED_STATUS *lt_status);
+
+extern int nok_ltp_encode (BsBitStream *bs, enum WINDOW_TYPE win_type, int num_of_sfb, 
+                           NOK_LT_PRED_STATUS *lt_status, int write_flag);
+
+#endif /* not defined _NOK_LTP_ENC_H */
+
--- /dev/null
+++ b/nok_pitch.c
@@ -1,0 +1,499 @@
+/**************************************************************************
+
+This software module was originally developed by
+Nokia in the course of development of the MPEG-2 AAC/MPEG-4 
+Audio standard ISO/IEC13818-7, 14496-1, 2 and 3.
+This software module is an implementation of a part
+of one or more MPEG-2 AAC/MPEG-4 Audio tools as specified by the
+MPEG-2 aac/MPEG-4 Audio standard. ISO/IEC  gives users of the
+MPEG-2aac/MPEG-4 Audio standards free license to this software module
+or modifications thereof for use in hardware or software products
+claiming conformance to the MPEG-2 aac/MPEG-4 Audio  standards. Those
+intending to use this software module in hardware or software products
+are advised that this use may infringe existing patents. The original
+developer of this software module, the subsequent
+editors and their companies, and ISO/IEC have no liability for use of
+this software module or modifications thereof in an
+implementation. Copyright is not released for non MPEG-2 aac/MPEG-4
+Audio conforming products. The original developer retains full right to
+use the code for the developer's own purpose, assign or donate the code to a
+third party and to inhibit third party from using the code for non
+MPEG-2 aac/MPEG-4 Audio conforming products. This copyright notice
+must be included in all copies or derivative works.
+Copyright (c)1997.  
+
+***************************************************************************/
+/**************************************************************************
+  nok_pitch.c  - Long Term Prediction (LTP) subroutines.
+ 
+  Author(s): Juha Ojanpera, Lin Yin
+  	     Nokia Research Center, Speech and Audio Systems.
+  *************************************************************************/
+
+
+/**************************************************************************
+  Version Control Information			Method: CVS
+  Identifiers:
+  $Revision: 1.1 $
+  $Date: 2000/01/05 21:41:37 $ (check in)
+  $Author: menno $
+  *************************************************************************/
+
+
+/**************************************************************************
+  External Objects Needed
+  *************************************************************************/
+/*
+  Standard library declarations.  */
+#include <math.h>
+#include <stdio.h>
+
+/*
+  Interface to related modules.  */
+#include "tf_main.h"
+#include "block.h"
+#include "interface.h"
+#include "nok_ltp_common.h"
+#include "nok_ltp_enc.h"
+
+
+/*************************************************************************
+  External Objects Provided
+  *************************************************************************/
+#include "nok_pitch.h"
+
+
+
+/**************************************************************************
+  Internal Objects
+  *************************************************************************/
+#include "nok_ltp_common_internal.h"
+
+static void lnqgj (double (*a)[LPC + 1]);
+
+static void w_quantize (double *freq, int *ltp_idx);
+
+/**************************************************************************
+  Title:	snr_pred
+
+  Purpose:	Determines is it feasible to employ long term prediction in 
+                the encoder.
+
+  Usage:        y = snr_pred(mdct_in, mdct_pred, sfb_flag, sfb_offset, 
+                             block_type, side_info, num_of_sfb)
+
+  Input:        mdct_in     -  spectral coefficients
+                mdct_pred   -  predicted spectral coefficients
+                sfb_flag    -  an array of flags indicating whether LTP is 
+                               switched on (1) /off (0). One bit is reseved 
+                               for each sfb.
+                sfb_offset  -  scalefactor boundaries
+                block_type  -  window sequence type
+		side_info   -  LTP side information
+		num_of_sfb  -  number of scalefactor bands
+
+  Output:	y - number of bits saved by using long term prediction     
+
+  References:	-
+
+  Explanation:	-
+
+  Author(s):	Juha Ojanpera, Lin Yin
+  *************************************************************************/
+
+double
+snr_pred (double *mdct_in,
+	  double *mdct_pred,
+	  int *sfb_flag,
+	  int *sfb_offset,
+	  enum WINDOW_TYPE block_type,
+	  int side_info,
+	  int num_of_sfb
+)
+{
+	int i, j, flen;
+	double snr_limit;
+	double num_bit, snr[NSFB_LONG];
+	double temp1, temp2;
+	double energy[BLOCK_LEN_LONG], snr_p[BLOCK_LEN_LONG];
+
+
+	if (block_type != ONLY_SHORT_WINDOW)
+    {
+		flen = BLOCK_LEN_LONG;
+		snr_limit = 1.e-30;
+    }
+	else
+    {
+		flen = BLOCK_LEN_SHORT;
+		snr_limit = 1.e-20;
+    }
+
+	for (i = 0; i < flen; i++)
+    {
+		energy[i] = mdct_in[i] * mdct_in[i];
+		snr_p[i] = (mdct_in[i] - mdct_pred[i]) * (mdct_in[i] - mdct_pred[i]);
+    }
+
+	for (i = 0; i < num_of_sfb; i++)
+    {
+		temp1 = 0.0;
+		temp2 = 0.0;
+		for (j = sfb_offset[i]; j < sfb_offset[i + 1]; j++)
+		{
+			temp1 += energy[j];
+			temp2 += snr_p[j];
+		}
+
+		if (temp2 < snr_limit)
+			temp2 = snr_limit;
+
+		if (temp1 > 1.e-20)
+			snr[i] = /*-*/10. * log10 (temp2 / temp1);
+		else
+			snr[i] = 0.0;
+
+		sfb_flag[i] = 1;
+
+		if (block_type != ONLY_SHORT_WINDOW)
+		{
+			if (snr[i] <= 0.0)
+			{
+				sfb_flag[i] = 0;
+				for (j = sfb_offset[i]; j < sfb_offset[i + 1]; j++)
+					mdct_pred[j] = 0.0;
+			}
+		}
+    }
+
+	num_bit = 0.0;
+	for (i = 0; i < num_of_sfb; i++)
+		num_bit += snr[i] / 6. * (sfb_offset[i + 1] - sfb_offset[i]);
+
+	if (num_bit < side_info)
+    {
+		num_bit = 0.0;
+		for (j = 0; j < flen; j++)
+			mdct_pred[j] = 0.0;
+		for (i = 0; i < num_of_sfb; i++)
+			sfb_flag[i] = 0;
+    }
+	else
+		num_bit -= side_info;
+
+
+	return (num_bit);
+}
+
+
+/**************************************************************************
+  Title:	prediction
+
+  Purpose:	Predicts current frame from the past output samples.
+
+  Usage:        prediction(buffer, predicted_samples, weight, delay, flen)
+
+  Input:	buffer	          - past reconstructed output samples
+                weight            - LTP gain (scaling factor)
+		delay             - LTP lag (optimal delay)
+                flen              - length of the frame
+
+  Output:	predicted_samples - predicted samples
+
+  References:	-
+
+  Explanation:  -
+
+  Author(s):	Juha Ojanpera, Lin Yin
+  *************************************************************************/
+
+void
+prediction (short *buffer,
+	    double *predicted_samples,
+	    double *weight,
+	    int delay,
+	    int flen
+)
+{
+
+  int i, j;
+  int offset;
+
+  for (i = 0; i < flen; i++)
+    {
+      offset = NOK_LT_BLEN - flen - delay + i - (LPC - 1) / 2;
+      predicted_samples[i] = 0.0;
+      for (j = 0; j < LPC; j++)
+	predicted_samples[i] += weight[j] * buffer[offset + j];
+    }
+
+}
+
+
+/**************************************************************************
+  Title:	estimate_delay
+
+  Purpose:      Estimates optimal delay between current frame and past
+                reconstructed output samples. The lag between 0...DELAY-1 
+                with the maximum correlation becomes the LTP lag (or delay).
+
+  Usage:        y = estimate_delay(sb_samples, x_buffer, flen)
+
+  Input:	sb_samples - current frame
+                x_buffer   - past reconstructed output samples
+                flen       - length of the frame
+
+  Output:	y  - LTP lag
+
+  References:	-
+
+  Explanation:	-
+
+  Author(s):	Juha Ojanpera, Lin Yin
+  *************************************************************************/
+
+int estimate_delay (double *sb_samples,
+					short *x_buffer,
+					int flen
+					)
+{
+#if 0
+	int i, j;
+	int delay;
+	double corr[DELAY];
+	double p_max, energy;
+
+	for (i = 0; i < DELAY; i++)
+	{
+		energy = 0.0;
+		corr[i] = 0.0;
+		for (j = 0; j < flen; j++)
+		{
+			corr[i] += x_buffer[NOK_LT_BLEN - i - j - 1] * sb_samples[flen - j - 1];
+			energy += x_buffer[NOK_LT_BLEN - i - j - 1] * x_buffer[NOK_LT_BLEN - i - j - 1];
+		}
+
+		if (energy != 0.0)
+			corr[i] = corr[i] / sqrt (energy);
+		else
+			corr[i] = 0.0;
+    }
+
+	p_max = 0.0;
+	delay = 0;
+	for (i = 0; i < DELAY; i++)
+		if (p_max < corr[i])
+		{
+			p_max = corr[i];
+			delay = i;
+		}
+
+	if (delay < (LPC - 1) / 2)
+		delay = (LPC - 1) / 2;
+
+	/*
+	fprintf(stdout, "delay : %d ... ", delay);
+	*/
+
+	return delay;
+#else
+	static int delay = 992;
+	delay += 32;
+	if (delay == 2048)
+		delay = 1024;
+	return delay;
+#endif
+}
+
+
+/**************************************************************************
+  Title:	pitch
+
+  Purpose:	Calculates LTP gains, quantizes them and finally scales
+                the past output samples (from 'delay' to 'delay'+'flen') 
+                to get the predicted frame.
+
+  Usage:	pitch(sb_samples, sb_samples_pred, x_buffer, ltp_coef, 
+                      delay, flen)
+
+  Input:        sb_samples      - current frame
+                x_buffer        - past reconstructed output samples
+                delay           - LTP lag
+                flen            - length of the frame
+
+  Output:	sb_samples_pred - predicted frame
+                ltp_coef        - indices of the LTP gains
+
+  References:	1.)  lnqgj
+                2.)  w_quantize
+                3.)  prediction
+
+  Explanation:	-
+
+  Author(s):	Juha Ojanpera, Lin Yin
+  *************************************************************************/
+
+void
+pitch (double *sb_samples,
+       double *sb_samples_pred,
+       short *x_buffer,
+       int *ltp_coef,
+       int delay,
+       int flen
+)
+{
+	int i, k, j;
+	int offset1, offset2;
+	double weight[LPC];
+	double r[LPC][LPC + 1];
+
+
+	for (i = 0; i < LPC; i++)
+		for (j = 0; j < LPC; j++)
+		{
+			offset1 = NOK_LT_BLEN - flen - delay + i - (LPC - 1) / 2;
+			offset2 = NOK_LT_BLEN - flen - delay + j - (LPC - 1) / 2;
+			r[i][j] = 0.0;
+			for (k = 0; k < flen; k++)
+				r[i][j] += x_buffer[offset1 + k] * x_buffer[offset2 + k];
+		}
+	for (i = 0; i < LPC; i++)
+		r[i][i] = 1.010 * r[i][i];
+
+	for (i = 0; i < LPC; i++)
+    {
+		offset1 = NOK_LT_BLEN - flen - delay + i - (LPC - 1) / 2; 
+		r[i][LPC] = 0.0;
+		for (k = 0; k < flen; k++)
+			r[i][LPC] += x_buffer[offset1 + k] * sb_samples[k];
+    }
+
+	for (i = 0; i < LPC; i++)
+		weight[i] = 0.0;
+
+	weight[(LPC - 1) / 2] = r[(LPC - 1) / 2][LPC] / r[(LPC - 1) / 2][(LPC - 1) / 2];
+
+	lnqgj (r);
+
+	for (i = 0; i < LPC; i++)
+		weight[i] = r[i][LPC];
+
+	/* 
+	fprintf(stdout, "unquantized weight : %f ... ", weight[0]);
+	*/
+
+	w_quantize (weight, ltp_coef);
+
+	/*
+	fprintf(stdout, "quantized weight : %f\n", weight[0]);    
+	*/
+
+	prediction (x_buffer, sb_samples_pred, weight, delay, flen);
+}
+
+
+/**************************************************************************
+  Title:	lnqgj
+
+  Purpose:	Calculates LTP gains.
+
+  Usage:	lnqgj(a)
+
+  Input:	a - auto-correlation matrix
+
+  Output:	a - LTP gains (in the last column)   
+
+  References:	-
+
+  Explanation:	-
+
+  Author(s):	Juha Ojanpera, Lin Yin
+  *************************************************************************/
+
+void
+lnqgj (double (*a)[LPC + 1])
+{
+	int nn, nnpls1, ip, i, j;
+	double p, w;
+
+	for (nn = 0; nn < LPC; nn++)
+    {
+		p = 0.0;
+		nnpls1 = nn + 1;
+		for (i = nn; i < LPC; i++)
+			if (p - fabs (a[i][nn]) < 0)
+			{
+				p = fabs (a[i][nn]);
+				ip = i;
+			}
+
+		if (p - 1.e-10 <= 0)
+			return;
+
+		if (p - 1.e-10 > 0)
+		{
+			for (j = nn; j <= LPC; j++)
+			{
+				w = a[nn][j];
+				a[nn][j] = a[ip][j];
+				a[ip][j] = w;
+			}
+
+			for (j = nnpls1; j <= LPC; j++)
+				a[nn][j] = a[nn][j] / a[nn][nn];
+
+			for (i = 0; i < LPC; i++)
+				if ((i - nn) != 0)
+					for (j = nnpls1; j <= LPC; j++)
+						a[i][j] = a[i][j] - a[i][nn] * a[nn][j];
+		}
+    }
+}
+
+
+/**************************************************************************
+  Title:	w_quantize
+
+  Purpose:	Quantizes LTP gain(s).
+
+  Usage:	w_quantize(freq, ltp_idx)
+
+  Input:	freq    - original LTP gain(s)
+
+  Output:	freq    - LTP gain(s) selected from the codebook
+                ltp_idx - corresponding indices of the LTP gain(s)
+
+  References:	-
+
+  Explanation:	The closest value from the codebook is selected to be the
+                quantized LTP gain.
+
+  Author(s):	Juha Ojanpera, Lin Yin
+  *************************************************************************/
+
+void
+w_quantize (double *freq, int *ltp_idx)
+{
+	int i, j;
+	double dist, low;
+
+
+	low = 1.0e+10;
+	for (i = 0; i < LPC; i++)
+    {
+		dist = 0.0;
+		for (j = 0; j < CODESIZE; j++)
+		{
+			dist = (freq[i] - codebook[j]) * (freq[i] - codebook[j]);
+			if (dist < low)
+			{
+				low = dist;
+				ltp_idx[i] = j;
+			}
+		}
+    }
+
+	for (j = 0; j < LPC; j++)
+		freq[j] = codebook[ltp_idx[j]];
+}
+
--- /dev/null
+++ b/nok_pitch.h
@@ -1,0 +1,43 @@
+/**************************************************************************
+
+This software module was originally developed by
+Nokia in the course of development of the MPEG-2 AAC/MPEG-4 
+Audio standard ISO/IEC13818-7, 14496-1, 2 and 3.
+This software module is an implementation of a part
+of one or more MPEG-2 AAC/MPEG-4 Audio tools as specified by the
+MPEG-2 aac/MPEG-4 Audio standard. ISO/IEC  gives users of the
+MPEG-2aac/MPEG-4 Audio standards free license to this software module
+or modifications thereof for use in hardware or software products
+claiming conformance to the MPEG-2 aac/MPEG-4 Audio  standards. Those
+intending to use this software module in hardware or software products
+are advised that this use may infringe existing patents. The original
+developer of this software module, the subsequent
+editors and their companies, and ISO/IEC have no liability for use of
+this software module or modifications thereof in an
+implementation. Copyright is not released for non MPEG-2 aac/MPEG-4
+Audio conforming products. The original developer retains full right to
+use the code for the developer's own purpose, assign or donate the code to a
+third party and to inhibit third party from using the code for non
+MPEG-2 aac/MPEG-4 Audio conforming products. This copyright notice
+must be included in all copies or derivative works.
+Copyright (c)1997.  
+
+***************************************************************************/
+#ifndef NOK_PITCH_H_
+#define NOK_PITCH_H_
+
+
+extern double snr_pred (double *mdct_in, double *mdct_pred, int *sfb_flag,
+                        int *sfb_offset, enum WINDOW_TYPE block_type, int side_info,
+                        int num_of_sfb);
+
+extern void prediction (short *buffer, double *predicted_samples, double *weight,
+                        int delay, int flen);
+
+extern int estimate_delay (double *sb_samples, short *x_buffer, int flen);
+
+extern void pitch (double *sb_samples, double *sb_samples_pred, short *x_buffer,
+                   int *ltp_coef, int delay, int flen);
+
+
+#endif
--- a/psych.h
+++ b/psych.h
@@ -64,14 +64,10 @@
 /* added by T. Araki (1997.07.10) */
 typedef struct {
   double hw[BLOCK_LEN_LONG*2];     /* Hann window table */
-  double st[BLOCK_LEN_LONG*5/2];   /* sin table*/
-  int    brt[BLOCK_LEN_LONG*2];    /* bit inverse table */
 } FFT_TABLE_LONG;
 
 typedef struct {
   double hw[BLOCK_LEN_SHORT*2];     /* Hann window table */
-  double st[BLOCK_LEN_SHORT*5/2];   /* sin table*/
-  int    brt[BLOCK_LEN_SHORT*2];    /* bit inverse table */
 } FFT_TABLE_SHORT;
  
 typedef struct {
--- a/tns.c
+++ b/tns.c
@@ -138,7 +138,7 @@
 #if 1
 	if (use_tns)
 	/* Doesn't work well on short windows. */
-	if (blockType != ONLY_SHORT_WINDOW)
+	if (blockType == ONLY_LONG_WINDOW)
 	/* Perform analysis and filtering for each window */
 	for (w=0;w<numberOfWindows;w++) {
 
--- a/transfo.c
+++ b/transfo.c
@@ -1,7 +1,10 @@
 #include <math.h>
 #include "transfo.h"
 
-fftw_complex_d FFTarray[2048];    /* the array for in-place FFT */
+fftw_complex_d FFTarray[2048];    /* the array for in-place FFT */
+extern int unscambled64[64];    /* the permutation array for FFT64*/
+extern int unscambled512[512];  /* the permutation array for FFT512*/
+
 #ifndef M_PI
 #define M_PI        3.14159265358979323846
 #endif
@@ -105,5 +108,362 @@
 		c = c * cfreq - s * sfreq;
 		s = s * cfreq + cold * sfreq;
     }
+}
+
+
+/*****************************
+  Fast MDCT Code
+*****************************/
+
+#define SWAP(a,b) tempr=a;a=b;b=tempr
+void CompFFT (double *data, int nn, int isign);
+void FFT (double *data, int nn, int isign);
+
+double *NewFloat (int N) {
+/* Allocate array of N Floats */
+
+    double *temp;
+
+    temp = (double *) malloc (N * sizeof (double));
+    if (!temp) {
+		exit (1);
+		}
+    return temp;
+}
+
+typedef double *pFloat;
+
+double **NewFloatMatrix (int N, int M) {
+/* Allocate NxM matrix of Floats */
+
+    double **temp;
+    int i;
+
+/* allocate N pointers to double arrays */
+    temp = (pFloat *) malloc (N * sizeof (pFloat));
+    if (!temp) {
+		exit (1);
+		}
+
+/* Allocate a double array M long for each of the N array pointers */
+
+    for (i = 0; i < N; i++) {
+		temp [i] = (double *) malloc (M * sizeof (double));
+		if (! temp [i]) {
+			exit (1);
+			}
+		}
+    return temp;
+}
+
+
+void IMDCT (double *data, int N) {
+
+    static double *FFTarray = 0;    /* the array for in-place FFT */
+    static double tempr, tempi, c, s, cold, cfreq, sfreq; /* temps for pre and post twiddle */
+    double freq = 2. * M_PI / N;
+    static double fac;
+    static int i, n;
+	int isign = -1;
+	int b = N/2;
+    int a = N - b;
+
+    /* Choosing to allocate 2/N factor to Inverse Xform! */
+    if (isign == 1) {
+      fac = 2.; /* 2 from MDCT inverse  to forward */
+    } 
+    else {
+      fac = 2. / N; /* remaining 2/N from 4/N IFFT factor */
+    }
+    
+    
+    {
+      static int oldN = 0;
+      if (N > oldN) {
+	if (FFTarray) {
+	  free (FFTarray);
+	  FFTarray = 0;
+	}
+	oldN = N;
+      }
+      
+      if (! FFTarray) 
+	FFTarray = NewFloat (N / 2);  /* holds N/4 complex values */
+    }
+    
+    /* prepare for recurrence relation in pre-twiddle */
+    cfreq = cos (freq);
+    sfreq = sin (freq);
+    c = cos (freq * 0.125);
+    s = sin (freq * 0.125);
+
+    for (i = 0; i < N / 4; i++) {
+      
+      /* calculate real and imaginary parts of g(n) or G(p) */
+      
+      if (isign == 1) { /* Forward Transform */
+	n = N / 2 - 1 - 2 * i;
+	if (i < b / 4) {
+	  tempr = data [a / 2 + n] + data [N + a / 2 - 1 - n]; /* use second form of e(n) for n = N / 2 - 1 - 2i */
+	} 
+	else {
+	  tempr = data [a / 2 + n] - data [a / 2 - 1 - n]; /* use first form of e(n) for n = N / 2 - 1 - 2i */
+	}
+	n = 2 * i;
+	if (i < a / 4) {
+	  tempi = data [a / 2 + n] - data [a / 2 - 1 - n]; /* use first form of e(n) for n=2i */
+	} 
+	else {
+	  tempi = data [a / 2 + n] + data [N + a / 2 - 1 - n]; /* use second form of e(n) for n=2i*/
+	}
+      } 
+      else { /* Inverse Transform */
+	tempr = -data [2 * i];
+	tempi = data [N / 2 - 1 - 2 * i];
+      }
+      
+      /* calculate pre-twiddled FFT input */
+      FFTarray [2 * i] = tempr * c + isign * tempi * s;
+      FFTarray [2 * i + 1] = tempi * c - isign * tempr * s;
+      
+      /* use recurrence to prepare cosine and sine for next value of i */
+      cold = c;
+      c = c * cfreq - s * sfreq;
+      s = s * cfreq + cold * sfreq;
+    }
+    
+    /* Perform in-place complex FFT (or IFFT) of length N/4 */
+    /* Note: FFT has physics (opposite) sign convention and doesn't do 1/N factor */
+    CompFFT (FFTarray, N / 4, -isign);
+    
+    /* prepare for recurrence relations in post-twiddle */
+    c = cos (freq * 0.125);
+    s = sin (freq * 0.125);
+    
+    /* post-twiddle FFT output and then get output data */
+    for (i = 0; i < N / 4; i++) {
+      
+      /* get post-twiddled FFT output  */
+      /* Note: fac allocates 4/N factor from IFFT to forward and inverse */
+      tempr = fac * (FFTarray [2 * i] * c + isign * FFTarray [2 * i + 1] * s);
+      tempi = fac * (FFTarray [2 * i + 1] * c - isign * FFTarray [2 * i] * s);
+      
+      /* fill in output values */
+      if (isign == 1) { /* Forward Transform */
+	data [2 * i] = -tempr;   /* first half even */
+	data [N / 2 - 1 - 2 * i] = tempi;  /* first half odd */
+	data [N / 2 + 2 * i] = -tempi;  /* second half even */
+	data [N - 1 - 2 * i] = tempr;  /* second half odd */
+      } 
+      else { /* Inverse Transform */
+	data [N / 2 + a / 2 - 1 - 2 * i] = tempr;
+	if (i < b / 4) {
+	  data [N / 2 + a / 2 + 2 * i] = tempr;
+	} 
+	else {
+	  data [2 * i - b / 2] = -tempr;
+	}
+	data [a / 2 + 2 * i] = tempi;
+	if (i < a / 4) {
+	  data [a / 2 - 1 - 2 * i] = -tempi;
+	} 
+	else {
+	  data [a / 2 + N - 1 - 2*i] = tempi;
+	}
+      }
+      
+      /* use recurrence to prepare cosine and sine for next value of i */
+      cold = c;
+      c = c * cfreq - s * sfreq;
+      s = s * cfreq + cold * sfreq;
+    }
+    
+    /*    DeleteFloat (FFTarray);   */
+    
+}
+
+
+void CompFFT (double *data, int nn, int isign) {
+
+    static int i, j, k, kk;
+    static int p1, q1;
+    static int m, n, logq1;
+    static double **intermed = 0;
+    static double ar, ai;
+    static double d2pn;
+    static double ca, sa, curcos, cursin, oldcos, oldsin;
+    static double ca1, sa1, curcos1, cursin1, oldcos1, oldsin1;
+
+
+/* Factorize n;  n = p1*q1 where q1 is a power of 2.
+    For n = 1152, p1 = 9, q1 = 128.		*/
+
+    n = nn;
+    logq1 = 0;
+
+    for (;;) {
+		m = n >> 1;  /* shift right by one*/
+		if ((m << 1) == n) {
+		    logq1++;
+		    n = m;
+			} 
+		else {
+			break;
+			}
+		}
+
+    p1 = n;
+    q1 = 1;
+    q1 <<= logq1;
+
+    d2pn = 2*M_PI / nn;
+
+{static int oldp1 = 0, oldq1 = 0;
+
+	if ((oldp1 < p1) || (oldq1 < q1)) {
+		if (intermed) {
+			free (intermed);
+			intermed = 0;
+			}
+		if (oldp1 < p1) oldp1 = p1;
+		if (oldq1 < q1) oldq1 = q1;;
+		}
+		
+	if (!intermed)
+		intermed = NewFloatMatrix (oldp1, 2 * oldq1);
+}
+
+/* Sort the p1 sequences */
+
+    for	(i = 0; i < p1; i++) {
+		for (j = 0; j < q1; j++){
+		    intermed [i] [2 * j] = data [2 * (p1 * j + i)];
+			intermed [i] [2 * j + 1] = data [2 * (p1 * j + i) + 1];
+			}
+		}
+
+/* compute the power of two fft of the p1 sequences of length q1 */
+
+    for (i = 0; i < p1; i++) {
+/* Forward FFT in place for n complex items */
+		FFT (intermed [i], q1, isign);
+		}
+
+/* combine the FFT results into one seqquence of length N */
+
+    ca1 = cos (d2pn);
+    sa1 = sin (d2pn);
+    curcos1 = 1.;
+    cursin1 = 0.;
+
+    for (k = 0; k < nn; k++) {
+		data [2 * k] = 0.;
+		data [2 * k + 1] = 0.;
+		kk = k % q1;
+
+		ca = curcos1;
+		sa = cursin1;
+		curcos = 1.;
+		cursin = 0.;
+
+		for (j = 0; j < p1; j++) {
+			ar = curcos;
+			ai = isign * cursin;
+			data [2 * k] += intermed [j] [2 * kk] * ar - intermed [j] [2 * kk + 1] * ai;
+			data [2 * k + 1] += intermed [j] [2 * kk] * ai + intermed [j] [2 * kk + 1] * ar;
+
+		    oldcos = curcos;
+		    oldsin = cursin;
+			curcos = oldcos * ca - oldsin * sa;
+			cursin = oldcos * sa + oldsin * ca;
+			}
+		oldcos1 = curcos1;
+		oldsin1 = cursin1;
+		curcos1 = oldcos1 * ca1 - oldsin1 * sa1;
+		cursin1 = oldcos1 * sa1 + oldsin1 * ca1;
+		}
+
+}
+
+
+void FFT (double *data, int nn, int isign)  {
+/* Varient of Numerical Recipes code from off the internet.  It takes nn
+interleaved complex input data samples in the array data and returns nn interleaved
+complex data samples in place where the output is the FFT of input if isign==1 and it
+is nn times the IFFT of the input if isign==-1. 
+(Note: it doesn't renormalize by 1/N when doing the inverse transform!!!)
+(Note: this follows physicists convention of +i, not EE of -j in forward
+transform!!!!)
+ */
+
+/* Press, Flannery, Teukolsky, Vettering "Numerical
+ * Recipes in C" tuned up ; Code works only when nn is
+ * a power of 2  */
+
+    static int n, mmax, m, j, i;
+    static double wtemp, wr, wpr, wpi, wi, theta, wpin;
+    static double tempr, tempi, datar, datai;
+    static double data1r, data1i;
+
+    n = nn * 2;
+
+/* bit reversal */
+
+    j = 0;
+    for (i = 0; i < n; i += 2) {
+		if (j > i) {  /* could use j>i+1 to help compiler analysis */
+			SWAP (data [j], data [i]);
+			SWAP (data [j + 1], data [i + 1]);
+			}
+		m = nn;
+		while (m >= 2 && j >= m) {
+			j -= m;
+			m >>= 1;
+			}
+		j += m;
+		}
+
+    theta = 3.141592653589795 * .5;
+    if (isign < 0)
+		theta = -theta;
+    wpin = 0;   /* sin(+-PI) */
+    for (mmax = 2; n > mmax; mmax *= 2) {
+		wpi = wpin;
+		wpin = sin (theta);
+		wpr = 1 - wpin * wpin - wpin * wpin; /* cos(theta*2) */
+		theta *= .5;
+		wr = 1;
+		wi = 0;
+		for (m = 0; m < mmax; m += 2) {
+			j = m + mmax;
+			tempr = (double) wr * (data1r = data [j]);
+			tempi = (double) wi * (data1i = data [j + 1]);
+			for (i = m; i < n - mmax * 2; i += mmax * 2) {
+/* mixed precision not significantly more
+ * accurate here; if removing double casts,
+ * tempr and tempi should be double */
+				tempr -= tempi;
+				tempi = (double) wr *data1i + (double) wi *data1r;
+/* don't expect compiler to analyze j > i + 1 */
+				data1r = data [j + mmax * 2];
+				data1i = data [j + mmax * 2 + 1];
+				data [i] = (datar = data [i]) + tempr;
+				data [i + 1] = (datai = data [i + 1]) + tempi;
+				data [j] = datar - tempr;
+				data [j + 1] = datai - tempi;
+				tempr = (double) wr *data1r;
+				tempi = (double) wi *data1i;
+				j += mmax * 2;
+				}
+			tempr -= tempi;
+			tempi = (double) wr * data1i + (double) wi *data1r;
+			data [i] = (datar = data [i]) + tempr;
+			data [i + 1] = (datai = data [i + 1]) + tempi;
+			data [j] = datar - tempr;
+			data [j + 1] = datai - tempi;
+			wr = (wtemp = wr) * wpr - wi * wpi;
+			wi = wtemp * wpi + wi * wpr;
+			}
+		}
 }
 
--- a/transfo.h
+++ b/transfo.h
@@ -2,6 +2,7 @@
 #define TRANSFORM_H 
 
 void MDCT(double* data, int N);
+void IMDCT(double* data, int N);
 void FFT(double *data, int nn, int isign);
 
 #define c_re(c)  ((c).re)