Friday, August 27, 2010

IIR Filtering Ideas

This post is going to be in pseudo-formal speak since I'm fleshing out ideas for my dissertation. Although that's tough without adding any equations.

The problem of constructing the optimal IIR controller for a given closed-loop plant and disturbance model can be solved as an LQR problem with a particular state-space system. However, the condition that the plant transfer matrix is commutable with the filter turns out to be overly restrictive. Suppose the transfer matrix can be factored into the product of minimum and non-minimnum phase matrices. We can define a new filter which is the product of the minimum phase component and the optimal filter F, this leaving the non-minimum phase component to be compensated. The requirement is now that this non-minimum phase component is commutable, which happens if its equal to a scalar transfer function times the identify matrix. If this is satisfied, we can perform the LQR problem to identify the controller, then multiply it by the inverse of the minimum phase component to recover the actual optimal filter that is implemented in the software.

For the adaptive optics experiment, the lack of significant DM dynamics simplifies the problem. If the non-minimum phase transfer matrix consists solely of n-step delays on the diagonal, then the optimal IIR filter is simply the n-step Kalman predictor for the disturbance model. The disturbance model is itself in innovations form, thus the Kalman predictor can be constructed directly from the state-space model generated by the subspace identification algorithm.

...some math showing how this is done....

The result is a filter which predicts the disturbance wavefronts n-steps ahead on the basis of the current wavefront measurement.

Beautiful.

Monday, August 23, 2010

Overly Ambitious

Now that the FIR filter is basically working for the single channel experiment, I think its a good time to set some targets for the next four weeks. The clear next step is to start working on a multichannel version using more modes. Theoretically there isn't that much difference here, but I suspect there will be practical problems with actuator saturation and other shittyness that will slow things down.

I think a reasonable goal is a 10 channel adaptive/optimal filter in 4 weeks. Here are some things that will need to happen:

1. Characterize closed-loop transfer matrix. How similar are the diagonal terms? How significant are the off diagonals? What's the best scalar transfer function approximation?

2. Get the adaptive controller working given an identified or ideal transfer matrix.

3. Write a script to compute the optimal multichannel FIR and IIR filters.

4. Get target camera working to compute Strehl ratios.

5. Find/write/steal a simulink block that can implement multichannel transfer functions. The pole at 1 in the pure integrator I'm using now might (will) lead to saturation when more modes are controlled.

6. Questions: How to the PSD's of the output channels compare? How does changing the number of modes alter steady-state performance?

A lot of this will require running many experiments or simulations, so there should be plenty of down time to pursue some theoretical stuff for the SISO case. Basically, I'd like to write m-files that do the following:

1. Compute the optimal IIR filter. How does performance compare to the FIR case?

2. Compute the optimal FIR gains using and RLS array algorithm.

3. Compute the optimal FIR gains using an RLS lattice filter.

4. Do something about implementing a Laguerre filter.

SInce these will all be m-files, I don't expect to implement these in the actual experiment right away. Mainly, I want to get some idea of the theoretical performance benefits by using increasingly complicated methods. Obviously I want to do all of this for the multichannel case one day, but the details of that are so complicated my head might explode first. We'll see.

Wednesday, August 18, 2010

Phucking Transpose

Yes, a freaking missing apostrophe was responsible for nearly a weeks delay. To try to narrow down the problem with the impulse response filter calculation, I was trying to match the results with the data driven m-file using made up disturbance models comprised of random matrices. It turned out that both filters were the same as long as the A matrix was diagonal. The only place in the code where this mattered were locations where I needed A transpose. The function call to Matlab's dlyap function was the freaking culprit. That piece of shit is one of those functions that's screwed me in the past, and of course the one time I wasn't careful it bits me in the ass.

Anyway, with the fix the impulse response method now returns basically the same FIR filter coefficients as using data. I also realized that the results from the adaptive loop I posted yesterday were crap. Somehow I was using the wrong model for the closed loop plant, chalk that up to shitty variable names. Here are the proper PSD's and modal outputs


My advisor thinks these results are stellar. Amazingly the PSD with the optimal FIR filter is pretty similar to the AO PSD, hopefully showing that my idea of dividing out the part of the disturbance cancelled by the classical loop is correct. In this case both the adaptive and fixed gain FIR filter are using 4 taps, so you might ask why the AO loop does better than the "optimal" filter at certain frequencies. The reason is that the adaptive loop can compensate somewhat for modeling error.

Right now I'm trying to run things without the shitty DM pause, which speeds things up to around 20Hz. I suspect the results won't be so peachy, but the time savings would be huge (8 minutes vs 40 for 10000 frames), and I wouldn't have to spend so much time quality time with youtube waiting for my experiments to run.

Also, I'd like to look at how the number of filter coefficients changes the steady-state performance, although I don't think adding many more taps will make much difference. Also I'd be cute to have the Strehl ratio performance to look at too.

Peace out.

Monday, August 16, 2010

Case of the Mondays

Here are the first results from using the "optimal" FIR filter, computed using the data-drive approach. Surprisingly it performs pretty well compared to the AO loop, the fact that its not spitting out absolute crap is a small miracle.



Here's a sampling of the modal output


I'm not sure what's going on with the AO loop. Looking at the commands with the adaptive loop closed shows lots of lower-end saturation going on about 1/3 of the time, so something is probably screwed up somewhere in the experiment. Hard to say at the moment since I'm running this remotely from home. The frustrating part is that it takes to friggin long to run an experiment, around 40 minutes for 10000 frames, that its easy to distract my already OCD mindset. I'm going to have to start doing this in simulations first.

