#include <malloc.h>
#include <math.h>
#include <stdio.h>
#include <stdlib.h>
#include "mex.h"

#define	SEQUENCE_LENGTH	33554432	/* 2^25 */
#define CODE_LENGTH		38400
#define BUFFER_LENGTH	25
#define ARRAY_LENGTH	26
#define LOG2(arg) (log( (double) arg)/log(2))

#define OFFSET		16777232
#define COUNT1		16777199
#define MAX_CODE_NUMBER	16777215
#define NDIMS			2



void scramble_long(unsigned long code_num,short int *scramble_real,short int *scramble_imag)
/*****************************************************************************************************
/void scramble_long(unsigned long code_num,short int *scramble_real,short int *scramble_imag)
/
/
/ Copyright 2002 The Mobile and Portable Radio Research Group
/
/This function generages the "long scramble code" by using the algorithm specified in 
/ETSI TS 125 213 V3.2.0 (2000-03).  The function generates this code by using the implementation 
/that is pictorially described in figure 5 of the aforementioned specificatoin.  
/This implemtation takes advantage of two shift registers, which this function emulates via 
/circular buffers.  The first circular buffer, X, is initiallized using code_num. The second circular 
/buffer Y is initiallized by setting all elements to 1.
/
/The function then performs the following operation using modulo 2 arithmetic
/	x(n+25) = x(n+3) + x(n+1)
/	y(n+25) = y(n+3) + y(n+2) + y(n+1) + y(n)
/	z(n) = x(n+4) + x(n+7) + x(n+18)
/	w(n) = y(n+4) + y(n+6) + y(n+17)
/	Clong1 = x(n) + y(n)
/	Clong2 = w(n) + z(n)
/	Scramble_long(n) = Clong1 * (1 + j*((-1)^n)*Clong2(2*floor(n/2)))
/
/Parameters
/	Input
/		code_num		unsigned value between 0 and 16777215 that determines the scramble code number
/
/	Output
/		*scramble_real	double valued array that stores the real part of the scramble code
/		*scramble_imag	double valued array that stores the imaginary part of the scramble code  
/******************************************************************************************************/

