## 2 Materials and Methods
### Normal Vector Estimation and the Coordinates of Interpolation Node  
Let \(\mathbf{Q}=\{\mathbf{X}_{t}=(x_{t},y_{t},z_{t})\,|\,t=1,\,\ldots,\,M\}\) represent the neighbor of \(\mathbf{p}_{1}\), in which \(M\) is the number of the neighbor point. The normal vector of the unknown isosurface at \(\mathbf{p}_{1}\) is associated with a tangent plane that is through the center of the point set \(\mathbf{Q}\). Let \(\mathbf{n}\) represent the unit normal vector and \(\mathbf{X}_{c}\) represent the center of \(\mathbf{Q}\). Then, \(\mathbf{X}_{c}\) can be computed from \(\mathbf{Q}\):  
\[\mathbf{X}_{c}=\frac{1}{M}{\sum_{t=1}^{M}\mathbf{X}_{t}}=(x_{c},y_{c},z_{c}), \tag{5}\]  
The unit normal vector \(\mathbf{n}\) can be computed by preceding a general PCA algorithm [4; 25]. The symmetric 3 \(\times\) 3 positive semi-definite covariance matrix \(\mathbf{CV}\) is formed using point set \(\mathbf{Q}\).  
\[\mathbf{CV}=\left(\begin{array}{c}\mathbf{X}_{1}-\mathbf{X}_{c}\\ \ldots\ldots\\ \mathbf{X}_{M}-\mathbf{X}_{c}\end{array}\right)^{T}\otimes\left(\begin{array} []{c}\mathbf{X}_{1}-\mathbf{X}_{c}\\ \ldots\ldots\\ \mathbf{X}_{M}-\mathbf{X}_{c}\end{array}\right), \tag{6}\]  
By a singular value decomposition (SVD) of covariance matrix \(\mathbf{CV}\), we get eigenvalues \(\lambda_{1}\geq\lambda_{2}\geq\lambda_{3}\), which are associated with eigenvectors \(\mathbf{v}_{1}\), \(\mathbf{v}_{2}\), \(\mathbf{v}_{3}\). According to the PCA algorithm, the eigenvector \(\mathbf{v}_{3}\) or \(-\mathbf{v}_{3}\) can be used as the normal vector of the plane. The unit normal vector \(\mathbf{n}=(n_{x},\,n_{y},\,n_{z})\) is computed by normalization of the eigenvector \(\mathbf{v}_{3}\) or \(-\mathbf{v}_{3}\).  
Once the tangent plane is calculated, the coordinates of the vertex \(\mathbf{p}=(x_{p},y_{p},\,z_{p})\) can be calculated by Equation (7) according to the given index of the cube \((i,j,k)\) and the cube's size \(L\).  
\[\left\{\begin{array}{l}x_{p}=i\times L\\ y_{p}=j\times L\\ z_{p}=k\times L\end{array}\right. \tag{7}\]  
The projection point \(\mathbf{p}_{1}\) in Section 2 is determined by the geometric relationship of tangent plane and \(\mathbf{p}\). By letting \(\mathbf{p}_{1}=(x,\,y,\,z)\), according to the geometric relation that \(\mathbf{p}_{1}\) is on the tangent plane, we can get the following expression:  
\[n_{x}(x-x_{c})+n_{y}(y-y_{c})+n_{z}(z-z_{c})=0, \tag{8}\]  
Consider that the projection direction is consistent with the normal direction of the tangent plane, the unit normal vector \(\mathbf{n}\) is proportional to the vector from \(\mathbf{p}\) to \(\mathbf{p}_{1}\).  
\[\frac{x-x_{p}}{n_{x}}=\frac{y-y_{p}}{n_{y}}=\frac{z-z_{p}}{n_{z}}, \tag{9}\]If Equation (9) is equal to a certain ratio \(r\), then the coordinates of projection point \(\mathbf{p}_{1}\) can be expressed by \(r\):  
\[\left\{\begin{array}{l}x=x_{p}+n_{x}r\\ y=y_{p}+n_{y}r\\ z=z_{p}+n_{z}r\end{array}\right., \tag{10}\]  
The ratio \(r\) can be calculated by replacing \(x\), \(y\), and \(z\) in Equation (8) with Equation (10)  
\[r=\frac{1}{n_{x}^{2}+n_{y}^{2}+n_{z}^{2}}(n_{x}(x_{c}-x_{p})+n_{y}(y_{c}-y_{p}) +n_{z}(z_{c}-z_{p})), \tag{11}\]  
Because \(\mathbf{n}\) is a unit normal vector, then  
\[n_{x}^{2}+n_{y}^{2}+n_{z}^{2}=1, \tag{12}\]  
Additionally, Equation (11) can be rewritten as:  
\[r=n_{x}(x_{c}-x_{p})+n_{y}(y_{c}-y_{p})+n_{z}(z_{c}-z_{p}), \tag{13}\]  
The coordinates of the projection point \(\mathbf{p}_{1}\) can be obtained by Equation (10)  
\[\left\{\begin{array}{l}x=x_{p}+n_{x}(n_{x}(x_{c}-x_{p})+n_{y}(y_{c}-y_{p})+ n_{z}(z_{c}-z_{p}))\\ y=y_{p}+n_{y}(n_{x}(x_{c}-x_{p})+n_{y}(y_{c}-y_{p})+n_{z}(z_{c}-z_{p}))\\ z=z_{p}+n_{z}(n_{x}(x_{c}-x_{p})+n_{y}(y_{c}-y_{p})+n_{z}(z_{c}-z_{p}))\end{array}\right., \tag{14}\]