174 Acta Polytechnica CTU Proceedings 2(1): 174–177, 2015 174 doi: 10.14311/APP.2015.02.0174 Application of Total Variation Minimization to Doppler Tomography M. Uemura1, T. Kato2, D. Nogami3, R. Mennickent4 1Hiroshima Astrophysical Science Center, Hiroshima University, Kagamiyama 1-3-1, Higashi-Hiroshima, Hiroshima 739-8526, Japan 2Department of Astronomy, Kyoto University, Kitashirakawa-Oiwake-cho, Sakyo-ku, Kyoto 606-8502, Japan 3Kwasan and Hida Observatories, Kyoto University, Yamashina-ku, Kyoto 607-8471, Japan 4Universidad de Concepción, Departamento de Astronomı́a, Casilla 160-C, Concepción, Chile Corresponding author: uemuram@hiroshima-u.ac.jp Abstract We have developed a new model of the Doppler tomography using total variation minimization (DTTVM). We demon- strated that this method can reconstruct localized and non-axisymmetric profiles possibly having sharp edges in the Doppler map. We apply this model to the real data of the dwarf nova, WZ Sge in superoutburst. DTTVM can reproduce the observed spectra with a high precision, while the previous models fail to reproduce localized sources. Keywords: cataclysmic variables - dwarf novae - optical - spectroscopy - individual: WZ Sge. 1 Introduction Doppler tomography (DT) is a widely-used method in the field of cataclysmic variables (CVs) [1]. DT is an ill-posed inverse problem in most cases. A filtered back- projection method provides a way to solve the problem. It masks high frequency components in the Doppler maps in order to uniquely determine the solution. The other way is to introduce regularization terms in ad- dition to the likelihood term. The maximum entropy method (MEM) has been used to regularize the DT problem. The DT with MEM determines the solution by maximizing the information entropy over the default image. The comparison between the back-projection method and MEM is reported in [2]. Those DT methods tend to smeared out the local- ized and sharp-edge profiles in the Doppler map. On the other hand, such profiles are expected to appear in accretion disks of CVs, such as spiral shocks and hot spot. It has been reported that there are significant residuals between the observed and model spectra that are obtained from DT (e.g. [3]). It is, however, un- clear if the residuals indicate real emission components that are not reproduced by the model. In this paper, we introduce a new DT by total variation minimization (TVM), which is described in [4]. TVM has received attention in the field of image reconstruction because it can reconstruct sharp-edge features by making image sparse in its gradient domain. 2 Doppler Tomography by Total Variation Minimization 2.1 Model Suppose there are M data given by d = {d0,d1, · · · ,dM}. The i-th element of d, di, is the flux at a radial velocity vR,i and phase φi. Also sup- pose n × n = N points in the Doppler map given by s = {s0,0,s0,1, · · · ,sn,n}. The j-th element in the map has the coordinate of (vx,j,vy,j) and intensity of sj. The goal of DT is to estimate ŝ from the observation d according to ŝ = arg min s { 1 M ‖d−As‖22 + λ N Φ(x) } . (1) The first term in the right side is a least-square term, and the second term is a regularization term which con- sists of a hyperparameter λ and a function Φ. A is a M × N matrix whose element depends on the phase, radial velocity, and instrumental response. For the re- sponse function, we assumed a Gaussian profile with a variance σ2, where σ corresponds to the velocity reso- lution of the observed spectra. An element of A, aij is then given by aij = 1 √ 2πσ2 exp { − (vR,i −vR,ij)2 2σ2 } . (2) TVM provides a method to solve ill-posed problems. 174 http://dx.doi.org/10.14311/APP.2015.02.0174 Application of Total Variation Minimization to Doppler Tomography The isotropic total variation ΦTVM,iso is defined as ΦTVM,iso(s) = ∑ i,j √ (si+1,j −si,j)2 + (si,j+1 −si,j)2. (3) ΦTVM,iso is used as Φ in equation (1) in our model. We use the TwIST algorithm to estimate ŝ in equa- tion (1) [5, 6].1 The calculation code of DTTVM can be available at our web site.2 2.2 Demonstrations using artificial data We show an example of the results of DTTVM using artificial data sets in figure 1. First, we made the as- sumed Doppler maps, then simulated observations of line-profile variations, and reconstructed the maps us- ing DTTVM. We calculated three Doppler maps from the simulated spectra with different values of λ in equa- tion (1). We assumed the original Doppler map con- taining three spot structures, as shown in the upper leftmost panel of figure 1. The upper two spots have sharp-edge and flat-top profiles of different size and in- tensity. The lower left spot has a Gaussian profile. DTTVM reproduced the original map in all three cases. In case (a), the edges of the flat-top spots are slightly smeared out compared to cases (b) and (c). In addi- tion, the peak intensity of the Gaussian spot is under- estimated in case (a). These smearing effects are less prominent for smaller values of λ. The model spectra and residuals are also shown in figure 1. The resid- uals are smaller for smaller λ. These results demon- strate that DTTVM can reconstruct both sharp-edge and smooth profiles in a Doppler map. 2.3 Application to the data of WZ Sge We apply DTTVM to the data of the dwarf nova, WZ Sge taken on the 10th day of the superoutburst in 2001 [7]. Figure 2 presents the Doppler maps and model spectra. The results were computed using log10 λ of (a) −2.2, (b) −3.2, and (c) −5.2. All Doppler maps exhibit a disk feature having non-axisymmetric distribution of intensity. The peak intensity of the localized bright re- gions is lower in case (a) than the others. As a result, the residuals are large in case (a), compared with the other cases. The model spectra reproduce the observed one without any prominent residuals in cases (b) and (c). Figure ?? shows the result of the Doppler tomog- raphy using the same data of WZ Sge calculated with MEM. We used the calculation code presented in [8].3 The calculation code defines α as a hyperparameter of the model. Cases (a), (b), and (c) in figure ?? were cal- culated with log α = −1.0, −2.0, and −3.0, respectively. In cases (a) and (b), smooth two-armed structures ap- pear on a weak sign of a disk feature. In case (c), patchy structures appear, but they are evidently artifacts. In all cases of (a), (b), and (c), there are significant resid- uals of the spectral data at −500 km s−1 for phase 1.75 and at +500 km s−1 for phase 1.60. Both features have counterparts at phases shifted by 0.5, which indicate they are rotating components with the orbital period. These rotating residuals definitely arise from localized structures in the Doppler map, which can be success- fully reconstructed by DTTVM: The former one cor- responds to the feature at (vx,vy) ∼ (0, 500) and the latter one to the lower right region in figure 2. 3 Summary As demonstrated above, DTTVM can reconstruct both sharp-edge and smooth profiles in a Doppler map. Hence, this model is advantageous for the reconstruc- tion of highly localized or edge structures in the Doppler map, such as the emission from the secondary star, hot spot, or shock wavefront. The Doppler map of WZ Sge obtained by DTTVM is surprisingly different from that built with the MEM method. The difference in the regularization term could cause the different maps es- pecially when the data sample is small. DTTVM can reproduce the observed spectra with a high accuracy, even in the case that the MEM method generate sig- nificant residuals in spectra. This indicates that TVM is probably a suitable regularization method for real Doppler map, compared with MEM. Acknowledgement This work was supported by JSPS KAKENHI Grant Number 22540252 and 25120007. We appreciate com- ments and suggestions from Dr. Shiro Ikeda. REM acknowledges support by Fondecyt grant 1110347 and the BASAL Centro de Astrofisica y Tecnologias Afines (CATA) PFB–06/2007. 1〈http://www.lx.it.pt/˜bioucas/TwIST/TwIST.htm〉 2〈http://home.hiroshima-u.ac.jp/uemuram/dttvm/〉 3〈http://www.mpa-garching.mpg.de/˜henk/〉 175 M. Uemura et al. References [1] Marsh, T. R. & Horne, K. 1988, MNRAS, 235, 269 doi:10.1093/mnras/235.1.269 [2] Marsh, T. R. 2001, in Astrotomography, Indirect Imaging Methods in Observational Astronomy, ed. H. M. J. Boffin, D. Steeghs, & J. Cuypers (Springer-Verlag), 1 [3] Tappert, C., Mennickent, R. E., Arenas, J., Matsumoto, K., & Hanuschik, R. W. 2003, A&A, 408, 651 [4] Uemura, M., Kato, T., Nogami, D., & Mennickent, R. E. 2013, PASJ, submitted [5] Bioucas-Dias, J. M. & Figueiredo, M. A. 2007b, in Image Processing, 2007. ICIP 2007. IEEE International Conference on Vol. 1 (IEEE), pp I–105 [6] Bioucas-Dias, J. M. & Figueiredo, M. A. 2007a, Image Processing, IEEE Transactions on, 16, 2992 [7] Nogami, D. & Iijima, T. 2004, PASJ, 56, 163 [8] Spruit, H. C. 1998, arXiv:astro-ph/9806141 DISCUSSION DMITRY KONONOV: 1. It’s nice that your method can restore flat features. 2. It’s in a sense incorrect to compare two results of inverse-problem solving methods if you have no physical basics under the structures you are trying to explain. MAKOTO UEMURA: As you pointed out, both TVM and MEM have no physical meaning. But, I think the comparison still makes a sense. We can select a better model based on residuals between observations and model. In the case of WZ Sge, the residuals ob- tained with our method are much smaller than those with MEM. The MEM method fails to reproduce the localized structures of the line profile, while our method can reproduce them. DMITRY BISIKALO: Is it possible to use TVM in 3D Doppler tomography? MAKOTO UEMURA: Yes, it is. The 3D Doppler tomography is also a linear problem, so it is easy to apply the method to it. Vx (km s−1) V y (k m s −1 ) −1000 0 +1000 − 1 0 0 0 0 + 1 0 0 0 Original − 1 5 1 0 −1000 0 +1000 − 1 0 0 0 0 + 1 0 0 0 (a)logλ=−1.0 −1000 0 +1000 (b) logλ=−2.0 −1000 0 +1000 (c) logλ=−3.0 RV (km s−1) O rb ita l p h a se −2000 0 +2000 0 .0 0 .2 0 .4 0 .6 0 .8 1 .0 −2000 0 +2000 0 .0 0 .2 0 .4 0 .6 0 .8 1 .0 0 +2000 0 +2000 0 +2000 0 +2000 0 +2000 Figure 1: Results of DTTVM using an artificial map having three spots. The artificial map is depicted in the upper left-most panel. The simulated spectra from the map is shown in the lower left-most panel. Three sets of results are shown with different λ, which is indicated by the top of the figure. For each set, the estimated Doppler map is shown in the upper panel. The lower left and right panels are the model spectra and residuals between the original and model spectra, respectively. 176 http://dx.doi.org/10.1093/mnras/235.1.269 Application of Total Variation Minimization to Doppler Tomography Vx (km s−1) V y (k m s −1 ) −1000 0 +1000 − 1 0 0 0 0 + 1 0 0 0 (a) logλ=−2.2 Vx (km s−1) V y (k m s −1 ) −1000 0 +1000 − 1 0 0 0 0 + 1 0 0 0 (b) logλ=−3.2 Vx (km s−1) V y (k m s −1 ) −1000 0 +1000 − 1 0 0 0 0 + 1 0 0 0 (c) logλ=−5.2 RV (km s−1) O rb ita l p h a se −2000 0 +2000 1 .0 1 .5 2 .0 2 .5 RV (km s−1) −2000 0 +2000 1 .0 1 .5 2 .0 2 .5 RV (km s−1) −2000 0 +2000 1 .0 1 .5 2 .0 2 .5 RV (km s−1) −2000 0 +2000 1 .0 1 .5 2 .0 2 .5 RV (km s−1) −2000 0 +2000 1 .0 1 .5 2 .0 2 .5 RV (km s−1) −2000 0 +2000 1 .0 1 .5 2 .0 2 .5 Figure 2: Results of DTTVM for the data of WZ Sge. Upper panels: the Doppler maps. Lower left and right panels: the model spectra and residuals between the model and observed spectra, respectively. The three sets of the panels were calculated with different λ, which are indicated in the figure. Vx (km s−1) V y (k m s −1 ) −1000 0 +1000 − 1 0 0 0 0 + 1 0 0 0 (a) logα=−1.0 Vx (km s−1) V y (k m s −1 ) −1000 0 +1000 − 1 0 0 0 0 + 1 0 0 0 (b) logα=−2.0 Vx (km s−1) V y (k m s −1 ) −1000 0 +1000 − 1 0 0 0 0 + 1 0 0 0 (c) logα=−3.0 RV (km s−1) O rb ita l p h a se −2000 0 +2000 1 .0 1 .5 2 .0 2 .5 RV (km s−1) −2000 0 +2000 1 .0 1 .5 2 .0 2 .5 RV (km s−1) −2000 0 +2000 1 .0 1 .5 2 .0 2 .5 RV (km s−1) −2000 0 +2000 1 .0 1 .5 2 .0 2 .5 RV (km s−1) −2000 0 +2000 1 .0 1 .5 2 .0 2 .5 RV (km s−1) −2000 0 +2000 1 .0 1 .5 2 .0 2 .5 Figure 3: As in figure 2, but for the MEM case. 177 Introduction Doppler Tomography by Total Variation Minimization Model Demonstrations using artificial data Application to the data of WZ Sge Summary