Turbulence modeling code.

[schlafly]

Here's a PR implementing the turbulence modeling code. It contains a few improvements on what was presented on Monday, though none of these were found to make a significant impact on performance and so are now dead code (using a rational quadratic kernel rather than a gaussian kernel led the code to find a very high alpha, in which limit the rational quadratic is Gaussian; adopting a gradient-of-wavefront covariance instead of an independent Gaussian covariance only led to a fraction of a micron RMS different turbulence solutions).

The expected entrance point is newx, newy = turbulence.correct(x_measure, y_measure, x_expect, y_expect) following Sergey's example. That takes a bunch of "default" options throughout the code.

The most expensive step is solving for the parameters of the covariance matrix (the measurement + positioner offset noise floor; the scale of the turbulence; the amplitude of the turbulence). That requires the covariance matrix be constructed a lot of times and cholesky-decomposed, and won't scale at present to 5k fibers efficiently. But since it's only important to get the 3 parameters defining the covariance matrix, I've just decided to send it only the central 500 fibers provided, which means that that step takes ~3 sec.

The application of the covariance matrix to solving for the location of all the fibers still requires one matrix inversion. There are tricks to avoid this, but I'm already using a trick to get the "leave one out" solutions for each fiber. I'm not sure that avoiding the inversion would help if I then have to do a separate solution for each fiber. For 5k fibers, that inversion takes 3 sec on a NERSC login node using 12 cores (mkl.set_num_threads(12)), or 12 sec on a single core.

[julienguy]

Thanks. Two comments before testing.

• optional arguments dx,dy of the function correct() are not used
• would you be able to use a prior knowledge on the uncertainty on the offsets we want to measure? this would be an easy way to include fiducials or positioners we know are not moving.

[schlafly]

I wanted to have the interface for incorporating dead fibers & fiducials, even though I don't use that information yet. Thus the unused dx & dy.

It would be straightforward to populate the diagonal of the covariance matrix with actual uncertainties. Currently it's populated with a constant. I'm less enthusiastic about it, but it would also be straightforward in principle to compute separate empirical uncertainties for the moving and non-moving positioners.

[julienguy]

Tested with success. Will now add as an option to desi_fvc_proc.

Please add a docstring with a few lines of explanation and a description of the arguments and output to the function correct(..), and at least print a 'not implemented' warning or error message if someone is setting dx=... dy=....

[julienguy]

Thanks. Merging.