Pacific-Design.com

    
Home Index

1. Caltech

2. Decode LFP cpp

Caltech / Decode LFP cpp /


/* 
   California Institute of Technology
   The Andersen Lab, Division of Biology

   File: Decode_LFP.cpp

   Description: Plexon neural recording system that counts spikes from strobed time 
   value until strobed stop time, then predicts the category and writes it to the 
   serial port.  Note, this package along with the Decode program is the most recent 
   single unit real time system.

    Begin          : June 1, 2001
    Written        : by Kevin Duraj
    Email          : kevin@vis.caltech.edu 
    Collaborator   : Daniella Meeker 
*/

#include "Decode_LFP.h"
int main(int argc, char **argv) {

    HANDLE hServerPollEvent;                    // Handle for plexon
    HANDLE hComm;                               // Handle for serial communications
    LPCTSTR portName = "COM2:";                 // Communications port
    PL_Event     *pServerTimeStampBuffer;
    double** Means;
    double** Stds;
    int     index;                               // Indicates which one is the most probable
    double  x, y;
    double Probability= 0.0;                     // the probability of the most probable direction
    int*  SpikeSum;                              // SpikeSum contains # of spikes in each neuron
    double  ReachProb[NUM_DIRECTIONS] = {0.0};   // Reach probability.
    int ReachDirection;                          // the true reach direction
    int TrialCounters[NUM_DIRECTIONS] = {0};     // number of db trials in each direction
    long StartTime = 0;
    long EndTime = 0;
    bool stopTrial = false, startTrial = false;
    int NumActiveUnits=0;
    int NumLFPChannels=1;
    //the number of local field potential channels being used.  this should maybe be an input, but...
    int NumUnitsPerDSP[MAX_NUM_UNITS];

    /******************* Setup Functions *******************/
    ReadParameters(argc, argv);

    if (!SetupPlexon(&hServerPollEvent))
        exit(PLEXON_ERROR);
    else
        printf("Setup Plexon\n");
    PL_GetNumUnits(NumUnitsPerDSP);
    if (!ReadNumNeurons(NumActiveUnits))
        exit (READ_DATABASE_ERROR);
    if (!CheckNumNeurons(NumActiveUnits, NumUnitsPerDSP))
        exit (WRONG_NUM_NEURONS_ERROR);
    Sleep(3000);
    if (!AllocateMemory(&pServerTimeStampBuffer, &Means, &Stds, &SpikeSum, NumActiveUnits))
        exit(NO_MEMORY_ERROR);
    if (ReadInDatabase(Means, Stds, NumActiveUnits, TrialCounters))
        printf("Using Existing Database...\n");
    else
        printf("Creating New Database...\n");
    if (!SetupSerialPort(hComm, portName))
        exit(SERIAL_PORT_ERROR);
    else
        printf("Setup Serial Port\n");
    /*******************    Main Loop    *******************/
    while (true)
    {
        int NumTimeStamps = MAX_TIMESTAMPS_PER_READ;
        PL_GetTimeStampStructures(&NumTimeStamps, pServerTimeStampBuffer);
        for (int i=0; i < NumTimeStamps; i++)
        {
            // Start Signal
            if (IsStartEvent(pServerTimeStampBuffer[i]))
            {
                CleanupTrial(startTrial, stopTrial, ReachProb, SpikeSum, NumActiveUnits);
                ProcessStartSignal(pServerTimeStampBuffer[i], StartTime, startTrial, stopTrial, ReachDirection);
            }
            // Data
            else if (startTrial && (! stopTrial) && IsDataEvent(pServerTimeStampBuffer[i]))
            {
                ReadDataEvent(pServerTimeStampBuffer[i], SpikeSum, NumUnitsPerDSP);
            }
            // Stop Signal
            else if (IsStopEvent(pServerTimeStampBuffer[i]))
            {
                ProcessStopSignal(pServerTimeStampBuffer[i], EndTime, startTrial, stopTrial);
                PrintSpikeSums(SpikeSum, NumActiveUnits);
                CalculateProbabilities(SpikeSum, LFP, ReachProb, Means, Stds, StartTime, EndTime, NumActiveUnits, TrialCounters);
                SelectDirections(Means,Stds,NumActiveUnits); //shiyan
                index = FindMostProbableDirection(ReachProb, Probability);
                CalculateXY(index, x, y);
                WriteOutDecodeData(hComm, index, Probability, x, y);
            }
            // Database Success
            else if (IsDatabaseSuccessEvent(pServerTimeStampBuffer[i]))
            {
                printf("Getting Database Success = %d\r\n", pServerTimeStampBuffer[i].Unit);
                if (! stopTrial) printf("\t GOT SUCCESS WITHOUT STOP!!!\n");
                CalculateDatabase(SpikeSum, Means, Stds, TrialCounters, StartTime, EndTime, ReachDirection,
                                  NumActiveUnits, NumLFPChannels);
            }
            // Failure
            else if (IsFailureEvent(pServerTimeStampBuffer[i]))
            {
                printf("Getting Failure Code = %d\r\n", pServerTimeStampBuffer[i].Unit);
            }
            // Cancel Signal
            else if (IsCancelEvent(pServerTimeStampBuffer[i]))
            {
                printf("Getting Cancel Code = %d\r\n", pServerTimeStampBuffer[i].Unit);
            }
        }
    }
    CleanupPlexon(hServerPollEvent, pServerTimeStampBuffer);
    DeallocateMemory(&pServerTimeStampBuffer, &Means, &Stds);
    CloseHandle(hComm);
    return 0;
}
/*****************************************************************************
WriteDatabase
*/
void WriteDatabase(double** Means, double** Stds, int NumActiveNeurons, int* TrialCounters)
{
    int i,j;
    ofstream fout(gFileName);
    fout << NumActiveNeurons << "\n";
    for (i = 0; i < NUM_DIRECTIONS; i++)
    {
        fout << TrialCounters[i] << " ";
    }
    fout << "\n";
    for (j = 0; j < NumActiveNeurons; j++)
    {
        for (i = 0; i < NUM_DIRECTIONS; i++)
        {
            fout << Means[j][i] << " " << Stds[j][i] << " \n";
        }
        fout << "\n";
    }
    fout.close();
}
/*****************************************************************************
PrintSpikeSums
*/
void PrintSpikeSums(int* Spikes, int NumActiveNeurons)
{
    printf("\n");
    for (int i = 0 ; i < NumActiveNeurons; i++)
    {
        printf("Neuron %d) Spikes = %d\n", i, Spikes[i]);
    }
    printf("\n");
}
/*****************************************************************************
CalculateDatabase
Calculating means and standard deviations for a new database (where we know
the total number of trials in each direction)
*/
void CalculateDatabase(int* Spikes, int* LFP, double** Means, double** Stds, int* TrialCounters,
                       long StartTime, long EndTime, int ReachDirection, int NumActiveNeurons, int NumLFPChannels)
{
    int i;
    double oldMean, deltaMean;

    double  length = (double)(EndTime-StartTime)/TIME_CONSTANT;
    printf("length=%g", length);
    printf("total time =%g\nStartTime=%g\nEndTime=%g\n", length, StartTime, EndTime);
    TrialCounters[ReachDirection]++;
    printf("Completed %d trials in direction %d\n", TrialCounters[ReachDirection], ReachDirection);


    for (i = 0; i < NumActiveNeurons; i++)
    {

        oldMean = Means[i][ReachDirection];
        Means[i][ReachDirection] = (oldMean*(TrialCounters[ReachDirection]-1)+Spikes[i]/length)/TrialCounters[ReachDirection];
        deltaMean = Means[i][ReachDirection]-oldMean;
        Stds[i][ReachDirection] = sqrt((pow(Stds[i][ReachDirection],2.0)*
                                        (TrialCounters[ReachDirection]-1.0)+pow(deltaMean,2.0)*
                                        TrialCounters[ReachDirection]+pow(Spikes[i]/length-Means[i][ReachDirection],2.0))/TrialCounters[ReachDirection]);
    }
    WriteDatabase(Means, Stds, NumActiveNeurons, TrialCounters);
}
/*****************************************************************************
CleanupTrial
*/
void CleanupTrial(bool& startTrial, bool& stopTrial, double* ReachProb, int* SpikeSum, int NumActiveUnits)
{
    startTrial = false;
    stopTrial = false;
    for (int j = 0; j < NUM_DIRECTIONS; j++)
    {
        ReachProb[j] = 0.0;
    }
    for (int k = 0; k < NumActiveUnits; k++)
    {
        SpikeSum[k] = 0;
    }
}
/*****************************************************************************
CheckNumNeurons
*/
bool CheckNumNeurons(int NumActiveUnits, int* NumUnitsPerDSP)
{
    int TotalNumUnitsPerDSP = 0;
    for (int i = 0; i < MAX_NUM_UNITS ; i++)
        TotalNumUnitsPerDSP+=NumUnitsPerDSP[i];
    if (NumActiveUnits != TotalNumUnitsPerDSP)
    {
        printf("Database and Plexon.GetNumUnits Disagree on Number of Neurons!!!");
        Sleep(5000);
        return false;
    }
    return true;
}
/*****************************************************************************
GetNeuronNum
*/
int GetNeuronNum(int* NumUnitsPerDSP, int channel, int unit)
{
    int neuronNum=0;
    for (int i = 0; i < (channel-1) ; i++)
        neuronNum+=NumUnitsPerDSP[i];
    neuronNum += (unit-1);
    return neuronNum;
}
/*****************************************************************************
ReadDecodeDataEvent
*/
void ReadDataEvent(PL_Event theEvent, int* SpikeSum, int* NumUnitsPerDSP)
{
    if (theEvent.Unit != 0)
        SpikeSum[GetNeuronNum(NumUnitsPerDSP,theEvent.Channel,theEvent.Unit)]++;
}
/*****************************************************************************
ProcessStopSignal
*/
void ProcessStopSignal(PL_Event theEvent, long& EndTime, bool startTrial, bool& stopTrial)
{
    printf("Got stop code!\n");
    EndTime = theEvent.TimeStamp;
    //printf("got EndTime: %g\n", (double)EndTime);
    if (!startTrial)
    {
        printf("\t BUT DIDN'T HAVE START!!!\n");
        Sleep(3000);
        exit(STOP_WITHOUT_START_ERROR);
    }
    else
    {

        stopTrial = true;
    }
}
/*****************************************************************************
ProcessStartSignal
*/
void ProcessStartSignal(PL_Event theEvent, long& StartTime, bool& startTrial, bool stopTrial, int& reachDirection)
{
    printf("Got start code for direction %d!\n", theEvent.Unit);
    StartTime = theEvent.TimeStamp + 200;
    // add 200 ms to shave off transient effect
    printf("NEW DECODER!!!: %g\n", (double)StartTime);
    if (stopTrial)
    {
        printf("\t BUT ALREADY GOT STOP!!!\n");
        Sleep(3000);
        exit(STOP_BEFORE_START_ERROR);
    }

    if (startTrial)
    {
        printf("\t BUT ALREADY HAD ONE!!!\n");
        Sleep(3000);
        exit(ALREADY_HAD_START_ERROR);
    }
    //StartTime = theEvent.TimeStamp;
    startTrial = true;
    reachDirection = theEvent.Unit-1; //plexon indicates start by sending 1-8
}
/*****************************************************************************
Checking Event Types
*/
bool IsCancelEvent(PL_Event theEvent)
{
    return (theEvent.Type == PL_ExtEventType
            && theEvent.Channel == PL_StrobedExtChannel);
}
bool IsDataEvent(PL_Event theEvent)
{
    return (theEvent.Type == PL_SingleWFType);
}
bool IsStopEvent(PL_Event theEvent)
{
    return ( theEvent.Type == PL_ExtEventType
             && theEvent.Channel == PL_StrobedExtChannel
             && theEvent.Unit == STOP_CODE);
}
bool IsStartEvent(PL_Event theEvent)
{
    return (theEvent.Type == PL_ExtEventType
            && theEvent.Channel == PL_StrobedExtChannel
            && theEvent.Unit >= MIN_START_CODE
            && theEvent.Unit <= MAX_START_CODE);
}
bool IsDatabaseSuccessEvent(PL_Event theEvent)
{
    return ( theEvent.Type == PL_ExtEventType
             && theEvent.Channel == PL_StrobedExtChannel
             && theEvent.Unit == DATA_SUCCESS_CODE);
}
bool IsDecodeSuccessEvent(PL_Event theEvent)
{
    return ( theEvent.Type == PL_ExtEventType
             && theEvent.Channel == PL_StrobedExtChannel
             && theEvent.Unit == DECODE_SUCCESS_CODE);
}
bool IsFailureEvent(PL_Event theEvent)
{
    return (theEvent.Type == PL_ExtEventType
            && theEvent.Channel == PL_StrobedExtChannel
            && theEvent.Unit >= MIN_FAILURE_CODE
            && theEvent.Unit <= MAX_FAILURE_CODE);

}
/*****************************************************************************
AllocateMemory
*/
bool AllocateMemory(PL_Event **pServerTimeStampBuffer, double** Means[], double** Stds[],
                    int* Spikes[], int NumActiveUnits)
{
    *pServerTimeStampBuffer = (PL_Event*) malloc(MAX_TIMESTAMPS_PER_READ*sizeof(PL_Event));
    if (pServerTimeStampBuffer == NULL)
    {
        printf("Couldn't allocate memory, I can't continue!\r\n");
        Sleep(3000); // pause before console window closes
        return false;
    }
    (*Spikes) = (int*) malloc(NumActiveUnits*sizeof(int));
    (*Means) = (double**) malloc(NumActiveUnits*sizeof(double*));
    (*Stds) = (double**) malloc(NumActiveUnits*sizeof(double*));
    for (int i = 0; i < NumActiveUnits; i++)
    {
        (*Means)[i]= (double*) malloc(NUM_DIRECTIONS*sizeof(double));
        (*Stds)[i]= (double*) malloc(NUM_DIRECTIONS*sizeof(double));
    }
    return true;
}
/*****************************************************************************
DeallocateMemory
*/
void DeallocateMemory(PL_Event** pServerTimeStampBuffer, double** Means[], double** Stds[])
{
    free(*pServerTimeStampBuffer);
    for (int i = 0; i < NUM_CHANNELS*NUM_NEURONS; i++)
    {
        free((*Means)[i]);
        free((*Stds)[i]);
    }
    free(*Means);
    free(*Stds);
}
/*****************************************************************************
SetupSerialPort
*/
bool SetupSerialPort(HANDLE& hComm, LPCTSTR portName)
{

    hComm = CreateFile( portName, GENERIC_WRITE, 0, 0, OPEN_EXISTING, 0, 0);
    if (hComm == INVALID_HANDLE_VALUE)
    {
        fprintf(stderr, "Problem opening serial port %s ...\r\n", portName);
        return false;
    }
    DCB dcb = {0};
    if (!GetCommState(hComm, &dcb))
    {
        printf("Couldn't get DCB!!!");
        Sleep(5000);
        return false;
    }
    dcb.BaudRate = LAB_VIEW_BAUD_RATE;
    dcb.StopBits = LAB_VIEW_STOP_BITS;
    dcb.fDtrControl = LAB_VIEW_DTR_CONTROL;
    dcb.fRtsControl = LAB_VIEW_RTS_CONTROL;
    dcb.Parity = LAB_VIEW_PARITY;
    dcb.ByteSize = LAB_VIEW_BYTE_SIZE;
    if (! SetCommState(hComm, &dcb))
    {
        fprintf(stderr, "Problem Setting Comm State...\n");
        Sleep(5000);
        return false;
    }
    return true;
}
/*****************************************************************************
WriteSerialData
*/
bool WriteSerialData(HANDLE hComm, char * stringBuffer)
{
    DWORD numBytesToWrite = strlen(stringBuffer);
    DWORD numBytesWritten;
    if (!WriteFile(hComm, stringBuffer, numBytesToWrite, &numBytesWritten, 0))
    {
        fprintf(stderr, "WriteFile failed: Error Code = %d, # Bytes Written=%d \r\n",
                GetLastError(), numBytesWritten);
        return false;
    }
    else
    {
        printf("Write file suceeded: %d bytes to write, %d bytes written : %c\n",
               numBytesToWrite, numBytesWritten, stringBuffer);
    }
    return true;
}
/*****************************************************************************
ReadParameters
*/
void ReadParameters(int argc, char** argv)
{
    if (argc >1)
    {
        for (int i=1; i <argc; i++)
        {
            if ( (argv[i][0] == '-') || (argv[i][0] == '/') )
            {
                switch (tolower(argv[i][1]))
                {
                case 'd':
                {
                    gFileName = argv[i+1];
                    i++;
                }
                break;
                default:
                {
                    Usage(argv[0]);
                }
                break;
                }
            }
            else
            {
                Usage(argv[0]);
            }
        }
    }
}
/*****************************************************************************
Usage
*/
void Usage(char *progname)
{
    fprintf(stderr,"Usage\n%s -d [database file name]\n", progname);
    Sleep(3000);
    exit(INCORRECT_PARAMETER_LIST_ERROR);
}
/*****************************************************************************
SetupPlexon
*/
bool SetupPlexon(HANDLE* hServerPollEvent)
{
    //** connect to the server
    PL_InitClientEx3(0, NULL, NULL);
    Sleep(20);
    //** check to make sure SortClient is running
    if (!PL_IsSortClientRunning())
    {
        printf("The SortClient is not running!...\r\n");
        return false;
    }
    //** open the event used to synchronize with the server
    *hServerPollEvent = OpenEvent(SYNCHRONIZE, FALSE, "PlexonServerEvent");
    return true;
}
/*****************************************************************************
CleanupPlexon
*/
void CleanupPlexon(HANDLE& hServerPollEvent, PL_Event* pServerTimeStampBuffer)
{
    CloseHandle(hServerPollEvent);
    PL_CloseClient();
    delete []pServerTimeStampBuffer;
}
/*****************************************************************************
ReadNumNeurons
*/
bool ReadNumNeurons(int& NumActiveUnits)
{
    ifstream input(gFileName);
    if (!input.is_open()) return false;
    input >> NumActiveUnits;
    input.close();
    printf("Using %d active neurons\n", NumActiveUnits);
    return true;
}
/*****************************************************************************
ReadInDatabase
*/
bool ReadInDatabase(double** Means, double** Stds, int NumActiveUnits, int* TrialCounters)
{
    int i,j;
    ifstream input(gFileName);
    if (!input.is_open()) return false;
    // skip NumActiveUnits on first line of file
    input.ignore(255,'\n');
    // Read in old trial counters
    for (i = 0; i < NUM_DIRECTIONS; i++)
    {
        input >> TrialCounters[i];
    }
    // Read in means and standard deviations
    for (j = 0; j < NumActiveUnits; j++)
    {
        for (i = 0; i < NUM_DIRECTIONS; i++)
        {
            input >> Means[j][i] >> Stds[j][i];
        }
    }
    input.close();
    printf("\nRead In Database:\n");
    for (i = 0; i < NUM_DIRECTIONS; i++)
        printf("%d_trials=%d\n", i, TrialCounters[i]);
    for (j = 0; j < NumActiveUnits; j++)
    {
        printf("Neuron (%d): \n");
        for (i = 0; i < NUM_DIRECTIONS; i++)
        {
            printf("%d_mean=%g, %d_std=%g, \n", i, Means[j][i],i, Stds[j][i]);
        }
        printf("\n");
    }
    printf("\n");
    for (j = 0; j < NumActiveUnits; j++)
    {
        for (i = 0; i < NUM_DIRECTIONS; i++)
        {
            Means[j][i] = 0.0;
            Stds[j][i] = 0.0;
        }
    }

    return true;
}
/*****************************************************************************
DisplayProbabilities
*/
void DisplayProbabilities(double* ReachProb)
{
    printf("\nReach Probabilities:\n");
    for (int j=0; j<NUM_DIRECTIONS; ++j)
    {
        printf("%g\n",ReachProb[j]);
    }
    printf("\n");
}
/*****************************************************************************
FindMostProbableDirection
Set index to j+1 because reach directions should be 1-8, not 0-7
*/
int FindMostProbableDirection(double* ReachProb, double& Probability)
{
    double  MostProb = -1e100;
    int index=-1;
    for (int i=0; i<NUM_DIRECTIONS; ++i)
    {
        MostProb = max(MostProb,ReachProb[i]);
    }
    for (int j=0; j<NUM_DIRECTIONS; ++j)
    {
        if (ReachProb[j] == MostProb)
        {
            index = j+1;
            break;
        }
    }
    printf("\nDirection index is %d \n",index);
    Probability=MostProb;
    return index;
}
/*****************************************************************************
SelectDirections
shiyan
*/
void SelectDirections(double** Means, double** Stds, int NumActiveUnits)
{
    double distance[NUM_DIRECTIONS][NUM_DIRECTIONS] = {0.0};
    int pairindex = -1;
    int i,j,k;
    double dummydistance = 0.0;
    for (k=0; k < NumActiveUnits; k++)
    {

        for (i=0; i<(NUM_DIRECTIONS-1); ++i)
        {
            for (j=(i+1); j<NUM_DIRECTIONS; ++j)

            {
                if (Stds[k][i]!=0.0 && Stds[k][j]!=0.0)
                {
                    distance[i][j] = pow(Means[k][i]-Means[k][j],2.0)/pow(Stds[k][i]-Stds[k][j],2.0 );
                }
            }
        }
        for (i=0; i<(NUM_DIRECTIONS-1); ++i)
        {
            for (j=(i+1); j<NUM_DIRECTIONS; ++j)
            {
                dummydistance = max(dummydistance,distance[i][j]);
            }
        }
        for (i=0; i<(NUM_DIRECTIONS-1); ++i)
        {
            for (j=(i+1); j<NUM_DIRECTIONS; ++j)
            {
                if (distance[i][j] == dummydistance)
                {
                    printf("\nDirection Pair are %d %d \n",i,j);
                }
            }
        }
    }
}

