Phase $ \ phi $ is similar to Delay
\hat{\phi}_k = \phi_j - \phi_i \tag{4.1}
However, there are two points to note.
Regarding 2., we will deal with it by selecting a point source as the phase calibration object. Regarding 1., in the case where the phase of each baseline is within $ \ pm \ pi / 2 $, it can be treated in the same way as the delay solution. Otherwise, it is generally better to solve with a complex solution.
The phase term of the antenna gain is
e^{i\hat{\phi}_k} = e^{i\phi_j} e^{-i\phi_i} \tag{4.2}
It is expressed as. It satisfies equation (4.1) and has no indefiniteness of $ 2n \ pi $. However, since it is a non-linear equation, iteration by successive approximation is required. In the following, $ E_j = e ^ {i \ phi_j} $ is used to simplify the notation. Suppose the initial value $ E ^ 0_i $ is close enough to the optimal solution. The baseline-based residual is $ r_k = E_k --E ^ 0_j E_i ^ {0 *} $. The correction vector $ \ boldsymbol {c} $ given to $ \ phi ^ 0_j $ to minimize the residual norm $ \ sum_k r_k r ^ * _k $ is
\left(
\begin{array}{c}
r_0 \\
r_1 \\
r_2 \\
r_3 \\
r_4 \\
r_5 \\
\vdots
\end{array}
\right) = \left( \begin{array}{cccc}
\frac{\partial {E_1 E^*_0}}{\partial \phi_1} & \frac{\partial {E_1 E^*_0}}{\partial \phi_2} & \frac{\partial {E_1 E^*_0}}{\partial \phi_3} & \cdots \\
\frac{\partial {E_2 E^*_0}}{\partial \phi_1} & \frac{\partial {E_2 E^*_0}}{\partial \phi_2} & \frac{\partial {E_2 E^*_0}}{\partial \phi_3} & \cdots \\
\frac{\partial {E_2 E^*_1}}{\partial \phi_1} & \frac{\partial {E_2 E^*_1}}{\partial \phi_2} & \frac{\partial {E_2 E^*_1}}{\partial \phi_3} & \cdots \\
\frac{\partial {E_3 E^*_0}}{\partial \phi_1} & \frac{\partial {E_3 E^*_0}}{\partial \phi_2} & \frac{\partial {E_3 E^*_0}}{\partial \phi_3} & \cdots \\
\frac{\partial {E_3 E^*_1}}{\partial \phi_1} & \frac{\partial {E_3 E^*_1}}{\partial \phi_2} & \frac{\partial {E_3 E^*_1}}{\partial \phi_3} & \cdots \\
\frac{\partial {E_3 E^*_2}}{\partial \phi_1} & \frac{\partial {E_3 E^*_2}}{\partial \phi_2} & \frac{\partial {E_3 E^*_2}}{\partial \phi_3} & \cdots \\
\vdots & \vdots & \vdots & \ddots
\end{array} \right)
\left(
\begin{array}{c}
c_1 \\
c_2 \\
c_3 \\
\vdots
\end{array}
\right) \tag{4.3}
It is shown by the equation. Note that the phase of refant is set to 0 and generality is not lost, so $ c_0 $ is always 0, not unknown. $ \ frac {\ partial} {\ partial \ phi_0} $ is also 0. When I write out the elements of the matrix $ P $ in the middle,
P = \left( \begin{array}{cccc}
i{E_1 E^*_0} & 0 & 0 & \cdots \\
0 & i{E_2 E^*_0} & 0 & \cdots \\
-i{E_2 E^*_0} & i{E_2 E^*_0} & 0\\
0 & 0 & i{E_3 E^*_0} & \cdots \\
-i{E_3 E^*_1} & 0 & i{E_3 E^*_1} & \cdots \\
0 & -i{E_3 E^*_2} & i{E_3 E^*_2} & \cdots \\
\vdots & \vdots & \vdots & \ddots
\end{array} \right) \tag{4.4}
It looks like. The gain of refant is $ E_0 = E ^ * _0 = 1 $. Equation (4.3) is
\boldsymbol{r} = P \boldsymbol{c}
It is given by the simultaneous equations, so if you solve this for $ \ boldsymbol {c} $
P^{*T}\boldsymbol{r} = P^{*T} P \boldsymbol{c} \\
\boldsymbol{c} = (P^{*T} P)^{-1} P^{*T}\boldsymbol{r}
It will be. Add the $ \ boldsymbol {c} $ obtained in this way to the initial phase value $ \ boldsymbol {\ phi} ^ 0 $, and add $ \ boldsymbol {\ phi} ^ 1 \ leftarrow \ boldsymbol {\ phi} ^ Improve the phase as 0 + \ boldsymbol {c} $. If you repeat this, it will eventually converge. The convergence test can be judged by the norm of the correction vector $ \ boldsymbol {c} $.
Calculate $ P ^ {* T} P $ by hand in advance.
P^{*T} P = \left( \begin{array}{ccccc}
N_a - 1 & -1 & -1 & -1 & \cdots \\
-1 & N_a - 1 & -1 & -1 & \cdots \\
-1 & -1 & N_a-1 & -1 & \cdots \\
-1 & -1 & -1 & N_a-1 & \cdots \\
\vdots & \vdots & \vdots & \vdots &\ddots
\end{array} \right)
And the same matrix as for the delay. Therefore, the inverse matrix
(P^TP)^{-1}= \frac{1}{N_a}\left( \begin{array}{ccccc}
2 & 1 & 1 & 1 & \cdots \\
1 & 2 & 1 & 1 & \cdots \\
1 & 1 & 2 & 1 & \cdots \\
1 & 1 & 1 & 2 & \cdots \\
\vdots & \vdots & \vdots & \vdots &\ddots
\end{array} \right)
Therefore, you can use Same matrix as when getting the solution of delay. The calculation of $ P ^ {* T} \ boldsymbol {r} $ is a bit tricky, but you can get it with the following Python code. It is unavoidable to end up in a loop with the number of antennas, but a loop with the number of baselines can be avoided.
PTY = np.zeros(antNum - 1, dtype=complex)
for ant_index in range(1,antNum):
index0, index1 = np.where(ant0 == ant_index)[0].tolist(), np.where(ant1 == ant_index)[0].tolist()
Y = np.zeros(blNum, dtype=complex)
Y[index0] = -1.0j* antGain[ant_index].conjugate()* antGain[ant1[index0]]
Y[index1] = 1.0j* antGain[ant_index]* antGain[ant0[index1]].conjugate()
PTY[ant_index - 1] += Y.dot(resid)
#
Summarizing the above, we show the function clphase_solve () that solves the antenna-based phase from the baseline-based visibility. The argument Vis is visibility, which is a vector of complex baseline numbers in Canonial order. IterNum is the number of iterations, which is set to 2 by default (if the S / N ratio is good, it will converge sufficiently). The return value is the antenna-based phase.
def clphase_solve(Vis, iterNum = 2):
# Initial setup
Vis = Vis / abs(Vis) # Normalization (discard amplitude)
blNum = len(Vis); antNum = Bl2Ant(blNum)[0]
ant0, ant1 = np.array(ANT0[0:blNum]), np.array(ANT1[0:blNum])
# Partial matrix
PTP_inv = (np.diag(np.ones(antNum - 1)) + 1.0) / antNum
antGain = np.append(1.0 + 0.0j, Vis[KERNEL_BL[0:antNum-1]])
#---- Loop for iteration
for iter_index in range(iterNum):
resid = Vis - (antGain[ant0] * antGain[ant1].conjugate())
PTY = np.zeros(antNum - 1, dtype=complex)
#---- Calculate P*T r
for ant_index in range(1,antNum):
index0 = np.where(ant0 == ant_index)[0].tolist()
index1 = np.where(ant1 == ant_index)[0].tolist()
Y = np.zeros(blNum, dtype=complex)
Y[index0] = -1.0j* antGain[ant_index].conjugate()* antGain[ant1[index0]]
Y[index1] = 1.0j* antGain[ant_index]* antGain[ant0[index1]].conjugate()
PTY[ant_index - 1] += Y.dot(resid)
#
#---- Correction for the antenna-based phase
antPhs = np.append(0, np.angle(antGain[1:antNum]) + PTP_inv.dot(PTY.real))
antGain = np.cos(antPhs) + 1.0j* np.sin(antPhs)
#
return antPhs
#
that's all. Return to the lottery
Recommended Posts