Introduction to Sampling

Mathematica definitions

Load needed packages for this notebook

<< Graphics`Colors`;
<< Miscellaneous`Audio`;

Set up defaults

Turn off annoying warning messages

Off[General::spell];
Off[General::spell1];

Set default options for specific functions

SetOptions[Plot, Frame -> True];
SetOptions[ListPlot, Frame -> True,  PlotJoined -> True,  
        Axes -> {True, False}];
SetOptions[Fourier, 
        FourierParameters -> {1, -1}];
SetOptions[InverseFourier, 
        FourierParameters -> {1, -1}];

Visualization functions

Fix time variable for sampled signals

FixTime[signal_, 
          sampleRate_] := 
Module[{time, 
            data}, 
time = 
            Table[i/sampleRate, {i, 1, 
                Length[signal]}]; 
data  =  {time, 
                signal} // Transpose
];

Generate periodic waveform from one period

PeriodicSignal[f, {dmin, dmax}, {t, min, max}] generates a periodic waveform for the interval [min, max]. It is assumed that one period of the function f(t) is defined over the interval [dmin, dmax].

PeriodicSignal[
          f_, {dmin_, dmax_}, {t_, min_, 
            max_}] := 
Module[{T, imin, 
            imax}, 
T  =  
            dmax  -  dmin; 
imin = 
            Floor[(min - dmin)/T]; 
imax = 
            Ceiling[(max - dmax)/T];  
Sum[
            f[t - i T] (UnitStep[t - dmin - i T]  -  
                  UnitStep[t - dmax - i T]), {i, imin, 
              imax}]
];

Plot periodic waveforms

PeriodicPlot[f, {dmin, dmax}, {t, min, max}, opts] plots a a periodic waveform for the interval [min, max]. It is assumed that one period of the function f(t) is defined over the interval [dmin, dmax]. The 'opts' argument passed plotting options to the Mathematica function Plot[ ].

PeriodicPlot[f_, {dmin_, dmax_}, {t_, min_, max_}, 
          opts___] := 
Module[{fnc}, 
fnc = PeriodicSignal[f, {dmin, dmax}, {t, min, max}]; 
Plot[
            fnc, {t, min, max}, opts]
];

Plot discrete-time signal as sequence of amplitude lines

DiscreteListPlot[signal_, color_:  Black, 
          opts___] := 