/*****************************************************************************
CalculateProbabilities
*/
void CalculateProbabilities(int* SpikeSum, double* LFP, double* ReachProb, double** Means, double** Stds,
                            long StartTime, long EndTime, int NumActiveUnits, int* TrialCounters)
{
    //variables for numerical error reduction:
    double  length = (double)(EndTime - StartTime)/40000.;
    //printf("start: %d, end %d\n", StartTime, EndTime);
    //this is the only place where we convert into timstamps into seconds
    printf("length %g\n", length);
    int i,j;
    int counter =0; //counts up the number of directions
    double Prior[NUM_DIRECTIONS]={0.0001};
    double Norm = 0;
    double Norm2 = 0.0;
    double minprob =-745;
    double maxprob =709;
    double correction = 0;
    double largest =0;
    double nextlargst = 0;
    double rpjmc=0;
    for (j=0; j<NUM_DIRECTIONS; ++j)if (TrialCounters[j]!=0)counter++;

    for (j=0; j<NUM_DIRECTIONS; ++j)

    {
        if (TrialCounters[j]!=0)Prior[j]=1.0/(double)counter;
        for (i=0; i < NumActiveUnits; i++)
        {
            //the log likelihood shifted by a constant

            if (Means[i][j]==1)Means[i][j]=1.000001;
            if (Means[i][j]==0)Means[i][j]=1.000001;
            //printf("SpikeSum[i] is %d", SpikeSum[i]);
            ReachProb[j] += (((SpikeSum[i]/length)*log(Means[i][j]))-Means[i][j]);
            //printf("sum: %d, length: %g, mean: %g\n", SpikeSum[i], length, Means[i][j]);
            //printf("undivided %g\n", ReachProb[j]);

        }
        printf("Direction %d has unnormed, unexp'd prob. %e\n", j, ReachProb[j]);
        if (ReachProb[j]>maxprob)
        {
            maxprob=ReachProb[j];
            correction = maxprob;
        }
        if (ReachProb[j]<minprob)
        {
            minprob=ReachProb[j];
            correction=minprob;
        }
        //the unlogged, correction shifted norm


    }
    printf("correction is: %g\n", correction);
    if (maxprob>709 && minprob<-745) printf("PROBABILITY RANGE ERROR");

    for (j=0; j<NUM_DIRECTIONS; j++)  if (fabs(ReachProb[j])>.01*fabs(correction))
            Norm+=exp(ReachProb[j]-correction);


    printf("norm is: %g \n", Norm);
    for (j=0; j<NUM_DIRECTIONS; j++)
    {

        ReachProb[j] = ReachProb[j] - log(Norm)-(correction);
        printf("logprob= %g\n", ReachProb[j]);
        printf("priorprob= %g\n", Prior[j]);
        ReachProb[j] = Prior[j]*exp(ReachProb[j]);
        Norm2 += ReachProb[j];
        printf("reachprob= %g\n", ReachProb[j]);
        //printf("reachprob is %g\n", ReachProb[j]);

    }
    for (j=0; j<NUM_DIRECTIONS; j++)
    {

        ReachProb[j] = ReachProb[j]/Norm2;

    }
    DisplayProbabilities(ReachProb);
    rpjmc=exp(710);
    printf("sanity: %g", rpjmc);
}
/*****************************************************************************
CalculateXY
*/
void CalculateXY(int index, double& x, double& y)
{
    x = 10.0*cos(((double) index)/((double)NUM_DIRECTIONS)*2.0*PI);
    y = 10.0*sin(((double) index)/((double)NUM_DIRECTIONS)*2.0*PI);
}
/*****************************************************************************
WriteOutData
        Not writing out XY data at this point.
*/
void WriteOutDecodeData(HANDLE hComm, int index, double Probability, double x, double y)
{
    char xBuffer[32];
    char yBuffer[32];
    char xyBuffer[32];
    // char indexBuffer[32];
    // char MostProbBuffer[32];
    char OutBuffer[64];
    int bufferlength;
    int prob;

    prob=(int)floor(Probability*100);
    bufferlength = sprintf(OutBuffer,"%d\n%d", index, prob);
    sprintf(xyBuffer,"%4.2f\n%4.2f\n",x,y);
    sprintf(xBuffer,"%4.2f\n",x);
    sprintf(yBuffer,"%4.2f\n",y);
    WriteSerialData(hComm, OutBuffer);

    printf("Wrote Serial data: %s (%d bytes)\n", OutBuffer, bufferlength);
}