Phase-space optics 2: Instantaneous correlation and pseudo-symbolic computation in MATLAB

As I noted in the previous post about FiO, I wish to share some of the simulations and basics about phase-space signal analysis. While doing simulations, I sometimes find it tricky to interpret analog equations and convert them to a sequence of programmed operations on discrete data. Also, one needs to be careful about not losing information when going from analog to discrete domain. This all-important step is known as A to D (analog to *discrete*) conversion and is responsible for Nyquist and Shannon being so famous.

I rely on two tools for my modeling needs: MATLAB (TM) most of the time and Mathematica (TM) sometimes. Mathematica is excellent for allowing quick visualization of the model due to its ability to work directly with symbols, whereas MATLAB is good for fast and maintainable implementation. It may be possible to write maintainable code in Mathematica but that requires going through steep learning curve, especially for folks who had learned classic ‘C’ as their first language. Is there a way of marrying symbolic ease of Mathematica with programming ease of MATLAB? I seem to have found one using anonymous functions. The code used in this post has been posted here on MATLAB file-exchange, which will be available within a day after moderation.

Let us consider the problem of calculating instantaneous correlation (which is a primary quantity required to compute and understand the Wigner distribution and Ambiguity function). For a signal it is defined as,

We know that the auto-correlation of a function is defined as Thus, we see that instantaneous correlation is in some sense the quantity which when integrated along gives us the auto-correlation function.

As we shall see in upcoming posts, Wigner distribution is the Fourier transform of the instantaneous correlation (Icorr) along the delay () dimension, whereas the Ambiguity function is Icorr’s Fourier transform along the normal () dimension. Therefore, there are two ways of interpreting what the instantaneous correlation is:

1. is the signal which is shifted by and The conjugate of the delayed version () when multiplied with the advanced version () gives us Icorr.
2. is the signal which is first dilated to give us This dilated signal’s conjugate is flipped to give The signal is advanced to give and the signal is delayed to give us

As is often the case – hindsight is 20/20. I could write above paragraph so clearly because I have already gone through simulating the problem with two tools. But things should be presented logically and not chronologically, so I explained the model first. Next are the simulations.

Let’s consider results from Mathematica first. Note that blue color in all figures corresponds to green to and red to their product.

• Define a chirp:

   Plot[{f[x] }, {x, -0.5, 1.5}, PlotRange -> 1]

  As shown below this function is a chirp that is ‘on’ between and The signal’s slow and fast ends help us identify the orientation and extent as it is modified.

  

• Treat the signal as a function of and see how the ‘early’ and ‘late’ signals change as is changed:

   Manipulate[
   Plot[{f[x + (c/2)], f[x - (c/2)], f[x + c/2]*f[x - c/2]}, {x, -2, 2},
   PlotStyle -> {{Blue, Dashed}, {Green, Dashed}, Red},
   PlotRange -> 1], {{c, 0.5}, -1, 1, 0.1}]

  This animation, when executed in Mathematica, immediately reveals the interpretation-1 given above. Following figure shows a snapshot of all three signals at which shows that the early signal starts at and the late signal at More importantly, we notice that the Icorr has support (i.e., it is non-zero) along the -axis (on either size of zero) up to the duration of the signal along

  

• Treat the signal as a function of and notice the effect of changing :

   Manipulate[
   Plot[{f[x + c/2], f[x - c/2], f[x + c/2]*f[x - c/2]}, {c, -2, 2},
   PlotStyle -> {{Blue, Dashed}, {Green, Dashed}, Red},
   PlotRange -> 1], {{x, 0.5}, 0, 1}]

  Above animation explains the interpretation-2 given above. Important thing to notice is that the dilated signal is shifted by rather than along the dilated axis. Following figure is a snapshot for and we notice that the signal starts earlier by whereas the flipped signal starts late by – leading to their perfect overlap.

  

Now let’s try our hand at MATLAB. Typical MATLAB computation involves constructing an independent variable as a vector (or a matrix) and using it to calculate dependent variables. In the above case, we are manipulating the independent variable () and (). Therefore, to calculate the modified signals, we must clearly understand what above manipulations actually do. This is the chicken-and-egg problem that I often encounter in absence of symbolic computation. Now let us see how anonymous functions can help us. Anonymous functions are defined as handles to an expression. Let’s walk through MATLAB code that produces the same results as above.

• Instead of defining independent and dependent variables, straight away write out an anonymous function. Use this anonymous function to compute and plot the signal:

   sigh=@(xi) (cos(10*pi*xi.^2).*(0<=xi & xi<=1));
   	 figure(1);
   	 subplot(311);
   	 plot(linspace(-0.5,1.5,1024),sigh(linspace(-0.5,1.5,1024)));
   	 title('Signal');

• Define an independent variable , use it to compute dependent signals and see how the ‘early’ and ‘late’ signals change as is changed:

   x=linspace(-2,2,1024);
   for c=-1:0.1:1
   figure(2); plot(x,sigh(x+c/2),'b--',x,sigh(x-c/2),'g--',x,sigh(x+c/2).*sigh(x-c/2),'r');
   legend({'Early','Late','Icorr'},'Location','West');
   title(['c=' num2str(c)]);
   pause(0.1);
   if(c==0.5)
   figure(1);
   subplot(312);
   hplot=plot(x,sigh(x+c/2),'b--',x,sigh(x-c/2),'g--',x,sigh(x+c/2).*sigh(x-c/2),'r');
   title(['Signal assumed a function of x: c=' num2str(c)]);
   legend({'Early','Late','Icorr'},'Location','West');
   end
   end

• Define an independent variable use it to compute dependent signals and see the effect of

   c=linspace(-2,2,1024);
   for x=0:0.1:1
   figure(2);
   plot(c,sigh(x+c/2),'b--',c,sigh(x-c/2),'g--',c,sigh(x+c/2).*sigh(x-c/2),'r');
   legend({'Early','Late','Icorr'},'Location','West');
   title(['x=' num2str(x)]);
   pause(0.1);
   if(x==0.5)
   figure(1);
   subplot(313);
   hplot=plot(c,sigh(x+c/2),'b--',c,sigh(x-c/2),'g--',c,sigh(x+c/2).*sigh(x-c/2),'r');
   title(['Signal assumed a function of c: x=' num2str(x)]);
   legend({'Early','Late','Icorr'},'Location','West');
   end
   end

  Following image shows the results of the above (which are the same as those obtained with Mathematica.)

  

From the above, we see that once anonymous function is defined, we can call it in a ‘nearly symbolic’ fashion and interpret the model. You might have noticed arcane (yet succinct and semantically accurate) style of coding in Mathematica vs. simpler and easy to interpret style of coding in MATLAB.

For sure, this is not the only way of obtaining discretized data from the model. Mathematica no doubt allows you to export your data in discrete format. MATLAB has its own symbolic toolbox that may allow discretization. But both of these approaches did not seem as straight-forward and efficient to me as using the anonymous functions.

---

Typeset using LaTeX2WP script by Luca Trevisan.