Module[{data, g1, 
            g2}, 
If[
            Length[signal[[1]]]  ==  2, 
data = 
              signal, 
data = 
              FixTime[signal, 
                1]
]; 
g1 = Line[{{#[[1]], 0}, {#[[1]], #[[2]]}}] &  /@  
              data; 
g2 = {PointSize[ .015], 
                  Point[#]} &  /@  data; 
Show[
            Graphics[{color, g1, g2}], 
            opts]
];

Filtering definitions

FourierRearrange[fft] rearranges the result of an FFT (Mathematica function Fourier[ ]), to a easier to visualize the spectrum represented by 'fft'. The length of the list 'fft' is assumed to be even; if it is not, the last term in 'fft' is dropped.

FourierRearrange[fft_,  
          samplingRate_] := 
Module[{fft2, n, 
            n2}, 
If[EvenQ[Length[fft]],  
            fft2  =  fft,  
            fft2  =  Drop[fft, -1]];  
n = 
            Length[fft]; 
          n2 = n/2; 
{Table[ 
                samplingRate i/n, {i, -n2 + 1, n2}], 
              Take[fft2, {n2 + 2, n2*2}]~Join~Take[fft2, {1, n2 + 1}]} // 
            Transpose
];

Filter[n, sampleRate, {lowFreq, highFreq}] generates an ideal discrete filter for an FFT list of length 'n' (sampled at 'sampleRate'), with a bandpass from 'lowFreq' to 'highFreq' (in Hz).

Filter[n_, 
          samplingRate_, {lowF_, 
            highF_}] := 
Module[{coeff, l, 
            h}, 
coeff = 
            Table[0, {n}]; 
l = Floor[lowF*n/samplingRate]; 
          h = Ceiling[highF*n/samplingRate]; 
If[
            l  ==  0,  l  =  1;  
            coeff = ReplacePart[coeff, 1, 1]]; 
coeff = 
            ReplacePart[coeff, 1, Table[{i}, {i, l + 1, h + 1}]]; 
          coeff = ReplacePart[coeff, 1, 
              Table[{i}, {i, n + 1 - h, 
                  n + 1 - l}]]
];

Signal reconstruction

Sinc[t_] := Sin[t]/t;

ReconstructSignal[t_, fs_, range_: {0, 5},  
          padding_:  1] := 
Plus  @@  
          Table[x[n, fs] Sinc[π fs (t - n/fs)], {n, 
              Floor[(range[[1]] - padding)*sampleFrequency], 
              Ceiling[(range[[2]] + padding)*
                  sampleFrequency]}];

Sampling visualization

Continuous-time function

xc[t_] := Cos[2  π t];

range = {0, 5};

gin = 
        Plot[xc[t], {t, range[[1]], range[[2]]},  
          PlotStyle -> {Blue,  Thickness[ .01]}];

g1 = 
        Plot[xc[t], {t, range[[1]], range[[2]]},  
          PlotStyle -> {Blue,  Dashing[{ .02}]}];

Sampler (CT to DT converter)

x[n_,  fs_] :=  xc[n/fs];

sampler[fs_, 
          range_: {0, 5}] := 
Table[{n/fs, x[n, fs]}, {n, 
            range[[1]]*fs, range[[2]]*fs}];

Sampled signal: How much is enough?  (sampling frequency = 100Hz)

sampleFrequency = 100;

xd  =  sampler[sampleFrequency,  range];
g2 = 
        DiscreteListPlot[xd, Red,  Frame -> True, 
          Axes -> {True, False}, PlotRange -> All];

Show[g1, g2];

Sampled signal: How much is enough?  (sampling frequency = 10Hz)

sampleFrequency = 10;

xd  =  sampler[sampleFrequency,  range];
g2 = 
        DiscreteListPlot[xd, Red,  Frame -> True, 
          Axes -> {True, False}, PlotRange -> All];

Show[g1, g2];

Sampled signal: How much is enough?  (sampling frequency = 5Hz)

sampleFrequency = 5;

xd  =  sampler[sampleFrequency,  range];
g2 = 
        DiscreteListPlot[xd, Red,  Frame -> True, 
          Axes -> {True, False}, PlotRange -> All];

Show[g1, g2];

Sampled signal: How much is enough?  (sampling frequency = 3Hz)

sampleFrequency = 3;

xd  =  sampler[sampleFrequency,  range];
g2 = 
        DiscreteListPlot[xd, Red,  Frame -> True, 
          Axes -> {True, False}, PlotRange -> All];

Show[g1, g2];

Sampled signal: How much is enough?  (sampling frequency = 2.1Hz)

sampleFrequency = 2.1;

xd  =  sampler[sampleFrequency,  range];
g2 = 
        DiscreteListPlot[xd, Red,  Frame -> True, 
          Axes -> {True, False}, PlotRange -> All];

Show[g1, g2];

Sampled signal: How much is enough?  (sampling frequency = 1.9Hz)

sampleFrequency = 1.9;

xd  =  sampler[sampleFrequency,  range];
g2 = 
        DiscreteListPlot[xd, Red,  Frame -> True, 
          Axes -> {True, False}, PlotRange -> All];

Show[g1, g2];

Sampled signal: How much is enough?  (sampling frequency = 1.4Hz)

sampleFrequency = 1.4;

xd  =  sampler[sampleFrequency,  range];
g2 = 
        DiscreteListPlot[xd, Red,  Frame -> True, 
          Axes -> {True, False}, PlotRange -> All];

Show[g1, g2];

Sampled signal: How much is enough?  (sampling frequency = 1.1Hz)

sampleFrequency = 1.1;

xd  =  sampler[sampleFrequency,  range];
g2 = 
        DiscreteListPlot[xd, Red,  Frame -> True, 
          Axes -> {True, False}, PlotRange -> All];

Show[g1, g2];

Ideal reconstruction from samples

Ideally reconstructed signal (sampling frequency = 10Hz)

sampleFrequency = 10;
reconstruction = 
        ReconstructSignal[t, sampleFrequency,  range,  10];
xd  =  sampler[sampleFrequency,  range];
g2 = 
        DiscreteListPlot[xd, Red,  Frame -> True, 
          Axes -> {True, False}, PlotRange -> All,  
          DisplayFunction -> Identity];
g3 = 
        Plot[reconstruction,  {t, 0, 5},  
          PlotStyle -> {Orange, Thickness[ .01]}, 
          PlotRange -> All,  
          DisplayFunction -> Identity];
ga = 
        Show[g3, 
          DisplayFunction -> Identity];
gb = 
        Show[g3, g2,  
          DisplayFunction -> Identity];
gc  = 
        Show[g1, g2, 
          DisplayFunction -> Identity];
gd = 
        Show[g1, g3,  
          DisplayFunction -> Identity];
Show[GraphicsArray[{{ga, gb}, {gc, gd}}], ImageSize -> 640, 
        DisplayFunction -> $DisplayFunction];

Ideally reconstructed signal (sampling frequency = 5Hz)

sampleFrequency = 5;
reconstruction = 
        ReconstructSignal[t, sampleFrequency,  range,  10];
xd  =  sampler[sampleFrequency,  range];
g2 = 
        DiscreteListPlot[xd, Red,  Frame -> True, 
          Axes -> {True, False}, PlotRange -> All,  
          DisplayFunction -> Identity];
g3 = 
        Plot[reconstruction,  {t, 0, 5},  
          PlotStyle -> {Orange, Thickness[ .01]}, 
          PlotRange -> All,  
          DisplayFunction -> Identity];
ga = 
        Show[g3, 
          DisplayFunction -> Identity];
gb = 
        Show[g3, g2,  
          DisplayFunction -> Identity];
gc  = 
        Show[g1, g2, 
          DisplayFunction -> Identity];
gd = 
        Show[g1, g3,  
          DisplayFunction -> Identity];
Show[GraphicsArray[{{ga, gb}, {gc, gd}}], ImageSize -> 640, 
        DisplayFunction -> $DisplayFunction];

Ideally reconstructed signal (sampling frequency = 3Hz)

sampleFrequency = 3. ;
reconstruction = 
        ReconstructSignal[t, sampleFrequency,  range,  10];
xd  =  sampler[sampleFrequency,  range];
g2 = 
        DiscreteListPlot[xd, Red,  Frame -> True, 
          Axes -> {True, False}, PlotRange -> All,  
          DisplayFunction -> Identity];
g3 = 
        Plot[reconstruction,  {t, 0, 5},  
          PlotStyle -> {Orange, Thickness[ .01]}, 
          PlotRange -> All,  
          DisplayFunction -> Identity];
ga = 
        Show[g3, 
          DisplayFunction -> Identity];
gb = 
        Show[g3, g2,  
          DisplayFunction -> Identity];
gc  = 
        Show[g1, g2, 
          DisplayFunction -> Identity];
gd = 
        Show[g1, g3,  
          DisplayFunction -> Identity];
Show[GraphicsArray[{{ga, gb}, {gc, gd}}], ImageSize -> 640, 
        DisplayFunction -> $DisplayFunction];

Ideally reconstructed signal (sampling frequency = 2.1Hz)

sampleFrequency = 2.1;
reconstruction = 
        ReconstructSignal[t, sampleFrequency,  range,  10];
xd  =  sampler[sampleFrequency,  range];
g2 = 
        DiscreteListPlot[xd, Red,  Frame -> True, 
          Axes -> {True, False}, PlotRange -> All,  
          DisplayFunction -> Identity];
g3 = 
        Plot[reconstruction,  {t, 0, 5},  
          PlotStyle -> {Orange, Thickness[ .01]}, 
          PlotRange -> All,  
          DisplayFunction -> Identity];
ga = 
        Show[g3, 
          DisplayFunction -> Identity];
gb = 
        Show[g3, g2,  
          DisplayFunction -> Identity];
gc  = 
        Show[g1, g2, 
          DisplayFunction -> Identity];
gd = 
        Show[g1, g3,  
          DisplayFunction -> Identity];
Show[GraphicsArray[{{ga, gb}, {gc, gd}}], ImageSize -> 640, 
        DisplayFunction -> $DisplayFunction];

Ideally reconstructed signal (sampling frequency = 1.9Hz)

sampleFrequency = 1.9;
reconstruction = 
        ReconstructSignal[t, sampleFrequency,  range,  10];
xd  =  sampler[sampleFrequency,  range];
g2 = 
        DiscreteListPlot[xd, Red,  Frame -> True, 
          Axes -> {True, False}, PlotRange -> All,  
          DisplayFunction -> Identity];
g3 = 
        Plot[reconstruction,  {t, 0, 5},  
          PlotStyle -> {Orange, Thickness[ .01]}, 
          PlotRange -> All,  
          DisplayFunction -> Identity];
ga = 
        Show[g3, 
          DisplayFunction -> Identity];
gb = 
        Show[g3, g2,  
          DisplayFunction -> Identity];
gc  = 
        Show[g1, g2, 
          DisplayFunction -> Identity];
gd = 
        Show[g1, g3,  
          DisplayFunction -> Identity];
Show[GraphicsArray[{{ga, gb}, {gc, gd}}], ImageSize -> 640, 
        DisplayFunction -> $DisplayFunction];

Ideally reconstructed signal (sampling frequency = 1.4Hz)

sampleFrequency = 1.4;
reconstruction = 
        ReconstructSignal[t, sampleFrequency,  range,  10];
xd  =  sampler[sampleFrequency,  range];
g2 = 
        DiscreteListPlot[xd, Red,  Frame -> True, 
          Axes -> {True, False}, PlotRange -> All,  
          DisplayFunction -> Identity];
g3 = 
        Plot[reconstruction,  {t, 0, 5},  
          PlotStyle -> {Orange, Thickness[ .01]}, 
          PlotRange -> All,  
          DisplayFunction -> Identity];
ga = 
        Show[g3, 
          DisplayFunction -> Identity];
gb = 
        Show[g3, g2,  
          DisplayFunction -> Identity];
gc  = 
        Show[g1, g2, 
          DisplayFunction -> Identity];
gd = 
        Show[g1, g3,  
          DisplayFunction -> Identity];
Show[GraphicsArray[{{ga, gb}, {gc, gd}}], ImageSize -> 640, 
        DisplayFunction -> $DisplayFunction];

Ideally reconstructed signal (sampling frequency = 1.1Hz)

sampleFrequency = 1.1;
reconstruction = 
        ReconstructSignal[t, sampleFrequency,  range,  10];
xd  =  sampler[sampleFrequency,  range];
g2 = 
        DiscreteListPlot[xd, Red,  Frame -> True, 
          Axes -> {True, False}, PlotRange -> All,  
          DisplayFunction -> Identity];
g3 = 
        Plot[reconstruction,  {t, 0, 5},  
          PlotStyle -> {Orange, Thickness[ .01]}, 
          PlotRange -> All,  
          DisplayFunction -> Identity];
ga = 
        Show[g3, 
          DisplayFunction -> Identity];
gb = 
        Show[g3, g2,  
          DisplayFunction -> Identity];
gc  = 
        Show[g1, g2, 
          DisplayFunction -> Identity];
gd = 
        Show[g1, g3,  
          DisplayFunction -> Identity];
Show[GraphicsArray[{{ga, gb}, {gc, gd}}], ImageSize -> 640, 
        DisplayFunction -> $DisplayFunction];