/* SigLib Analytical Gignal By Applying A Hilbert Transform
Hilbert transform applied to a real input signal on the imaginary
output path and applies a corresponding delay to the real
output. The flow diagram is as follows :

	  +--     Delay     ----> R
	--|
	  +-- Hilbert xform ----> I

Analytic signals give positive output frequencies but no
negative output frequencies.

The Hilbert Transform uses an N tap FIR filter to phase shift every
component in a signal by 90 degrees (N is odd ordered).

The defining equations for the Hilbert transform are :

          2        2 ( n * PI )
h(n) = ______ * sin  (________)		for n = +-1, +-2, . . +- N/2
       n * PI        (	  2   )

and h(0) = 0, for n = 0
*/

#include <stdio.h>
#include <siglib.h>
#include "GraphFunctions.h"

/* Define constants */
#define	SAMPLE_LENGTH	((SLArrayIndex_t)512)
#define	FFT_SIZE		SAMPLE_LENGTH
#define	LOG_FFT_SIZE	((SLArrayIndex_t)9)

#define	FILTER_LENGTH	31L
#define	DELAY_LENGTH	((FILTER_LENGTH + 1L) / 2L)

SLData_t	SinePhase;
SLData_t	*pImagData, *pMagnitude, *pDelay, *pFFTCoeffs, *pTempDelay;
SLData_t	*pFilterTaps, *pFilterTapsStart, *pFilterState;
SLArrayIndex_t	FilterIndex;
SLData_t	*pSrc1, *pSrc2;

void	main(void);