{
	unsigned long temp_code_num;
	unsigned short *x_buf_front,*x_buf_end;		//Points to the begining and end of the X circular buffer
	unsigned short *y_buf_front,*y_buf_end;		//Points to the begining and end of the Y circular buffer
	unsigned short *x0;		//Points to x(n)
	unsigned short *x3;		//Points to x(n+3)
	unsigned short *x4;		//Points to x(n+4)
	unsigned short *x7;		//Points to x(n+7)
	unsigned short *x18;	//Points to x(n+18)
	unsigned short *x25;	//Points to x(n+25)
	unsigned short *x_temp;	//Temporary Pointer
	unsigned short *y0;		//Points to y(n)
	unsigned short *y1;		//Points to y(n+1)
	unsigned short *y2;		//Points to y(n+2)
	unsigned short *y3;		//Points to y(n+3)
	unsigned short *y4;		//Points to y(n+4)
	unsigned short *y6;		//Points to y(n+6)
	unsigned short *y17;	//Points to y(n+17)
	unsigned short *y25;	//Points to y(n+25)
	unsigned short *y_temp;	//Temporary Pointer	
	unsigned short x1_hold; //Stores x(n)
	unsigned short y1_hold; //Stores y(n)
	unsigned short x2_hold;	//Stores the result of x(n+4) + x(n+7) + x(n+18)
	unsigned short y2_hold;	//Stores the result of y(n+4) + y(n+6) + y(n+17)
	short int *c_long_1,*c_long_1_temp;	//Pointers and temporary pointers to Clong1
	short int *c_long_2,*c_long_2_temp;	//Pointers and temporary pointers to Clong2
	short int *temp_real,*temp_imag;	//Temporary Pointers for *scramble_real and *scramble_imag
	int k,num_bits;
//	FILE *fp;


	/* Will shift registers using two circular buffers: 
	   one for the "x" sequence and one for the "y" sequence
	*/

	//Allocate space for code sequences
	if ((c_long_1 = (short int *) mxCalloc(CODE_LENGTH,sizeof(short int)))==NULL)
	{
		printf("\nc_long_1 array not allocated!--exiting\n");
		exit(-2);
	}
	if ( (c_long_2 = (short int *) mxCalloc(CODE_LENGTH,sizeof(short int)))==NULL)
	{
		printf("\nc_long_2 array not allocated!--exiting\n");
		exit(-2);
	}
	c_long_1_temp = c_long_1;
	c_long_2_temp = c_long_2;

	//Allocate buffers
	if ((x_buf_front = (unsigned short *) mxCalloc(ARRAY_LENGTH,sizeof(unsigned short)))==NULL)
	{
		printf("\nx_buf_front array not allocated!--exiting\n");
		exit(-2);
	}
	if ((y_buf_front = (unsigned short *) mxCalloc(ARRAY_LENGTH,sizeof(unsigned short)))==NULL)
	{
		printf("\n y_buf_front array not allocated!--exiting\n");
		exit(-2);
	}

	//Assign pointers to the end of buffer
	x_buf_end=x_buf_front+BUFFER_LENGTH;
	y_buf_end=y_buf_front+BUFFER_LENGTH;

	//Initialize buffers
	/* the code number (code_num) determines the inital condition for
	   the x seqeunce generator.  We therefore need to determine the 
	   binary representation of code_num and store it in the first 24
	   elements of x */
	
	x_temp=x_buf_front;
	temp_code_num=code_num;
	num_bits = (int) (LOG2(temp_code_num)+1);
	for (k=0; k<num_bits; k++)
	{
		*x_temp++ = 0x00000001 & temp_code_num;
		temp_code_num >>= 1;		
	}
	//Set the 24th element in the array equal to 1
	*(x_buf_front+24) = 1;

	//Initialize the "y" buffer by setting all 25 elements equal to 1
	y_temp=y_buf_front;
	for (k=0;k<BUFFER_LENGTH;k++) *y_temp++ = 1;

	//We now implement the sequence computation.
	//First initialize the locations of all pointers
	//Pointers for "x" buffer
	x0=x_buf_front;			//x(n)
	x3=x_buf_front+3;		//x(n+3)
	x4=x_buf_front+4;		//x(n+4)
	x7=x_buf_front+7;		//x(n+7)
	x18=x_buf_front+18;		//x(n+18)
	x25=x_buf_front;		//x(n+25)
	//Pointers for "y" buffer
	y0=y_buf_front;			//y(n)
	y1=y_buf_front+1;		//y(n+1)
	y2=y_buf_front+2;		//y(n+2)
	y3=y_buf_front+3;		//y(n+3)
	y4=y_buf_front+4;		//y(n+4)
	y6=y_buf_front+6;		//y(n+6)
	y17=y_buf_front+17;		//y(n+17)
	y25=y_buf_front;		//y(n+25)

	temp_real=scramble_real;	//assign temporary pointers to output parameter
	temp_imag=scramble_imag;		//assign temporary pointers to output parameter

	for (k=0;k<CODE_LENGTH;k++)
	{
		x1_hold = *x0;	//Shift register output for "x" sequence
		y1_hold = *y0;	//Shift register output for "x" sequence
		
		//Perform Shift register operations
		*x25 = *x0++ ^ *x3++;		//x(n+25) = x(n) + x(n+3)
		*y25 = (((*y0++ ^ *y1++) ^ *y2++) ^ *y3++);	//y(n+25) = y(n) + y(n+1) + y(n+2) + y(n+3)
		x2_hold = ((*x4++ ^ *x7++) ^ *x18++);		//z(n) = x(n+4) + x(n+7) + x(n+18)
		y2_hold = ((*y4++ ^ *y6++) ^ *y17++);		//w(n) = y(n+4) + y(n+6) + y(n+17)
		x25++;
		y25++;
				
		//Create output for c_long_1 and c_long_2 sequences
		// convert a "0" to +1 and a "1" to -1
		*c_long_1_temp++ =  (1 - 2 * (short int) (x1_hold ^ y1_hold));
		*c_long_2_temp++ =  (1 - 2 * (short int) (x2_hold ^ y2_hold));

		//Check for pointer overflow
		if (x18 == x_buf_end) x18 = x_buf_front;
		else if (x7 == x_buf_end) x7 = x_buf_front;
		else if (x4 == x_buf_end) x4 = x_buf_front;
		else if (x3 == x_buf_end) x3 = x_buf_front;
		else if (x0 == x_buf_end) 
		{
			x0 = x_buf_front;
			x25 = x_buf_front;
		}
		if (y17 == y_buf_end) y17 = y_buf_front;
		else if (y6 == y_buf_end) y6 = y_buf_front;
		else if (y4 == y_buf_end) y4 = y_buf_front;
		else if (y3 == y_buf_end) y3 = y_buf_front;
		else if (y2 == y_buf_end) y2 = y_buf_front;
		else if (y1 == y_buf_end) y1 = y_buf_front;
		else if (y0 == y_buf_end) 
		{
			y0 = y_buf_front;
			y25 = y_buf_front;
		}
			
	}

//	Used for debugging purposes
/*
	fp=fopen("sequence_new.txt","w");
	c_long_2_temp=c_long_2;
	c_long_1_temp=c_long_1;
	

	fprintf(fp,"index c1 c2 index c1 c2 index c1 c2\n");
	for (k=0;k<101;k+=2) 
	{
		fprintf( fp,"%d \t %d \t %d \n",k,*c_long_1_temp++,*c_long_2_temp);
		fprintf( fp,"%d \t %d \t %d \n",k,*c_long_1_temp++,*c_long_2_temp);
		c_long_2_temp+=2;
	}
	fclose(fp);
*/	

	//Fill real part of the complex scramble code
	temp_real=scramble_real;
	c_long_1_temp=c_long_1;
	for (k=0;k<CODE_LENGTH;k++) *temp_real++ = *c_long_1_temp++;

	//Fill imaginary part part of the complex scramble code
	c_long_1_temp=c_long_1;
	c_long_2_temp=c_long_2;
	temp_imag=scramble_imag;
	for (k=0;k<CODE_LENGTH;k+=2)
	{
		*temp_imag++ = (*c_long_1_temp++) * (*c_long_2_temp);
		*temp_imag++ = -(*c_long_1_temp++) * (*c_long_2_temp);
		c_long_2_temp+=2;
	}
	mxFree(c_long_1);
	mxFree(c_long_2);
	mxFree(x_buf_front);
	mxFree(y_buf_front);
}