This is good for a first step, but there are still some outstanding questions I'd like to look at this week. The first few deal with this optimal FIR filter calculation:

1. Determine what's really causing the difference between the data and impulse response drive methods to finding the optimal gains.

2. What's the real disturbance model that should be used in the calculation? What's the difference between computing it and identifying it from i/o data? Should it be SISO or MISO?

3. How does the filter order affect performance?

4. Write a script to determine the optimal IIR filter by solving an LQR problem.

Also, all of this stuff so far has been for the first focus mode. Sooner or later I'm going to have to do everything over again al MIMO, so I'd be nice to have some heads up if there are potential problems in the road. The first step is to look at the transfer matrix for multiple modes with the classical loop closed. Everything depends on this being diagonal with the same SISO tf on the diagonals. If this doesn't hold to a reasonable extent then there could be serious limitations. With that in mind:

1. How similar are the diagonal transfer functions for each mode? Models identified with significant saturation are garbage.

2. If they're all of the same form, but with a different gain, can the transfer matrix be factored into a single transfer function times a static gain matrix? If so, can this matrix just be incorporated into the poke matrix?

3. What's the difference between doing a MIMO subspace ID and multiple SISO id's?

4. Does the simulation even have enough accuracy to identify the model for multiple channels?

All this will be much faster if I just suck it up and do it in silico first.

8.16.10

Still having trouble getting decent results with the FIR filter. I still can't get the impulse response method to agree with what the least-squares solution spits out, even when I include the noise covariance in the state-space model. I might try a third method, there the filter coefficients are spit out from a finite-time LQR problem.

One thing I noticed is that the disturbance model I use isn't exactly the state-space system identified directly from the disturbances. Instead, its the part of the disturbance left over after going through the classical control loop. After untangling the block diagram, you end up dividing the original disturbance system by some transfer function involving the plant model. The problem I have is that this transfer function might not be exactly minimum phase, so you can end up with an unstable disturbance model to put into your FIR calculation.

I think its time to visit my advisor.

Tuesday, August 10, 2010

8.10.10 [2]

I feel like I'm loosing the script here, so a quick sitrep on what's going on: I currently have 2 different methods written to calculate the optimal FIR disturbance rejection filter. One uses state-space models of the disturbance and plant to generate impulse response sequences. These can be used to form a Weiner-Hopf problem and solved for the optimal coefficients a la earlier work we did on jitter control. The second method uses the models to simulate data directly. Given enough samples, the solution to another (similar) linear equation yields the coefficients.

Theoretically, both of these should give similar results (I think). But so far no luck. Here's a comparison of the output PSD with filters calculated using each method compared to using no filter (F=1).



I used the actual plant disturbance models I identified from the experiment. Using the data approach works pretty well, although its pretty cumbersome, and would be stupid with multiple modes. The method using the impulse responses, however, is just crap, clearly making things worse.

Confusingly, both methods crap out the identical filter with less complicated disturbance models. I think the problem is that the input noise covariance matrix isn't really accounted for in the state-space model of the disturbance. It comes into play in the data driven case since I have to use it to generate the input, but it doesn't show up directly in the impulse response as the moment. It should be easily, however, to incorporate it into the state-space model by multiplying the input matrix.

Thats the plan for the afternoon, as soon as I finish blogging here in a coffee shop.

8.10.10

After doing a system ID on the new disturbance model, I was finally able to get my "optimal" FIR filter working. Or rather, functional since in the experiment I get saturation and in the simulation I get crap. Right now I'm trying to write a similar script that finds the optimal coefficients by using simulated data directly in the least-squares problem, instead of just impulse response terms. We'll see if I get the same results as with the other method.

I can already see the next stop on this pain train. The adaptive filter basically solves this least-squares problem recursively using either LMS or RLS filter, so the obvious next step is to write my own adaptive code. The only piece missing is my ignorance about writing S-functions for Simulink, but there are hints that I might have to learn that eventually anyway.

Tuesday, August 03, 2010

One Thing Leads To Another

I think my plan to do a subspace ID on the actual disturbances as measured by the WFS is a good one, but it turned out that the code we have to do it requires the measurements to be on a square grid. Stupidly, I've been using a rectangular 19x22 measurement area on the WFS, with the DM area (found by adding all the influence functions) shifted to the right by around 5 columns. Putting disturbances from a state space model on this odd configuration required scaling the image, shifting it in an attempt to align it with the DM area, an only then performing a least-squares fit with the poke matrix. Naturally the results were shit, as shown in the video in 6.18.10.

My slope calculation code is ludicrously cumbersome, so it took a day to properly rewrite it to use a 19x19 grid (actually a down-sampled 38x38 grid) and make sure it was bug-free. Realigning the beam and removing tilt so that the area for both DMs was centered took another. My advisor commented that he doesn't know how I keep everything straight in my head with so much shit going on. I assured him I have no idea what I'm doing.

Anyway, its good to periodically tinker with the experiment and realign everything anyway to keep my monkey skills sharp. Today I put some state-space generated disturbances on the DM, which now only requires rescaling the state-space output to a slightly larger grid. The resulting phases look much better. Here's a comparison between (L to R) the state-space output (1 phase screen model on a 17x17 grid), the DM phase predicted using the poke matrix (on a 19x19 grid), and the actual phase measured by the WFS (with the bias removed of course).



If you squint you can actually see the similarity in the phases as it flows across the aperture. Either way it much cleaner than the random flatulence I was getting before.

Tomorrow I'll look at how this affects the performance of the adaptive loop. Now that the grid is square I can also go ahead with my original plan and compute an optimal FIR filter.