void main(void)
{
	GraphObject *h2DGraph;							/* Declare graph object */

	h2DGraph =										/* Initialize graph */
		Create2DGraph ("Analytic Signal",			/* Graph title */
					   "Time",						/* X-Axis label */
					   "Magnitude",					/* Y-Axis label */
					   SV_AUTO_SCALE,				/* Scaling mode */
					   SV_SIGNED,					/* Sign mode */
					   SV_GRAPH_LINE,				/* Graph type */
					   "localhost");				/* Graph server */
	if (h2DGraph == NULL)							/* Graph creation failed - e.g is server running ? */
	{
		printf ("\nGraph creation failure. Please check that the server is running\n");
		exit (1);
	}

	printf ("Hilbert transform filter length   => %d", FILTER_LENGTH);

	pFilterTaps = SUF_VectorArrayAllocate (FILTER_LENGTH);
	pFilterState = SUF_VectorArrayAllocate (FILTER_LENGTH);
	pSrc1 = SUF_VectorArrayAllocate (SAMPLE_LENGTH);
	pSrc2 = SUF_VectorArrayAllocate (SAMPLE_LENGTH);
	pImagData = SUF_VectorArrayAllocate (SAMPLE_LENGTH);
	pMagnitude = SUF_VectorArrayAllocate (SAMPLE_LENGTH);
	pFFTCoeffs = SUF_FftCoefficientAllocate (FFT_SIZE);
	pTempDelay = SUF_VectorArrayAllocate (DELAY_LENGTH);
	pDelay = SUF_VectorArrayAllocate (DELAY_LENGTH);

	pFilterTapsStart = pFilterTaps;
													/* Initialise Hilbert transformer coefficients */
	SIF_HilbertTransformer (pFilterTaps,			/* Pointer to filter coefficients */
							FILTER_LENGTH);			/* Filter length */
													/* Initialise FIR filter for Hilbert transformer */
	SIF_Fir (pFilterState,							/* Pointer to filter state array */
			 &FilterIndex,							/* Pointer to filter index register */
			 FILTER_LENGTH);						/* Filter length */
													/* Initialise the delay function */
	SIF_Delay (pDelay,								/* Pointer to filter index register */
			   SIGLIB_NULL_ARRAY_INDEX_PTR,			/* Unused */
			   DELAY_LENGTH);						/* Filter length */

													/* Initialise FFT */
	SIF_Fft (pFFTCoeffs,							/* Pointer to FFT coefficients */
			 SIGLIB_NULL_ARRAY_INDEX_PTR,			/* Pointer to bit reverse address table - NOT USED */
			 FFT_SIZE);								/* FFT Size */

	SinePhase = SIGLIB_ZERO;
	SDA_SignalGenerate (pSrc1,						/* Pointer to destination array */
						SIGLIB_SINE_WAVE,			/* Signal type - Sine wave */
						SIGLIB_ONE,					/* Signal peak level */
						SIGLIB_FILL,				/* Fill (overwrite) or add to existing array contents */
						((SLData_t)0.05),			/* Signal frequency */
						SIGLIB_ZERO,				/* D.C. Offset */
						SIGLIB_ZERO,				/* Unused */
						SIGLIB_ZERO,				/* Signal end value - Unused */
						&SinePhase,					/* Signal phase - maintained across array boundaries */
						SIGLIB_NULL_DATA_PTR,		/* Unused */
						SAMPLE_LENGTH);				/* Output array length */

	SinePhase = SIGLIB_ZERO;
	SDA_SignalGenerate (pSrc2,						/* Pointer to destination array */
						SIGLIB_SINE_WAVE,			/* Signal type - Sine wave */
						SIGLIB_ONE,					/* Signal peak level */
						SIGLIB_FILL,				/* Fill (overwrite) or add to existing array contents */
						((SLData_t)0.005),			/* Signal frequency */
						SIGLIB_ZERO,				/* D.C. Offset */
						SIGLIB_ZERO,				/* Unused */
						SIGLIB_ZERO,				/* Signal end value - Unused */
						&SinePhase,					/* Signal phase - maintained across array boundaries */
						SIGLIB_NULL_DATA_PTR,		/* Unused */
						SAMPLE_LENGTH);				/* Output array length */

	SDA_Multiply2 (pSrc1, pSrc2, pSrc1, SAMPLE_LENGTH);

	Display2DGraph (h2DGraph,						/* Graph handle */
				    "Source signal",				/* Title of the dataset */
				    pSrc1,							/* Array of Double dataset */
				    SAMPLE_LENGTH,					/* Number of data points */
					SV_GRAPH_LINE,					/* Graph type */
				    SV_MAGENTA,						/* Colour */
					SV_HIDE_MARKERS,				/* Marker enable / disable */
					SV_GRAPH_NEW);					/* New graph */
	printf ("\nSource signal\nPlease hit <Carriage Return> to continue . . ."); getchar ();

	SDA_Copy (pSrc1,								/* Pointer to source array */
			  pSrc2,								/* Pointer to destination array */
			  SAMPLE_LENGTH);						/* Array length */

													/* Perform real FFT */
	SDA_Rfft (pSrc2,								/* Pointer to real array */
			  pImagData,							/* Pointer to imaginary array */
			  pFFTCoeffs,							/* Pointer to FFT coefficients */
			  SIGLIB_NULL_ARRAY_INDEX_PTR,			/* Pointer to bit reverse address table - NOT USED */
			  FFT_SIZE,								/* FFT size */
			  LOG_FFT_SIZE);						/* log2 FFT size */

													/* Calculate real magnitude from complex */
	SDA_Magnitude (pSrc2,							/* Pointer to real source array */
				   pImagData,						/* Pointer to imaginary source array */
				   pMagnitude,						/* Pointer to magnitude destination array */
				   SAMPLE_LENGTH);					/* Array length */

	SDA_10Log10 (pMagnitude,						/* Pointer to source array */
				 pMagnitude,						/* Pointer to destination array */
				 SAMPLE_LENGTH);					/* Array length */
	Display2DGraph (h2DGraph,						/* Graph handle */
				    "Spectrum of source signal",				/* Title of the dataset */
				    pMagnitude,						/* Array of Double dataset */
				    SAMPLE_LENGTH,					/* Number of data points */
					SV_GRAPH_LINE,					/* Graph type */
				    SV_MAGENTA,						/* Colour */
					SV_HIDE_MARKERS,				/* Marker enable / disable */
					SV_GRAPH_NEW);					/* New graph */
	printf ("\nSpectrum of source signal\nPlease hit <Carriage Return> to continue . . ."); getchar ();

													/* Apply Hilbert transformerer */
	SDA_Fir (pSrc1,									/* Input array to be filtered */
			 pSrc2,									/* Filtered output array */
			 pFilterState,							/* Pointer to filter state array */
			 pFilterTaps,							/* Pointer to filter coefficients */
			 &FilterIndex,							/* Pointer to filter index register */
			 FILTER_LENGTH,							/* Filter length */
			 SAMPLE_LENGTH);						/* Array length */

													/* Delay real component to generate an analytical signal */
	SDA_ShortDelay (pSrc1,							/* Pointer to source array */
					pSrc1,							/* Pointer to destination array */
					pDelay,							/* Pointer to temporary delayed array */
					pTempDelay,						/* Temporary destination array pointer */
					DELAY_LENGTH,					/* Sample delay count */
					SAMPLE_LENGTH);					/* Array length */

	Display2DGraph (h2DGraph,						/* Graph handle */
				    "Real output",					/* Title of the dataset */
				    pSrc1,							/* Array of Double dataset */
				    SAMPLE_LENGTH,					/* Number of data points */
					SV_GRAPH_LINE,					/* Graph type */
				    SV_MAGENTA,						/* Colour */
					SV_HIDE_MARKERS,				/* Marker enable / disable */
					SV_GRAPH_NEW);					/* New graph */
	Display2DGraph (h2DGraph,						/* Graph handle */
				    "Imaginary output",				/* Title of the dataset */
				    pSrc2,							/* Array of Double dataset */
				    SAMPLE_LENGTH,					/* Number of data points */
					SV_GRAPH_LINE,					/* Graph type */
				    SV_BLUE,						/* Colour */
					SV_HIDE_MARKERS,				/* Marker enable / disable */
					SV_GRAPH_ADD);					/* New graph */

	SUF_MemoryFree (pFilterTaps);					/* Free memory */
	SUF_MemoryFree (pFilterState);
	SUF_MemoryFree (pSrc1);
	SUF_MemoryFree (pSrc2);
	SUF_MemoryFree (pImagData);
	SUF_MemoryFree (pMagnitude);
	SUF_MemoryFree (pDelay);
	SUF_MemoryFree (pFFTCoeffs);
	SUF_MemoryFree (pTempDelay);

}


