Читать книгу Muography - Группа авторов - Страница 38
2.2 LINEAR INVERSION
ОглавлениеIn this section, we introduce the linear inversion method used by Nishiyama et al. (2014a). The method is based on the theory described in Tarantola (2005). By assuming that the observed number of muons and density probability distribution is Gaussian, the problem can be solved with a least‐squares method by introducing a priori information about the density structure as elements of the covariance matrix.
If the volume of interest is subdivided into voxels with index j and density ρ j (j = 1, 2, …), then the relationship between density length d i and ρ j can be represented as:
where i is the index of the muon pathway for any direction, A ij is the length across the j th voxel, and d i is the density length derived from the muon attenuation along the i th muon pathway. d i can represent the number of observed muons rather than the density length, but in this paper, we first treat d i as the density length. Equation 2.1 can be represented by vectors ,, and the matrix A:
In many cases, multi‐directional muography based on equation 2.2 is an under‐determined system that limits the spatial and density resolution. If is the observational value of , then and form a covariance matrix of and , respectively. In addition, we also define as initial values, set the model expectation as , and denote the posterior covariance matrix as . The relationship between these parameters is (e.g., Nishiyama et al., 2014a; Tarantola, 2005):
We now consider the meaning of the covariance matrices used in equations 2.3 and 2.4 in the field of probability theory. If is a set of probability variables, then the covariance matrix is defined as follows:
where V(x i ) is the deviation of x i and Cov(x i , x j ) = E[(x i − E[x i ])(x j − E[x j ])] (E[x i ] is the expected value of the probability variable x i ) and qualitatively reflects the correlation between x i and x j .
When solving for using equations 2.3 and 2.4, the systemis typically under‐determined, and thus every correlation Cov(ρ i , ρ j ) cannot be fixed. As such, we introduce another assumption into the density model. For the simplicity and effectivity for an under‐determined system, we may constrain the relationship between the i th and j th voxels as follows:
where σ ρ is the density contrast deviation, l(i, j) is the distance between the i th and j th voxel, and L 0 is the correlation length. Equation 2.6 assumes that the internal density is continuous on a spatial scale L 0, and that the density contrast is typically within σ ρ . However, equation 2.6 is just one possible example, and it is possible to assume different constraints in the model based on the expected structure. In this case, ρ 0, σ ρ , and L 0 are a priori parameters.
We now consider the matrix elements of in the simplest case. This matrix is also represented like equation 2.5, but the covariance between d i and d j (i ≠ j) is zero, if we assume there is no contamination of low momentums due to deflection in the mountain, and no cross talk between the adjacent angular bins in all detectors for simplicity. The elements can be written as:
(2.7)
where is the accidental error of di.
In equations 2.1 and 2.2, we described d i and as the density length; this can be replaced by the number of muons N i such that:
The elements of matrix A ij are different between equations 2.1 and 2.8. In equation 2.1, the elements of Aij can be calculated from the topology, size, and shape of the voxels that are defined. To calculate the elements of matrix , we need to consider the relationship between the density length X i and number of observed muons N i . This can be written as X i = f i (N i ) or N i = g i (X i ) by using the function f i and its inverse function g i . The functions f i and g i are not linear. We then linearize N i around with small variation δX i :
(2.9)
Given that N i can be written as and , δ N i becomes:
Introducing a small variation from the initial density value , δN i and δX i can be given as:
From equations 2.10 and 2.12 we obtain:
Finally, by comparing equations 2.11 and 2.13, the relationship between and is derived as:
(2.14)
There are some advantages to using equation 2.8 instead of equation 2.1. For example, we can define the relationship between the i th and k th bins (Fig. 2.1). If N' k is the number of muons after merging angular bins, then, by comparing N′ k = ∑ i N i , , and using equation 2.8, . Hence, when we merge the bins in angular space, we sum in the same way as the muon numbers. In the case that the mountain ridge line crosses the k th bin, we sum the smaller bins without bins above the ridge line. Another advantage is the increase in the calculation accuracy. The muon flux varies with pathway and elevation angle, and the accuracy of the calculation decreases as the bin size increases. In Jourde et al. (2015), this type of calculation was described using a form of integration.
Figure 2.1 Schematic of merging angular bins in generating a muographic image. Horizontal and vertical axes represent relative azimuth and elevation angle. Minimum angular bins constructed from the dotted gray lines (index i) belong to larger bins defined by the solid black squares (index k). The dashed black line shows the mountain ridge line.
The reconstructed density image obtained with equation 2.3 depends on the values of the parameters σ ρ and L 0 in equation 2.6. One solution to find the best‐fit parameters was described in Nishiyama et al. (2014a), which is to search for the parameters in a forward simulation by minimizing the discrepancy Φ between the assumed input density model and evaluated :
One clear advantage of this linear inversion method compared to the method described in the next section is the higher spatial resolution that can be achieved when combined with gravity observation data. However, the linear inversion result depends on the constraint function and parameters in equation 2.6.
