Loading [MathJax]/extensions/tex2jax.js
mxlib
c++ tools for analyzing astronomical data and other tasks by Jared R. Males. [git repo]
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Modules Pages
influenceFunctions.hpp
Go to the documentation of this file.
1/** \file influenceFunctions.hpp
2 * \author Jared R. Males (jaredmales@gmail.com)
3 * \brief Utilities for generating and analyzing deformable mirror influence functions.
4 * \ingroup mxAO_files
5 *
6 */
7
8#ifndef __influenceFunctions_hpp__
9#define __influenceFunctions_hpp__
10
11#include "../improc/eigenCube.hpp"
12#include "../math/func/gaussian.hpp"
13#include "../ioutils/fits/fitsFile.hpp"
14
15#include <mx/math/eigenLapack.hpp>
16// #include <mx/gnuPlot.hpp>
17#include <mx/ioutils/pout.hpp>
18#include <mx/math/randomT.hpp>
19
20#include <vector>
21
22#include "aoPaths.hpp"
23
24namespace mx
25{
26
27namespace AO
28{
29
30template <typename realT>
31struct influenceFunctionGaussianSpec
32{
33 std::string dmName;
34 realT dmSz;
35 int linNAct;
36 realT diameter;
37 realT coupling;
38 realT pupilSz;
39 bool offsetOdd;
40
41 std::vector<int> badActi;
42 std::vector<int> badActj;
43 std::vector<int> badActType;
44
45 influenceFunctionGaussianSpec()
46 {
47 pupilSz = 0;
48 offsetOdd = false;
49 }
50
52 {
53 if( pupilSz <= 0 )
54 pupilSz = dmSz;
55
56 realT pupilScale = pupilSz / dmSz;
57
58 realT act_space = ( dmSz - 1.0 ) / linNAct;
59
60 std::vector<realT> xcoords, ycoords;
61
62 realT xc, yc, xcen, ycen;
63 realT _rad;
64
65 // The center of the pupil
66 xcen = 0.5 * ( pupilSz - 1 ); // 0.5*(linNAct*pupilScale - 1);
67 ycen = 0.5 * ( pupilSz - 1 ); // 0.5*(linNAct*pupilScale - 1);
68
69 realT xoff = 0.0;
70
71 int q = 0;
72 for( int i = 0; i < linNAct; ++i )
73 {
74 xc = ( i + 0.5 ) * act_space - 0.5 * ( dmSz - 1 ); // coordinate in pixels w.r.t. dm, not pupil
75
76 for( int j = 0; j < linNAct; ++j )
77 {
78 int skip = 0;
79
80 if( badActi.size() > 0 && badActi.size() == badActj.size() )
81 {
82 for( int k = 0; k < badActi.size(); ++k )
83 {
84 if( badActi[k] == i && badActj[k] == j )
85 {
86 skip = 1;
87 }
88 }
89 }
90
91 if( skip )
92 continue;
93
94 yc = ( j + 0.5 ) * act_space - 0.5 * ( dmSz - 1 ); // coordinate in pixels w.r.t dm, not pupil
95
96 if( offsetOdd )
97 {
98 if( j % 2 )
99 xoff = 0.0;
100 else
101 xoff = 0.5 * act_space;
102 }
103
104 _rad = sqrt( pow( xc + xoff, 2 ) + pow( yc, 2 ) );
105
106 if( _rad / act_space <= 0.5 * diameter + .01 || diameter == 0 )
107 {
108 xcoords.push_back( xcen + ( xc + xoff ) + 0.0 * ( dmSz - pupilSz ) );
109 ycoords.push_back( ycen + yc + 0.0 * ( dmSz - pupilSz ) );
110 ++q;
111 }
112 }
113 }
114
115 std::cout << "Number of actuators: " << q << " " << xcoords.size() << "\n";
116
117 /*gnuPlot gp;
118 gp.command("set size square");
119 gp.plot(xcoords, ycoords, "w p ps 0.75");
120 gp.point(xcen, ycen);
121 gp.circle(xcen, ycen, 0.5*(pupilSz-1), "w l", "Pupil");
122 */
123 realT fwhm = act_space * 2. * sqrt( log2( 2 ) ) / sqrt( -log2( coupling ) );
124
125 mx::improc::eigenCube<realT> act_inf( pupilSz, pupilSz, xcoords.size() );
126
127 for( int i = 0; i < xcoords.size(); ++i )
128 {
129 xc = xcoords[i];
130 yc = ycoords[i];
131
133 act_inf.image( i ).data(), pupilSz, pupilSz, 0.0, 1.0, xc, yc, fwhm2sigma( fwhm ) );
134 }
135
136 if( badActi.size() == badActj.size() && badActi.size() == badActType.size() )
137 {
138 Eigen::Array<realT, -1, -1> im;
139
140 std::cerr << "ok trying bad acts\n";
141 std::cerr << badActi.size() << "\n";
142 for( int k = 0; k < badActi.size(); ++k )
143 {
144 if( badActType[k] == 0 )
145 continue;
146
147 im.resize( pupilSz, pupilSz );
148
149 xc = ( badActi[k] + 0.5 ) * act_space - 0.5 * ( dmSz - 1 ); // coordinate in pixels w.r.t. dm, not pupil
150 yc = ( badActj[k] + 0.5 ) * act_space - 0.5 * ( dmSz - 1 ); // coordinate in pixels w.r.t dm, not pupil
151
152 if( offsetOdd )
153 {
154 if( badActj[k] % 2 )
155 xoff = 0.0;
156 else
157 xoff = 0.5 * act_space;
158 }
159
160 _rad = sqrt( pow( xc + xoff, 2 ) + pow( yc, 2 ) );
161
162 std::cerr << xc << " " << yc << "\n";
163
164 if( _rad / act_space <= 0.5 * diameter + .01 || diameter == 0 )
165 {
166 xc = ( xcen + ( xc + xoff ) + 0.0 * ( dmSz - pupilSz ) );
167 yc = ( ycen + yc + 0.0 * ( dmSz - pupilSz ) );
168
169 std::cerr << xc << " " << yc << "\n";
170
171 mx::math::func::gaussian2D( im.data(), pupilSz, pupilSz, 0.0, -150.0, xc, yc, fwhm2sigma( fwhm ) );
172
173 for( int i = 0; i < xcoords.size(); ++i )
174 {
175 act_inf.image( i ) += im;
176 }
177 }
178 }
179 }
180 else
181 {
182 std::cerr << "Didn't do bad atuators!\n";
183 }
184
187
188 head.append( "LINNACT", linNAct, "Linear number of actuators" );
189 head.append( "DIAMACT", diameter, "Diameter of DM in actuators" );
190 head.append( "DIAMPIX", diameter * act_space, "Diameter of DM in pixels" );
191 head.append( "ACTSPACE", act_space, "Actuator spacing, in pixels" );
192 head.append( "IFTYPE", "GAUSSIAN", "Type of influence function" );
193 head.append( "COUPLING", coupling, "Coupling at 1 actuator spacing" );
194 head.append( "FWHM", fwhm / act_space, "Resultant FWHM in actuators" );
195 head.append( "OFFSTODD", offsetOdd, "Whether odd rows are offset by 1/2" );
196
197 std::string fName = mx::AO::path::dm::influenceFunctions( dmName, true );
198 ff.write( fName, act_inf, head );
199
200 //--------------------------------------------
201 //-- Write coords file --
202 //--------------------------------------------
203 Eigen::Array<realT, -1, -1> coords( 2, xcoords.size() );
204 for( int i = 0; i < xcoords.size(); ++i )
205 {
206 coords( 0, i ) = xcoords[i];
207 coords( 1, i ) = ycoords[i];
208 }
209
210 fName = mx::AO::path::dm::actuatorPositions( dmName, true );
211 ff.write( fName, coords );
212 }
213};
214
215/// Create a set of Gaussian influence functions (IFs) on a square grid.
216/** The width of the IFs is specified by the coupling, which is the height of the IF at 1 actuator separation. For a
217 * FWHM = 1 actuator, coupling = 0.0625. The IFs are normalized to a height of 1.
218 *
219 * \param [in] dmName the name of the deformable mirror (DM) (the mx::AO path will be constructed from this).
220 * \param [in] dmSz the linear size of the DM in pixels.
221 * \param [in] linNAct the linear number of actuators across the DM.
222 * \param [in] diameter is the diameter of the DM, in actuators. If 0 then the DM is a square.
223 * \param [in] coupling is the height of the IF at 1 actuator separation. E.G., for a FWHM = 1 actuator, set this to
224 * 0.0625. \param [in] pupilSz is the size of the pupil, and therefore the size of the maps (pupilSz x pupilSz), Is set
225 * to dmSz if 0, can be pupilSz < dmSx for an oversized DM. The pupil is assumed to be centered.
226 *
227 * \tparam realT is the real floating point type in which to do calculations.
228 */
229template <typename realT>
231 const std::string &dmName,
232 realT dmSz,
233 int linNAct,
234 realT diameter,
235 realT coupling,
236 realT couplingRange, ///< Uniformly distributed uncertainy in coupling, in fraction of the coupling.
237 realT pupilSz = 0,
238 bool offsetOdd = false )
239{
240 if( pupilSz <= 0 )
241 pupilSz = dmSz;
242
243 realT pupilScale = pupilSz / dmSz;
244
245 realT act_space = ( dmSz - 1.0 ) / linNAct;
246
247 std::vector<realT> xcoords, ycoords;
248
249 realT xc, yc, xcen, ycen;
250 realT _rad;
251
252 // The center of the pupil
253 xcen = 0.5 * ( pupilSz - 1 ); // 0.5*(linNAct*pupilScale - 1);
254 ycen = 0.5 * ( pupilSz - 1 ); // 0.5*(linNAct*pupilScale - 1);
255
256 realT xoff = 0.0;
257
258 mx::math::uniDistT<realT> uniVar; ///< Uniform deviate, used in shiftRandom.
259
260 if( couplingRange > 0 )
261 {
262 uniVar.seed();
263 }
264
265 for( int i = 0; i < linNAct; ++i )
266 {
267 xc = ( i + 0.5 ) * act_space - 0.5 * ( dmSz - 1 ); // coordinate in pixels w.r.t. dm, not pupil
268
269 for( int j = 0; j < linNAct; ++j )
270 {
271 yc = ( j + 0.5 ) * act_space - 0.5 * ( dmSz - 1 ); // coordinate in pixels w.r.t dm, not pupil
272
273 if( offsetOdd )
274 {
275 if( j % 2 )
276 xoff = 0.0;
277 else
278 xoff = 0.5 * act_space;
279 }
280
281 _rad = sqrt( pow( xc + xoff, 2 ) + pow( yc, 2 ) );
282
283 if( _rad / act_space <= 0.5 * diameter + .01 || diameter == 0 )
284 {
285 xcoords.push_back( xcen + ( xc + xoff ) + 0.0 * ( dmSz - pupilSz ) );
286 ycoords.push_back( ycen + yc + 0.0 * ( dmSz - pupilSz ) );
287 }
288 }
289 }
290
291 std::cout << "Number of actuators: " << xcoords.size() << "\n";
292
293 /*gnuPlot gp;
294 gp.command("set size square");
295 gp.plot(xcoords, ycoords, "w p ps 0.75");
296 gp.point(xcen, ycen);
297 gp.circle(xcen, ycen, 0.5*(pupilSz-1), "w l", "Pupil");
298 */
299
300 realT fwhm = act_space * 2. * sqrt( log2( 2 ) ) / sqrt( -log2( coupling ) );
301
302 mx::improc::eigenCube<realT> act_inf( pupilSz, pupilSz, xcoords.size() );
303
304 realT fw;
305
306 for( int i = 0; i < xcoords.size(); ++i )
307 {
308 xc = xcoords[i];
309 yc = ycoords[i];
310
311 if( couplingRange > 0 )
312 {
313 fw = fwhm + fwhm * ( 1.0 - 2.0 * uniVar ) * couplingRange;
314 }
315 else
316 {
317 fw = fwhm;
318 }
319
320 mx::math::func::gaussian2D<realT>(
321 act_inf.image( i ).data(), pupilSz, pupilSz, 0.0, 1.0, xc, yc, math::func::fwhm2sigma( fw ) );
322 }
323
326
327 head.append( "LINNACT", linNAct, "Linear number of actuators" );
328 head.append( "DIAMACT", diameter, "Diameter of DM in actuators" );
329 head.append( "DIAMPIX", diameter * act_space, "Diameter of DM in pixels" );
330 head.append( "ACTSPACE", act_space, "Actuator spacing, in pixels" );
331 head.append( "IFTYPE", "GAUSSIAN", "Type of influence function" );
332 head.append( "COUPLING", coupling, "Coupling at 1 actuator spacing" );
333 head.append( "FWHM", fwhm / act_space, "Resultant FWHM in actuators" );
334 head.append( "OFFSTODD", offsetOdd, "Whether odd rows are offset by 1/2" );
335
336 std::string fName = mx::AO::path::dm::influenceFunctions( dmName, true );
337 ff.write( fName, act_inf, head );
338
339 //--------------------------------------------
340 //-- Write coords file --
341 //--------------------------------------------
342 Eigen::Array<realT, -1, -1> coords( 2, xcoords.size() );
343 for( int i = 0; i < xcoords.size(); ++i )
344 {
345 coords( 0, i ) = xcoords[i];
346 coords( 1, i ) = ycoords[i];
347 }
348
349 fName = mx::AO::path::dm::actuatorPositions( dmName, true );
350 ff.write( fName, coords );
351}
352
353/// Calculate the pseudo-inverse and mirror modes for a set of influence functions.
354/** The influence functions are a cube of 2D maps of mirror deformations, and the matrix \f$ A \f$ is formed by
355 * vectorizing the influence functions. The pseudo inverse is then calculated from the SVD as follows: \f$ A = U S V^T
356 * \f$ which gives the Moore-Penrose pseudo-inverse: \f$ A^+ = V S^+ U^T\f$. The columns of U contain the orthogonal
357 * mirror modes, and S contains the singular values. These are both written to disk along with the pseudo-inverse.
358 *
359 * \param[in] dmName is the name of the deformable mirror.
360 * \param[in] maxCondition [optional] the maximum condition number to accept for the pseudo-inverse. If < 0 [default]
361 * then the user is interactively asked what to use.
362 *
363 * \tparam realT is the real floating point type in which to do calculations.
364 */
365template <typename realT>
366void ifPInv( const std::string &dmName, realT maxCondition = -1 )
367{
368 std::string pinvName = mx::AO::path::dm::pseudoInverse( dmName );
369 std::string mmodesName = mx::AO::path::dm::mirrorModes( dmName, true );
370 std::string svsName = mx::AO::path::dm::singularValues( dmName );
371 std::string ifName = mx::AO::path::dm::influenceFunctions( dmName );
372
373 std::cout << "mmodesName: " << mmodesName << "\n";
374
376
378 ff.read( ifName, actInf );
379
380 int nrows = actInf.rows();
381 int ncols = actInf.cols();
382 int nplanes = actInf.planes();
383
384 Eigen::Array<realT, -1, -1> rowInf = actInf.asVectors();
385 actInf.resize( 0, 0, 0 );
386
387 Eigen::Array<realT, -1, -1> PInv;
388
389 realT condition;
390 int nRejected;
391
392 Eigen::Array<realT, -1, -1> U, S, VT;
393
394 int interact = MX_PINV_PLOT | MX_PINV_ASK;
395 if( maxCondition >= 0 )
396 interact = MX_PINV_NO_INTERACT;
397
398 mx::math::eigenPseudoInverse( PInv, condition, nRejected, U, S, VT, rowInf, maxCondition, interact );
399
401 head.append( "MAXCONDN", maxCondition, "Max. condition no. in pseudo inverse" );
402 head.append( "CONDN", condition, "Actual condition number of pseudo inverse" );
403 head.append( "NREJECT", nRejected, "Number of s.v.s rejected" );
404
405 // std::string pinvFile = pinvName;//ifName + "_pinv.fits";
406
407 ff.write( pinvName, PInv, head );
408
409 // This maps to a cube, but does not copy or take ownership.
410 mx::improc::eigenCube<realT> mmodes( U.data(), nrows, ncols, nplanes );
411
412 ff.write( mmodesName, mmodes, head );
413
414 std::ofstream fout;
415 fout.open( svsName );
416
417 for( int i = 0; i < S.rows(); ++i )
418 {
419 fout << S( i ) << "\n";
420 }
421
422 fout.close();
423}
424
425/// Calculate the modes-to-commands matrix for a set of modes
426/**
427 * Given the pseudo-inverse \f$ A^+ \f$ of the influence functions, the commands to the actuators to take a mirror
428 shape \f$ \vec{s} \f$ are calculated as
429 * \f[
430 \vec{c} = \mathbf{A^+} \vec{s}.
431 * \f]
432 * Now given a basis set \f$ \mathbf{M} \f$, the mirror shape is specified as
433 \f[
434 \vec{s}= \sum_i m_i M_i
435 \f]
436 * where \f$ \vec{m} = m_0, m_1, \cdot m_i \cdot \f$ are the modal coefficients. If the mirror-to-commands matrix, \f$
437 \mathbf{M2c} \f$ gives commands as
438 * \f[
439 \vec{c} = \mathbf{M2c} \vec{m}
440 * \f]
441 * then we can calculate \f$ \mathbf{M2c} \f$ as:
442 * \f[
443 \mathbf{M2c} = \mathbf{A}^{+T} \mathbf{M}
444 \f]
445 *
446 * \param[out] M2c is the M2c matrix, allocated during the calculationg
447 * \param[in] Ap is the pseudo-inverse of the influence function basis set
448 * \param[in] M is the modal basis set, a cube of shapes
449 *
450 * \tparam realT is the real floating point type used for all calculations.
451 */
452template <typename realT>
453void m2cMatrix( Eigen::Array<realT, -1, -1> &M2c, Eigen::Array<realT, -1, -1> &Ap, mx::improc::eigenCube<realT> &M )
454{
455 M2c = Ap.matrix().transpose() * M.asVectors().matrix();
456}
457
458/// Calculate the modes-to-commands matrix for a set of modes
459/**
460 * A wrapper for \ref m2cMatrix, where the various matrices are here specified with names, which
461 * in turn generate mx::AO paths
462 *
463 * \param[in] dmName is the name of the deformable mirror
464 * \param[in] basisName is the name of the modal basis
465 * \param[in] pupilName is the name of the pupil
466 *
467 * \tparam realT is the real floating point type in which to do calculations.
468 */
469template <typename realT>
470void m2cMatrix( const std::string &dmName, const std::string &basisName )
471{
472 std::string M2cFName;
474
475 Eigen::Array<realT, -1, -1> Ap; // The pseudo-inverse
476 std::string pinvFName = mx::AO::path::dm::pseudoInverse( dmName );
477 ff.read( pinvFName, Ap ); // read from the file
478
479 mx::improc::eigenCube<realT> M; // The modal basis, an image cube of shapes
480 std::string basisFName;
481
482 basisFName = mx::AO::path::basis::modes( basisName );
483 ff.read( basisFName, M ); // read from file
484
485 Eigen::Array<realT, -1, -1> M2c; // The M2c matrix
486 m2cMatrix( M2c, Ap, M );
487
488 M2cFName = mx::AO::path::dm::M2c( dmName, basisName, true );
489
490 ff.write( M2cFName, M2c );
491}
492
493/// Generate the modes-to-commands matrix for a set of modes on the modal DM
494/**
495 * Generates an Identity matrix of the appropriate size.
496 *
497 * \param[in] basisName is the name of the modal basis
498 *
499 * \tparam realT is the real floating point type in which to do calculations.
500 */
501template <typename realT>
502void modalDMM2cMatrix( const std::string &basisName )
503{
505
506 mx::improc::eigenCube<realT> M; // The modal basis, an image cube of shapes
507
508 std::string basisFName;
509
510 basisFName = mx::AO::path::basis::modes( basisName );
511
512 ff.read( basisFName, M ); // read from file
513
514 Eigen::Array<realT, -1, -1> M2c; // The M2c matrix
515
516 M2c.resize( M.planes(), M.planes() );
517 M2c.setZero();
518
519 for( int i = 0; i < M2c.rows(); ++i )
520 {
521 M2c( i, i ) = 1;
522 }
523
524 std::string M2cFName;
525
526 M2cFName = mx::AO::path::dm::M2c( "modalDM", basisName, true );
527
528 ff.write( M2cFName, M2c );
529}
530
531/// Calculate the basis set as it will be reproduced by the DM.
532/** The projected modes are written to the path specified by mx::AO::path::dm::projectedModes().
533 *
534 * \param [in] dmName is the name of the DM.
535 * \param [in] basisName is the name of the basis.
536 *
537 * \tparam realT is the real floating point type in which to do calculations.
538 */
539template <typename realT>
540void modesOnDM( const std::string &dmName, const std::string &basisName )
541{
542 Eigen::Array<realT, -1, -1> m2c; //, m, c;
543
544 std::string m2cName;
545
546 m2cName = mx::AO::path::dm::M2c( dmName, basisName );
547
549 ff.read( m2cName, m2c );
550
551 std::cout << m2c.rows() << "\n";
552 std::cout << m2c.cols() << "\n";
553
554 int nmodes = m2c.cols();
555 int nact = m2c.rows();
556
557 // m.resize(nmodes,1);
558
559 mx::improc::eigenCube<realT> act_inf, projModes;
560
561 std::string fName = mx::AO::path::dm::influenceFunctions( dmName, true );
562 ff.read( fName, act_inf );
563
564 projModes.resize( act_inf.rows(), act_inf.cols(), nmodes );
565 projModes.setZero();
566
567#pragma omp parallel for
568 for( int i = 0; i < nmodes; ++i )
569 {
570 Eigen::Array<realT, -1, -1> m, c;
571
572 std::cout << i + 1 << "/" << nmodes << "\n";
573 m.resize( nmodes, 1 );
574
575 m.setZero();
576 m( i, 0 ) = 1.0;
577 c = m2c.matrix() * m.matrix();
578
579 // #pragma omp parallel
580 //{
581 // Eigen::Array<realT, -1, -1> tmp(act_inf.rows(), act_inf.cols());
582 // tmp.setZero();
583
584 // #pragma omp for
585 for( int j = 0; j < nact; ++j )
586 {
587 projModes.image( i ) += c( j, 0 ) * act_inf.image( j );
588 // tmp += c(j,0) * act_inf.image(j);
589 }
590 // #pragma omp critical
591 // projModes.image(i) += tmp;
592 //}
593 }
594
595 std::string oName;
596 oName = mx::AO::path::dm::projectedModes( dmName, basisName, true );
597
598 ff.write( oName, projModes );
599}
600
601} // namespace AO
602} // namespace mx
603
604#endif //__influenceFunctions_hpp__
Standardized paths for the mx::AO system.
std::string actuatorPositions(const std::string &dmName, bool create=false)
The path for the deformable mirror (DM) actuator positions.
Definition aoPaths.hpp:118
std::string influenceFunctions(const std::string &dmName, bool create=false)
The path for the deformable mirror (DM) influence functions.
Definition aoPaths.hpp:104
Class to manage interactions with a FITS file.
Definition fitsFile.hpp:52
int write(const dataT *im, int d1, int d2, int d3, fitsHeader *head)
Write the contents of a raw array to the FITS file.
int read(dataT *data)
Read the contents of the FITS file into an array.
Definition fitsFile.hpp:825
Class to manage a complete fits header.
void append(fitsHeaderCard card)
Append a fitsHeaderCard to the end of the header.
An image cube with an Eigen-like API.
Definition eigenCube.hpp:30
Eigen::Map< Eigen::Array< dataT, Eigen::Dynamic, Eigen::Dynamic > > asVectors()
Return an Eigen::Eigen::Map of the cube where each image is a vector.
Eigen::Map< Eigen::Array< dataT, Eigen::Dynamic, Eigen::Dynamic > > image(Index n)
Returns a 2D Eigen::Eigen::Map pointed at the specified image.
A random number type, which functions like any other arithmetic type.
Definition randomT.hpp:64
void seed(typename ranengT::result_type seedval)
Set the seed of the random engine.
Definition randomT.hpp:101
constexpr units::realT k()
Boltzmann Constant.
Definition constants.hpp:69
int eigenPseudoInverse(Eigen::Array< dataT, -1, -1 > &PInv, dataT &condition, int &nRejected, Eigen::Array< dataT, -1, -1 > &U, Eigen::Array< dataT, -1, -1 > &S, Eigen::Array< dataT, -1, -1 > &VT, Eigen::Array< dataT, -1, -1 > &A, dataT &maxCondition, dataT alpha=0, int interact=MX_PINV_NO_INTERACT)
Calculate the pseudo-inverse of a patrix using the SVD.
realT gaussian2D(const realT x, const realT y, const realT G0, const realT G, const realT x0, const realT y0, const realT sigma)
Find value at position (x,y) of the 2D arbitrarily-centered symmetric Gaussian.
Definition gaussian.hpp:127
floatT fwhm2sigma(floatT fw)
Convert from FWHM to the Gaussian width parameter.
Definition gaussian.hpp:64
void modalDMM2cMatrix(const std::string &basisName)
Generate the modes-to-commands matrix for a set of modes on the modal DM.
void modesOnDM(const std::string &dmName, const std::string &basisName)
Calculate the basis set as it will be reproduced by the DM.
void m2cMatrix(Eigen::Array< realT, -1, -1 > &M2c, Eigen::Array< realT, -1, -1 > &Ap, mx::improc::eigenCube< realT > &M)
Calculate the modes-to-commands matrix for a set of modes.
void ifPInv(const std::string &dmName, realT maxCondition=-1)
Calculate the pseudo-inverse and mirror modes for a set of influence functions.
void influenceFunctionsGaussian(const std::string &dmName, realT dmSz, int linNAct, realT diameter, realT coupling, realT couplingRange, realT pupilSz=0, bool offsetOdd=false)
Create a set of Gaussian influence functions (IFs) on a square grid.
The mxlib c++ namespace.
Definition mxError.hpp:106