Communications Eye Diagram and Q-Factor procedures

/////////////////////////////////////////////////////////////////////////////////////////////
// Eye Diagram generation procedure for Igor Pro.               //
// Author: Paul Anthony Haigh, Northumbria University           //
// Contact: paul.haigh@northumbria.ac.uk                    //
/////////////////////////////////////////////////////////////////////////////////////////////

#pragma rtGlobals=1                     // Use modern global access method.

Function EyeDiagram(amplitude_vector,time_vector,nsamp)
//Declare Parameters
wave/Z amplitude_vector
wave time_vector
variable nsamp

// KillWave to avoid cross contamination of waves
killwaves/Z tvect_EyeDiagram
killwaves/Z amplitude_EyeDiagram

//Declare local variables
variable osamp
variable niterations
variable zeroed_iterations
variable i,j
variable amp_size
variable req_size
variable time_size
variable tsamp
// Generate new waves
Make/D/N=(numpnts(amplitude_vector))/O amplitude_EyeDiagram

amplitude_EyeDiagram=amplitude_vector                           // Copy amplitude_vector to avoid modification of original signal
DeletePoints 0, ceil(nsamp/2), amplitude_EyeDiagram             // Offset to centre of the bit - ensures a full 1-0 and 0-1 transition is displayed

osamp=2*nsamp;                                                  // Oversamp (osamp) by two (hence why offsetting by half a bit)
niterations=numpnts(time_vector)/osamp                         
Make/D/N=(niterations+nsamp)/O tvect_EyeDiagram

tsamp=time_vector[2]-time_vector[1]                             // Find sampling interval for time vector
for(i=0;i<=niterations;i+=1)   
    for(j=0;j!=osamp+1;j+=1)
        tvect_EyeDiagram[i*osamp+j]=j+1                         // Generate the time vector
    endfor
endfor
tvect_EyeDiagram=tvect_EyeDiagram*tsamp                         // Convert to bit rate relative time

// Insert zeros and convert them to NaNs to ensure no overlapping of bits
zeroed_iterations=niterations+niterations/osamp
for (i=0;i<=zeroed_iterations;i+=osamp+1)
    InsertPoints i, 1, tvect_EyeDiagram
    InsertPoints i, 1, amplitude_EyeDiagram
    if (tvect_EyeDiagram[i]==0)
        tvect_EyeDiagram[i]=NaN
        if (amplitude_EyeDiagram[i]==0)
            amplitude_EyeDiagram[i]=NaN
        endif
    endif
endfor

// Create conditions for variable sizes
req_size=zeroed_iterations
time_size=numpnts(tvect_EyeDiagram)
amp_size=numpnts(amplitude_EyeDiagram)

// Ensure waves are the same size using variable size conditions
if (amp_size!=req_size)
    DeletePoints req_size, amp_size-req_size, amplitude_EyeDiagram
    if (time_size!=req_size)
        DeletePoints req_size, time_size-req_size, tvect_EyeDiagram
    endif
endif

Display amplitude_EyeDiagram vs tvect_EyeDiagram
ModifyGraph tick=2,mirror=1,fStyle=1,fSize=14,standoff=0,width=425.197,height=283.465,lsize=2,rgb=(0,0,0)
Qfactor(amplitude_vector,nsamp)
end

/////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
/////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

/////////////////////////////////////////////////////////////////////////////////////////////
// Qfactor calculation procedure for Igor Pro.                  //
// Author: Paul Anthony Haigh, Northumbria University           //
// Contact: paul.haigh@northumbria.ac.uk                    //
/////////////////////////////////////////////////////////////////////////////////////////////

// Declare function and dependent variables
Function Qfactor(amplitude_vector,nsamp)
wave/Z amplitude_vector
variable nsamp

// KillWave to avoid cross contamination of waves
killwaves/Z amplitude_Qfactor
killwaves/Z amplitude_Qfactor_downsamp
killwaves/Z upper_vect
killwaves/Z lower_vect
killwaves/Z Q_factor
// Generate new waves
Make/D/N=(numpnts(amplitude_vector))/O amplitude_Qfactor
Make/D/N=(numpnts(amplitude_vector)/nsamp)/O amplitude_Qfactor_downsamp

// Procedure variables
variable i,j,k,q,ii,jj
variable mean_upper
variable mean_lower
variable std_upper
variable std_lower
variable Q_factor

amplitude_Qfactor=amplitude_vector                  // Copy amplitude_vector to avoid modification of original signal
WaveStats /M=1/Q/Z amplitude_Qfactor                    // Get WaveStats to find V_avg
amplitude_Qfactor=amplitude_Qfactor-V_avg               // Make data bi-polar

// Downsample in the centre of the signal - data must start from the 0-1 or 1-0 transition. Change the first condition to select a different bit.
for (j=ceil(nsamp/2);j<=floor(numpnts(amplitude_Qfactor)/nsamp);j+=nsamp)
    amplitude_Qfactor_downsamp[k]=amplitude_Qfactor[j]
    k+=1
endfor

WaveStats /M=1/Q/Z amplitude_Qfactor_downsamp       // Update WaveStats to find new V_avg

// Calculate mean and standard deviation of both levels
ii=0
jj=0
Make/D/N=(numpnts(amplitude_Qfactor_downsamp))/O upper_vect     // Generate new vectors for each level
Make/D/N=(numpnts(amplitude_Qfactor_downsamp))/O lower_vect

// Sort downsampled data into the two new vectors
for (q=0;q<=numpnts(amplitude_Qfactor_downsamp);q+=1)
    if (amplitude_Qfactor_downsamp[q]>=V_avg && amplitude_Qfactor_downsamp[q]!=0)
        upper_vect[ii]=amplitude_Qfactor_downsamp[q]
        ii+=1
    elseif (amplitude_Qfactor_downsamp[q]<V_avg && amplitude_Qfactor_downsamp[q]!=0)
        lower_vect[jj]=amplitude_Qfactor_downsamp[q]
        jj+=1
    endif
endfor

// Remove all points that follow the data to avoid corruption of the WaveStats
DeletePoints ii, numpnts(amplitude_Qfactor_downsamp), upper_vect
DeletePoints jj, numpnts(amplitude_Qfactor_downsamp), lower_vect

WaveStats /M=1/Q/Z upper_vect   // Find V_avg of upper_vect
mean_upper=V_avg                // Assign to a variable
WaveStats /M=2/Q/Z upper_vect   // Find V_sdev of upper_vect
std_upper=V_sdev                // Assign to a variable
KillVariables/Z V_avg               // Remove variables to ensure no cross contamination
KillVariables/Z V_sdev 

WaveStats /M=1/Q/Z lower_vect   // Find V_avg of lower_vect
mean_lower=V_avg                // Assign to a variable
WaveStats /M=2/Q/Z lower_vect   // Find V_sdev of lower_vect
std_lower=V_sdev                // Assign to a variable
KillVariables/Z V_avg               // Remove variables to ensure no cross contamination
KillVariables/Z V_sdev

// Calculate Q-Factor and print
Q_factor=(mean_upper-mean_lower)/(std_upper+std_lower)
printf "Q-Factor = %g\r", Q_factor

// KillWaves to ensure no cross contamination
KillWaves/Z upper_vect
KillWaves/Z lower_vect
End