49 from typing
import List
54 modalIMat: np.ndarray,
55 modeSelect: np.ndarray =
None,
56 modeGains: np.ndarray =
None,
58 """ Generic numpy modal interaction matrix inversion function
62 M2V: (nActu x nModes) : modal basis matrix
64 modalIMat: (nSlopes x nModes) : modal interaction matrix
66 modeSelect: (nModes, dtype=bool): (Optional):
67 mode selection, mode at False is filtered
69 modeGains: (nModes, dtype=bool): (Optional):
70 modal gains to apply. These are gain in the reconstruction sens, ie
71 they are applied multiplicatively on the command matrix
73 if modeSelect
is None:
74 modeSelect = np.ones(modalIMat.shape[1], dtype=bool)
76 modeGains = np.ones(modalIMat.shape[1], dtype=np.float32)
78 return M2V.dot(modeGains[:,
None] * np.linalg.inv(modalIMat[:, modeSelect].T.dot(
79 modalIMat[:, modeSelect])).dot(modalIMat[:, modeSelect].T))
82 def cmat_init(ncontrol: int, rtc: Rtc, p_controller: conf.Param_controller,
83 p_wfss: List[conf.Param_wfs], p_atmos: conf.Param_atmos,
84 p_tel: conf.Param_tel, p_dms: List[conf.Param_dm],
85 nmodes: int = 0) ->
None:
86 """ Compute the command matrix on the GPU
94 p_controller: (Param_controller) : controller settings
96 p_wfss: (list of Param_wfs) : wfs settings
98 p_atmos: (Param_atmos) : atmos settings
100 p_tel : (Param_tel) : telescope settings
102 p_dms: (list of Param_dm) : dms settings
104 M2V : (np.ndarray[ndim=2, dtype=np.float32]): (optional) KL to volts matrix (for KL cmat)
106 nmodes: (int) : (optional) number of kl modes
108 if (p_controller.type == scons.ControllerType.LS):
109 print(
"Doing imat svd...")
111 rtc.d_control[ncontrol].svdec_imat()
112 print(
"svd done in %f s" % (time.time() - t0))
113 eigenv = np.array(rtc.d_control[ncontrol].d_eigenvals)
114 imat = np.array(rtc.d_control[ncontrol].d_imat)
115 maxcond = p_controller.maxcond
116 if (eigenv[0] < eigenv[eigenv.shape[0] - 1]):
117 mfilt = np.where((eigenv / eigenv[eigenv.shape[0] - 3]) < 1. / maxcond)[0]
119 mfilt = np.where((1. / (eigenv / eigenv[2])) > maxcond)[0]
120 nfilt = mfilt.shape[0]
122 print(
"Building cmat...")
124 if not p_controller.do_kl_imat:
125 print(
"Filtering ", nfilt,
" modes")
126 rtc.d_control[ncontrol].build_cmat(nfilt)
131 Dp_filt = np.linalg.inv(D_filt.T.dot(D_filt)).dot(D_filt.T)
132 if (p_controller.klgain
is not None):
133 Dp_filt *= p_controller.klgain[
None, :]
134 cmat_filt = p_controller._M2V.dot(Dp_filt)
135 rtc.d_control[ncontrol].set_cmat(cmat_filt)
137 print(
"cmat done in %f s" % (time.time() - t0))
139 if (p_controller.type == scons.ControllerType.MV):
140 Cn = np.zeros(p_controller._imat.shape[0], dtype=np.float32)
142 for k
in p_controller.nwfs:
143 Cn[ind:ind + 2 * p_wfss[k]._nvalid] = noise_cov(k, p_wfss[k], p_atmos, p_tel)
144 ind += 2 * p_wfss[k]._nvalid
146 rtc.d_control[ncontrol].load_noisemat(Cn)
147 print(
"Building cmat...")
148 rtc.d_control[ncontrol].build_cmat(p_controller.maxcond)
150 if (p_controller.TTcond ==
None):
151 p_controller.set_TTcond(p_controller.maxcond)
153 if (
"tt" in [dm.type
for dm
in p_dms]):
154 rtc.d_control[ncontrol].filter_cmat(p_controller.TTcond)
156 p_controller.set_cmat(np.array(rtc.d_control[ncontrol].d_cmat))
161 return np.linalg.svd(DtD)
165 """ Compute a command matrix in Btt modal basis (see error breakdown) and set
166 it on the sutra_rtc. It computes by itself the volts to Btt matrix.
170 rtc: (Rtc) : rtc object
172 dms: (Dms): dms object
174 p_dms: (list of Param_dm): dms settings
176 p_geom: (Param_geom): geometry settings
180 IFs = basis.compute_IFsparse(dms, p_dms, p_geom).T
182 IFtt = IFs[:, -2:].toarray()
183 IFpzt = IFs[:, :n - 2]
185 Btt, P = basis.compute_btt(IFpzt, IFtt)
189 def get_cmat(D, nfilt, Btt=None, rtc=None, svd=None):
190 """Compute a command matrix from an interaction matrix 'D'
194 get_cmat(D,nfilt,Btt=BTT,rtc=RTC)
195 get_cmat(D,nfilt,svd=SVD)
198 D: (np.ndarray[ndim=2, dtype=np.float32]): interaction matrix
200 nfilt: (int): number of element to filter
202 Btt: (np.ndarray[ndim=2, dtype=np.float32]): Btt modal basis
206 svd: (tuple of np.ndarray[ndim=1, dtype=np.float32): svd of D.T*D (obtained from np.linalg.svd)
208 nfilt = max(nfilt, 0)
209 if (Btt
is not None):
210 if (svd
is not None):
211 raise ValueError(
"Btt and SVD cannt be used together")
213 raise ValueError(
"Btt cannot be used without rtc")
215 index = np.concatenate((np.arange(n - nfilt - 2), np.array([n - 2,
216 n - 1]))).astype(int)
217 Btt_filt = Btt[:, index]
221 Dmp = np.linalg.inv(Dm.T.dot(Dm)).dot(Dm.T)
223 cmat = Btt_filt.dot(Dmp)
225 if (svd
is not None):
234 DtDx = v.T.dot(np.diag(s_filt)).dot(u.T)
237 return cmat.astype(np.float32)