International Journal of Biomedical Imaging

International Journal of Biomedical Imaging / 2011 / Article
Special Issue

Modern Mathematics in Biomedical Imaging

View this Special Issue

Research Article | Open Access

Volume 2011 |Article ID 285130 | https://doi.org/10.1155/2011/285130

Yuchuan Wei, Hengyong Yu, Ge Wang, "Inverse Fourier Transform in the Gamma Coordinate System", International Journal of Biomedical Imaging, vol. 2011, Article ID 285130, 16 pages, 2011. https://doi.org/10.1155/2011/285130

Inverse Fourier Transform in the Gamma Coordinate System

Academic Editor: Yangbo Ye
Received28 May 2010
Accepted29 Jul 2010
Published26 Oct 2010

Abstract

This paper provides auxiliary results for our general scheme of computed tomography. In 3D parallel-beam geometry, we first demonstrate that the inverse Fourier transform in different coordinate systems leads to different reconstruction formulas and explain why the Radon formula cannot directly work with truncated projection data. Also, we introduce a gamma coordinate system, analyze its properties, compute the Jacobian of the coordinate transform, and define weight functions for the inverse Fourier transform assuming a simple scanning model. Then, we generate Orlov's theorem and a weighted Radon formula from the inverse Fourier transform in the new system. Furthermore, we present the motion equation of the frequency plane and the conditions for sharp points of the instantaneous rotation axis. Our analysis on the motion of the frequency plane is related to the Frenet-Serret theorem in the differential geometry.

1. Introduction

In [1], we developed a general scheme for 2D and 3D, parallel- and divergent-beam computed tomography (CT). Different from traditional Radon's or Tuy's formulas starting with the inverse Fourier transform in the spherical coordinate system, our new framework is based on an instantaneous cylindrical coordinate system. However, in [1], the coordinate system was not formally defined, and the Jacobian of the transform was not explicitly calculated. Given the basic role of our general scheme in CT theory and some confusion in the literature [2–5], here, we treat the Ξ“ coordinate system in a mathematically strict way.

The organization of this paper is as follows. In Section 2, we demonstrate how the inverse Fourier transform formula generates different reconstruction formulas, with a special attention on different properties of projection truncation. In Section 3, we define the Ξ“ coordinate system and deduce its Jacobian factor. In Section 4, we study the motion of the frequency plane and present the condition for sharp points of the instantaneous rotation axis. In Section 5, a simplified model for helical CT reconstruction is given. In the appendixes, we present a weighted Ξ› reconstruction formula, a weighted Radon's formula, a deduction of Orlov's formula, and an alternative way to study the motion of the frequency plane. For convenience, we use the same notations as in [1].

2. Inverse Fourier Transform in Commonly Used Coordinate Systems

In this section, we write the inverse Fourier transform in commonly used coordinate systems, derive corresponding reconstruction formulas, and demonstrate their difference on truncation of projection data.

As shown in Figure 1, in the 3D space 𝑅3,Ξ¨(βƒ—π‘Ÿ) is an object function to be reconstructed, whose Fourier transform is ⃗Ψ(π‘˜). According to Fourier analysis, we have Ξ¨ξ€·ξ€Έ=ξ€œβƒ—π‘Ÿβˆžβˆ’βˆžξ€œβˆžβˆ’βˆžξ€œβˆžβˆ’βˆžξΞ¨ξ‚€βƒ—π‘˜ξ‚ξ‚€βƒ—ξ‚exp2πœ‹π‘–π‘˜β‹…βƒ—π‘Ÿπ‘‘π‘˜1π‘‘π‘˜2π‘‘π‘˜3,(1) where βƒ—π‘Ÿ=(π‘₯,𝑦,𝑧) and βƒ—π‘˜=(π‘˜1,π‘˜2,π‘˜3) are 3D vectors in the real space and Fourier domain, respectively. The unit vectors along the three axes are βƒ—βƒ—i,j, and βƒ—k. Note that the frequency vector βƒ—π‘˜ and the unit vector βƒ—k are different. This formula states that to reconstruct an 3D object Ξ¨(βƒ—π‘Ÿ), we need its Fourier transform ⃗Ψ(π‘˜) in the whole 3D space. In X-ray CT, the Fourier transform ⃗Ψ(π‘˜) can be calculated when sufficiently many parallel projections are measured.

Now, let us specify the simplest 3D reconstruction problem. In Figure 1(a), 𝑆(πœ™)=𝑆(sinπœ™,βˆ’cosπœ™,0) with πœ™βˆˆ[0,πœ‹] is a point on a great semicircle of the unit sphere. In the 3D space, we define the unit vectors ⃗𝑒=𝑂𝑆=(sinπœ™,βˆ’cosπœ™,0) and ξ‹π‘’βŸ‚=(cosπœ™,sinπœ™,0). Then, ⃗𝑒,ξ‹π‘’βŸ‚ and βƒ—k are another set of orthogonal unit vectors. As shown in Figure 1(c), the new components of the vector βƒ—π‘Ÿ can be denoted as 𝑒=βƒ—π‘Ÿβ‹…βƒ—π‘’,𝜌=βƒ—π‘Ÿβ‹…ξ‹π‘’βŸ‚βƒ—,𝑧=βƒ—π‘Ÿβ‹…k.(2) The projection at point 𝑆 is defined by 𝑃𝑆(𝜌,𝑧)=π‘ƒπœ™ξ€œ(𝜌,𝑧)=βˆžβˆ’βˆžΞ¨ξ‚€π‘’βƒ—π‘’+πœŒξ‹π‘’βŸ‚βƒ—k+𝑧𝑑𝑒.(3) Our purpose is to reconstruct Ξ¨(βƒ—π‘Ÿ) when π‘ƒπœ™(𝜌,𝑧) is known for all πœ™βˆˆ[0,πœ‹]. Let us see how the selection of a coordinate system will determine the reconstruction formula, paying an attention to longitudinal data truncation.

2.1. Signed Cylindrical Coordinate System

In a cylindrical coordinate system, the inverse Fourier transform can be rewritten as Ξ¨ξ€·ξ€Έ=ξ€œβƒ—π‘Ÿπœ‹0ξ€œβˆžβˆ’βˆžξ€œβˆžβˆ’βˆžξΞ¨ξ‚€βƒ—π‘˜ξ‚ξ‚€βƒ—ξ‚||π‘˜exp2πœ‹π‘–π‘˜β‹…βƒ—π‘ŸπœŒ||π‘‘π‘˜πœŒπ‘‘π‘˜3π‘‘πœ™,(4) where (π‘˜πœŒ,π‘˜3,πœ™) with βˆ’βˆž<π‘˜πœŒ<∞,βˆ’βˆž<π‘˜3<∞,0β‰€πœ™<πœ‹ is the signed cylindrical coordinates of βƒ—π‘˜. As shown in Figure 1(b), we have βƒ—βƒ—π‘˜ξ€·π‘˜π‘˜=𝜌,π‘˜3ξ€Έ,πœ™=π‘˜πœŒξ‹π‘’βŸ‚(πœ™)+π‘˜3βƒ—k.(5) When πœ™0 increases from 0 to Ο€, the plane πœ™=πœ™0 passes through every frequency point exactly once (except for the points on the axis π‘‚π‘˜3 which have a zero measure), and the normal vector 𝑂𝑆 of the plane πœ™=πœ™0 moves along a great semicircle on the unit sphere with the starting point A(0,βˆ’1,0) and the end point B(0,1,0). Note that the coordinate transform from (π‘˜1,π‘˜2,π‘˜3) to (π‘˜πœŒ,π‘˜3,πœ™) is one to one (except for the points on the axis π‘‚π‘˜3 with zero measure), and the absolute value of Jacobian is |π‘˜πœŒ|.

From the parallel projection π‘ƒπœ™(𝜌,𝑧), we can calculate the Fourier transform of the object ξΞ¨ξ€·π‘˜πœŒ,π‘˜3ξ€Έ=𝑃,πœ™πœ™ξ€·π‘˜πœŒ,π‘˜3ξ€Έ=ξ€œβˆžβˆ’βˆžξ€œβˆžβˆ’βˆžπ‘ƒπœ™ξ€Ίξ€·π‘˜(𝜌,𝑧)expβˆ’2πœ‹π‘–πœŒπœŒ+π‘˜3π‘§ξ€Έξ€»π‘‘πœŒπ‘‘π‘§.(6) Substituting (6) into (4), we have Ξ¨ξ€·ξ€Έ=ξ€œβƒ—π‘Ÿπœ‹0ξ€œβˆžβˆ’βˆžξ€œβˆžβˆ’βˆžξπ‘ƒπœ™ξ€·π‘˜πœŒ,π‘˜3⃗||π‘˜exp2πœ‹π‘–π‘˜β‹…βƒ—π‘ŸπœŒ||π‘‘π‘˜πœŒπ‘‘π‘˜3π‘‘πœ™.(7) The reconstruction scheme behind this formula can be denoted as π‘ƒπœ™ξΞ¨ξ‚€βƒ—π‘˜ξ‚ξ€·ξ€Έ(𝜌,𝑧)βŸΆβŸΆΞ¨βƒ—π‘Ÿ,(8) which indicates that we first calculate all the Fourier components from projections and then reconstruct the object using the inverse Fourier transform. In this process, the object is reconstructed as a whole body.

In fact, Formula (7) can be written in the filtered backprojection form as Ξ¨ξ€·ξ€Έ=ξ€œβƒ—π‘Ÿπœ‹0ξ‚π‘ƒπœ™(𝜌,𝑧)π‘‘πœ™,(9) where the filtered projection is given by ξ‚π‘ƒπœ™=ξ€œ(𝜌,𝑧)βˆžβˆ’βˆžξ€œβˆžβˆ’βˆžξπ‘ƒπœ™ξ€·π‘˜πœŒ,π‘˜3⃗||π‘˜exp2πœ‹π‘–π‘˜β‹…βƒ—π‘ŸπœŒ||π‘‘π‘˜πœŒπ‘‘π‘˜3=ξ€œβˆžβˆ’βˆžξ€œβˆžβˆžξπ‘ƒπœ™ξ€·π‘˜πœŒ,π‘˜3ξ€Έξ€·exp2πœ‹π‘–π‘˜πœŒπœŒξ€Έξ€·exp2πœ‹π‘–π‘˜3𝑧||π‘˜πœŒ||π‘‘π‘˜πœŒπ‘‘π‘˜3=π‘ƒπœ™1(𝜌,𝑧)βˆ—βˆ’2πœ‹2𝜌2.(10) The symbol βˆ— denotes the 1D convolution operation about variable 𝜌. Since the filtration is along the horizontal direction only, the projection can be longitudinally truncated, that is, a plane with a specified 𝑧 value can be reconstructed independently. The reconstruction scheme behind (9) with (10) can be denoted by ξΞ¨ξ‚€βƒ—π‘˜ξ‚βŸΆπ‘ƒπœ™ξ€·ξ€Έ(𝜌,𝑧)βŸΆΞ¨βƒ—π‘Ÿ,(11) which states that certain frequency components form the projection, and the object can be reconstructed directly from projections.

When projections are longitudinally truncated, the Fourier transform cannot be calculated. How and why can we reconstruct a part of the object?

The answer is as follows. What we want to reconstruct is the object Ξ¨(βƒ—π‘Ÿ) instead of the Fourier transform ⃗Ψ(π‘˜). Though we cannot calculate every single frequency component, we can calculate the contribution of the frequency components on the whole frequency plane (as well as its neighborhood) to a point to be reconstructed, see (10). In the aforementioned case, projections contain as all frequency components as the object does, but the ratio of the frequency components is not suitable. Then, the only thing we need to do is to perform a 1D filtration to adjust this ratio. An analog is that the Fourier components can be considered as raw material, projections a semifinished product, and the reconstructed object the final product. If we use the semifinished product to make a product, we do not need change it back to the raw material. The role of Fourier analysis is to tell us how far from the semifinished product to the final product. Simply speaking, the essence of image reconstruction from projections is filtration. Therefore, a suitable coordinate system should show the character of the filtration in an easy way.

2.2. Signed Spherical Coordinate System

In the signed spherical coordinate system shown in Figure 2, the inverse Fourier transform can be rewritten as Ξ¨ξ€·ξ€Έ=ξ€œβƒ—π‘Ÿπœ‹0ξ€œπœ‹/2βˆ’πœ‹/2ξ€œβˆžβˆ’βˆžξΞ¨ξ‚€βƒ—π‘˜ξ‚ξ‚€βƒ—ξ‚π‘˜exp2πœ‹π‘–π‘˜β‹…βƒ—π‘Ÿ2cosπœ—π‘‘π‘˜π‘‘πœ—π‘‘πœ™,(12) where (π‘˜,πœ—,πœ™) with βˆ’βˆž<π‘˜<∞,βˆ’πœ‹/2<πœ—<πœ‹/2,0β‰€πœ™<πœ‹ is the signed spherical coordinates of βƒ—π‘˜.

For convenience, we denote βƒ—βƒ—π‘˜=π‘˜(π‘˜,πœ—,πœ™)=π‘˜βƒ—π‘›, where ⃗𝑛=⃗𝑛(πœ—,πœ™)=ξ‹π‘’βŸ‚βƒ—(πœ™)cosπœ—+ksinπœ— is a unit vector.

For a unit vector ⃗𝑛=⃗𝑛(πœ—,πœ™), and a real number π‘™βˆˆ(βˆ’βˆž,∞), the Radon transform of the object is defined as an integral 𝑅Ψ𝑙,⃗𝑛=𝑅Ψ(𝑙,πœ—,πœ™)=βƒ—π‘Ÿβ‹…βƒ—π‘›=π‘™Ξ¨ξ€·ξ€Έπ‘‘βƒ—π‘Ÿ2βƒ—π‘Ÿ,(13) on the plane described by βƒ—π‘Ÿβ‹…βƒ—π‘›=𝑙.

The Fourier transform can be calculated from the Radon transform by ξΞ¨ξ‚€βƒ—π‘˜ξ‚=Ψ=ξ€œπ‘˜βƒ—π‘›βˆžβˆ’βˆžξ€·ξ€Έπ‘…Ξ¨π‘™,⃗𝑛exp(βˆ’2πœ‹π‘–π‘˜π‘™)𝑑𝑙.(14) Thus, the inverse Fourier transform (12) becomes Radon's formulaΞ¨ξ€·ξ€Έ1βƒ—π‘Ÿ=βˆ’4πœ‹2ξ€œπœ‹0ξ€œπœ‹/2βˆ’πœ‹/2πœ•2πœ•π‘™2||||𝑙=βƒ—π‘Ÿβ‹…βƒ—π‘›ξ€·ξ€Έ1𝑅Ψ𝑙,⃗𝑛cosπœ—π‘‘πœ—π‘‘πœ™=βˆ’4πœ‹2βˆ‡2ξ€œπœ‹0ξ€œπœ‹/2βˆ’πœ‹/2ξ€·ξ€Έ||𝑅Ψ𝑙,⃗𝑛𝑙=βƒ—π‘Ÿβ‹…βƒ—π‘›1cosπœ—π‘‘πœ—π‘‘πœ™=βˆ’8πœ‹2βˆ‡2ξ€œΞ©ξ€·ξ€Έ||𝑅Ψ𝑙,⃗𝑛𝑙=βƒ—π‘Ÿβ‹…βƒ—π‘›π‘‘βƒ—π‘›,(15) where βˆ‡2 is the Laplace operator, and Ξ© is the unit sphere.

Radon's formula tells us that to reconstruct the object at a point, what we need is the Radon transform of the planes through the point and its neighborhood. However, what we measure by the detector is line integrals along X-rays instead of the Radon transform [6] but it can be calculated from projections by ξ€·ξ€Έ=ξ€œπ‘…Ξ¨π‘™,⃗𝑛=𝑅Ψ(𝑙,πœ—,πœ™)βˆžβˆ’βˆžπ‘ƒπœ™(𝑙cosπœ—βˆ’π‘£sinπœ—,𝑙sinπœ—+𝑣cosπœ—)𝑑𝑣.(16) We call (16) the first relation between the Radon transform and parallel projection.

Formulas (15) and (16) form a scheme to reconstruct the object from its parallel projections. Here, the projection is not allowed to be longitudinally truncated. From Section 2.1, however, we have seen that the object can be reconstructed slice by slice from parallel projections. Why is the truncation not allowed when using the Radon formula?

For comparison, we rewrite Radon’s formula in the filtered backprojection format Ξ¨ξ€·ξ€Έ=ξ€œβƒ—π‘Ÿπœ‹0ξ‚π‘ƒπœ™(𝜌,𝑧)π‘‘πœ™.(17) with ξ‚π‘ƒπœ™ξ€œ(𝜌,𝑧)=βˆžβˆ’βˆžξ€œπœ‹/2βˆ’πœ‹/2ξΞ¨ξ‚€βƒ—π‘˜ξ‚ξ‚€βƒ—ξ‚π‘˜exp2πœ‹π‘–π‘˜β‹…βƒ—π‘Ÿ2=1cosπœ—π‘‘π‘˜π‘‘πœ—βˆ’4πœ‹2ξ€œπœ‹/2βˆ’πœ‹/2πœ•2πœ•π‘™2||||𝑙=βƒ—π‘Ÿβ‹…βƒ—π‘›=𝜌cosπœ—+𝑧sinπœ—ξ€·ξ€Έπ‘…Ξ¨π‘™,⃗𝑛cosπœ—π‘‘πœ—.(18) Comparing (18) with (10), we have the following relation between the Radon transform and parallel-beam projections 1βˆ’4πœ‹2ξ€œπœ‹/2βˆ’πœ‹/2πœ•2πœ•π‘™2||||𝑙=βƒ—π‘Ÿβ‹…βƒ—π‘›=𝜌cosπœ—+𝑧sinπœ—ξ€·ξ€Έπ‘…Ξ¨π‘™,⃗𝑛(πœ—,πœ™)cosπœ—π‘‘πœ—=π‘ƒπœ™1(𝜌,𝑧)βˆ—βˆ’2πœ‹2𝜌2.(19) We call this the second relation between the Radon transform and parallel-beam projections. Based on the relationship between fan-beam and parallel-beam projections [1], the right hand side of (19) can be expressed in terms of cone-beam projections as well.

Substituting (16) into (19), we get an identity of parallel projections 1βˆ’4πœ‹2ξ€œπœ‹/2βˆ’πœ‹/2πœ•2πœ•π‘™2||||𝑙=𝜌cosπœ—+𝑧sinπœ—Γ—ξ€œβˆžβˆ’βˆžπ‘ƒπœ™(𝑙cosπœ—βˆ’π‘£sinπœ—,𝑙sinπœ—+𝑣cosπœ—)𝑑𝑣cosπœ—π‘‘πœ—=π‘ƒπœ™(𝜌,𝑧)βˆ—1βˆ’2πœ‹2𝜌2.(20) Now, we see that reconstruction of an object from parallel-beam projections does allow truncation of projection data but to verify the identity (20) needs the information on the whole projection domain. Reconstruction from parallel-beam projections using the Radon formula means numerically verifying this identity. We point out that if the projection π‘ƒπœ™(𝜌,𝑧) is replaced by a general function with two variables, this identity holds true as well; see Appendix A.

Similarly, in cone-beam CT, reconstruction from cone-beam projections using the Radon formula, that is, the Grangeat's framework, means numerically verifying the cone-beam version of this identity. This is the reason why the Grangeat framework cannot deal with the long object problem.

We see that the Radon formula comes from the inverse Fourier transform in the signed spherical coordinate system. In this simplest reconstruction problem, there are two shortcomings with the Radon formula: (1) what we measure with X-ray is line integrals rather than planar integrals, (2) the projection cannot be longitudinally truncated.

2.3. Standard Spherical Coordinate System

In a standard spherical coordinate system, the inverse Fourier transform is written as Ξ¨ξ€·ξ€Έ=ξ€œβƒ—π‘Ÿ02πœ‹ξ€œπœ‹/2βˆ’πœ‹/2ξ€œβˆž0ξΞ¨ξ‚€βƒ—π‘˜ξ‚ξ‚€βƒ—ξ‚π‘˜exp2πœ‹π‘–π‘˜β‹…βƒ—π‘Ÿ2cosπœ—π‘‘π‘˜π‘‘πœ—π‘‘πœ™,(21) where (π‘˜,πœ—,πœ™) with  0β‰€π‘˜<∞,βˆ’πœ‹/2<πœ—<πœ‹/2,0β‰€πœ™<2πœ‹ is the standard spherical coordinates of βƒ—π‘˜, as shown in Figure 2.

From the well-known generalized function relations [7] ξ€œπ›Ώ(𝑙)=βˆžβˆ’βˆžπ‘–exp(2πœ‹π‘–π‘˜π‘™)π‘‘π‘˜,=ξ€œπœ‹π‘™βˆžβˆ’βˆžsgn(π‘˜)exp(2πœ‹π‘–π‘˜π‘™)π‘‘π‘˜,(22) we have 12𝑖𝛿(𝑙)+=ξ€œ2πœ‹π‘™βˆž01exp(2πœ‹π‘–π‘˜π‘™)π‘‘π‘˜,2π›Ώξ…žξ…žπ‘–(𝑙)+πœ‹π‘™3=βˆ’4πœ‹2ξ€œβˆž0π‘˜2exp(2πœ‹π‘–π‘˜π‘™)π‘‘π‘˜.(23) Based on the convolution theorem, we have ξ€œβˆž0ξΞ¨ξ‚€βƒ—π‘˜ξ‚ξ‚€βƒ—ξ‚π‘˜exp2πœ‹π‘–π‘˜β‹…βƒ—π‘Ÿ2=ξ€œπ‘‘π‘˜βˆž0ξΞ¨ξ€·ξ€Έπ‘˜βƒ—π‘›exp(2πœ‹π‘–π‘˜π‘™)π‘˜21π‘‘π‘˜=βˆ’8πœ‹2πœ•2πœ•2π‘™ξ€·ξ€Έβˆ’π‘–π‘…Ξ¨π‘™,⃗𝑛4πœ‹3𝑙3ξ€·ξ€Έ.βˆ—π‘…Ξ¨π‘™,⃗𝑛(24) Therefore, the inverse Fourier transform becomes Ξ¨ξ€·ξ€Έ1βƒ—π‘Ÿ=βˆ’8πœ‹2ξ€œ02πœ‹ξ€œπœ‹/2βˆ’πœ‹/2πœ•2πœ•π‘™2ξ€·ξ€Έ||||𝑅Ψ𝑙,⃗𝑛𝑙=βƒ—π‘Ÿβ‹…βƒ—π‘›βˆ’π‘–cosπœ—π‘‘πœ—π‘‘πœ™4πœ‹3ξ€œ02πœ‹ξ€œπœ‹/2βˆ’πœ‹/2ξ‚Έξ€·ξ€Έβˆ—1𝑅Ψ𝑙,⃗𝑛𝑙3ξ‚Ή||||𝑙=βƒ—π‘Ÿβ‹…βƒ—π‘›cosπœ—π‘‘πœ—π‘‘πœ™.(25) Based on the odd and even symmetry, we have the Radon formula again [6] Ξ¨ξ€·ξ€Έ1βƒ—π‘Ÿ=βˆ’4πœ‹2ξ€œπœ‹0ξ€œπœ‹/2βˆ’πœ‹/2πœ•2πœ•π‘™2ξ€·ξ€Έ||||𝑅Ψ𝑙,⃗𝑛𝑙=βƒ—π‘Ÿβ‹…βƒ—π‘›cosπœƒπ‘‘πœƒπ‘‘πœ™.(26) The difference between the signed and standard spherical coordinate systems is little, and the signed spherical coordinate system leads to Radon's formula more easily.

It is important to underline that the inverse Fourier transform may yield various reconstruct formulas in different coordinate systems and allow different degrees of data truncation. Among the commonly used coordinate systems, the signed cylindrical coordinate system is the most convenient one to solve the simplest 3D reconstruction problem. To handle the general 3D reconstruction problem, in the following we will introduce a variant of the signed cylindrical coordinate system, which is referred to as the Ξ“ coordinate system.

3. Ξ“ Coordinate System and Its Jacobian Factor

What plays an important role in the 3D reconstruction field is a variant of the cylindrical coordinate system, which can be named the Ξ“ coordinate system.

As shown in Figure 3(b),Ξ“AB is a three-times differentiable curve connecting points A(0,βˆ’1,0) and B(0,1,0) on the unit sphere, whose length is πœƒ0β‰₯πœ‹.𝑆(πœƒ) is a point on Ξ“AB parameterized by πœƒβˆˆ[0,πœƒ0] the length of the segment AS. Let us introduce three orthogonal unit vectors ⃗𝑒3=𝑂𝑆,⃗𝑒1=π‘‘π‘‘πœƒβƒ—π‘’3,⃗𝑒2=⃗𝑒3×⃗𝑒1,(27) where ⃗𝑒1 represents the tangential direction, ⃗𝑒2 the instantaneous rotation axis of the frequency plane. Note that they are functions of βˆ«πœƒ=𝑆𝐴|𝑑⃗𝑒3| (see Figures 3(a) and 3(b)).

In the Fourier domain, we call the plane through the origin 𝑂 and orthogonal to the vector 𝑂𝑆 the frequency plane, denoted as ∏(𝑂𝑆). Every point βƒ—βˆ(ξ‚”π‘˜βˆˆπ‘‚π‘†) can be expressed as βƒ—π‘˜=πœ”1⃗𝑒1+πœ”2⃗𝑒2,(28) with πœ”1=βƒ—π‘˜β‹…βƒ—π‘’1,πœ”2=βƒ—π‘˜β‹…βƒ—π‘’2.

More generally, for every given triplet (πœ”1,πœ”2,πœƒ) with πœ”1∈(βˆ’βˆž,∞),πœ”2∈(βˆ’βˆž,∞),πœƒβˆˆ[0,πœƒ0], we can define a 3D vector βƒ—βƒ—π‘˜ξ€·πœ”π‘˜=1,πœ”2ξ€Έ,πœƒ=πœ”1⃗𝑒1(πœƒ)+πœ”2⃗𝑒2(πœƒ),(29) and the triplet (πœ”1,πœ”2,πœƒ) is a set of Ξ“ coordinates of the point βƒ—π‘˜.

However, a given point βƒ—π‘˜βˆˆπ‘…3 may be represented by several sets of Ξ“ coordinates, because when 𝑆(πœƒ) moves along Ξ“AB, the rotation axis ⃗𝑒2 keeps changing, and hence the frequency plane ∏(𝑂𝑆) may scan the given point βƒ—π‘˜ more than once. The number of times of the point being scanned, denoted as ⃗𝐽(π‘˜), equals the number of intersections between the curve Ξ“AB and the great circle orthogonal to the vector βƒ—π‘˜ (Figure 3(c)). To find the all intersections, one can solve the equation βƒ—π‘˜β‹…βƒ—π‘’3(πœƒ)=0,(30) to obtain the solutions 1πœƒ,2πœƒ,…,⃗𝐽(π‘˜)πœƒ. For every solution, such as 1πœƒ, one can further have 1πœ”1=βƒ—π‘˜β‹…βƒ—π‘’1(1πœƒ),1πœ”2=βƒ—π‘˜β‹…βƒ—π‘’2(1πœƒ). Thus, one gets the ⃗𝐽(π‘˜) groups of coordinates βƒ—βƒ—π‘˜ξ€·π‘˜=1πœ”1,1πœ”2,1πœƒξ€Έ=βƒ—π‘˜ξ€·2πœ”1,2πœ”2,2πœƒξ€Έβƒ—π‘˜ξ‚΅=β‹―=π½ξ‚€βƒ—π‘˜ξ‚πœ”1,π½ξ‚€βƒ—π‘˜ξ‚πœ”2,π½ξ‚€βƒ—π‘˜ξ‚πœƒξ‚Ά.(31) We call this new system the Ξ“ coordinate system, since the system is based on the curve Ξ“AB. This coordinate system is a variant of the cylindrical system. For a differential arc π‘†π‘†ξ…ž, the new system is similar to a cylindrical coordinate system if the scanned volume is considered. When Γ𝐴𝐡 is a great semicircle, the Ξ“ coordinate system becomes the cylindrical coordinate system.

In the Ξ“ coordinate system, the inverse Fourier transform can be expressed as Ξ¨ξ€·ξ€Έ=ξ€œβƒ—π‘Ÿπœƒ00ξ€œβˆžβˆ’βˆžξ€œβˆžβˆ’βˆžπ‘Šξ€·πœ”1,πœ”2ξ€Έξ‚€βƒ—ξ‚Γ—ξΞ¨ξ‚€βƒ—π‘˜ξ‚||πœ”,πœƒexp2πœ‹π‘–π‘˜β‹…βƒ—π‘Ÿ1||π‘‘πœ”1π‘‘πœ”2π‘‘πœƒ,(32) with βƒ—π‘˜=πœ”1⃗𝑒1(πœƒ)+πœ”2⃗𝑒2(πœƒ). In [1], we obtained the absolute value of the Jacobian, |πœ”1| for the coordinate transform from (π‘˜1,π‘˜2,π‘˜3) to (πœ”1,πœ”2,πœƒ), by similarity between the Ξ“ coordinate system and the cylindrical system.

Based on (29) in the Ξ“ coordinate system, the absolute value of the Jacobian factor can be calculated as follows: ||||ξ‚΅πœ•βƒ—π‘˜πœ•πœ”1Γ—πœ•βƒ—π‘˜πœ•πœ”2ξ‚Άβ‹…πœ•βƒ—π‘˜||||=|||ξ€·πœ•πœƒβƒ—π‘’1×⃗𝑒2ξ€Έβ‹…ξ‚€πœ”1π‘‘π‘‘πœƒβƒ—π‘’1+πœ”2π‘‘π‘‘πœƒβƒ—π‘’2|||=|||⃗𝑒3β‹…ξ‚€πœ”1π‘‘π‘‘πœƒβƒ—π‘’1+πœ”2π‘‘π‘‘πœƒβƒ—π‘’2|||=||πœ”1||,(33) where we have used the relations discussed in the next section⃗𝑒3β‹…π‘‘π‘‘πœƒβƒ—π‘’1=βˆ’1,⃗𝑒3β‹…π‘‘π‘‘πœƒβƒ—π‘’2=0.(34) Then, a weight function π‘Š(πœ”1,πœ”2,πœƒ) is defined for all the coordinates (πœ”1,πœ”2,πœƒ), such that the summation of the weight function on the ⃗𝐽(π‘˜) groups of coordinates of the same βƒ—π‘˜βˆˆπ‘…3 remains one, that is, π½ξ‚€βƒ—π‘˜ξ‚ξ“π‘—=1π‘Šξ€·π‘—πœ”1,π‘—πœ”2,π‘—πœƒξ€Έ=1.(35) An example of the weight functions is [1, 8] π‘Šξ€·πœ”1,πœ”2ξ€Έξ€·π‘˜,πœƒ=sgn2ξ€Έξ€·πœ”sgn1ξ€Έ,(36) with π‘˜2=βƒ—βƒ—π‘˜β‹…π‘—=[πœ”1⃗𝑒1(πœƒ)+πœ”2⃗𝑒2βƒ—(πœƒ)]⋅𝑗.

With this weight function, (32) becomes a parallel-beam reconstruction formula from the Fourier slice theorem, and further Ye-Wang’s cone-beam formula [8] via the relationship between parallel-beam and cone-beam projections [1]. With the weight function, the filtration is no longer along the tangential direction generally.

Clearly, the weight function is not unique. For instance, for frequency points with ⃗𝐽(π‘˜)=3, in addition to the weight function defined by (36), (1βˆ’1+1), we can also choose (1/3+ 1/3+1/3), (0+1+0), or a more complicated one, (1+1βˆ’1) for some points and (βˆ’1+1+1) for the rest points, leading to different parallel-beam reconstruction formulas. Furthermore, using the approach described in [1], we can arrive at Katsevich’s helical cone-beam formulas [9–11], and so forth.

If one takes the weight function π‘Š(πœ”1,πœ”2,πœƒ)=1 on the right side of (32), what we will reconstruct is actually another function Ξ¦ξ€·ξ€Έ=ξ€œβƒ—π‘Ÿπœƒ00ξ€œβˆžβˆ’βˆžξ€œβˆžβˆ’βˆžξ‚€βƒ—ξ‚ξΞ¨ξ‚€βƒ—π‘˜ξ‚||πœ”exp2πœ‹π‘–π‘˜β‹…βƒ—π‘Ÿ1||π‘‘πœ”1π‘‘πœ”2=ξ€œπ‘‘πœƒπ‘…3ξ‚€βƒ—ξ‚ξΞ¨ξ‚€βƒ—π‘˜ξ‚π½ξ‚€βƒ—π‘˜ξ‚π‘‘exp2πœ‹π‘–π‘˜β‹…βƒ—π‘Ÿ3βƒ—π‘˜.(37) From (37), we observe that the reconstructed image Ξ¦(βƒ—π‘Ÿ) equals the true object Ξ¨(βƒ—π‘Ÿ) if and only if ⃗𝐽(π‘˜)=1, that is, the trajectory Ξ“AB is a great semicircle. The condition for Ξ¦(βƒ—π‘Ÿ) to be an approximation of the object function, that is, Ξ¦(βƒ—π‘Ÿ)β‰ˆΞ¨(βƒ—π‘Ÿ), is that ⃗𝐽(π‘˜)=1 in a major portion of the Fourier domain and ⃗𝐽(π‘˜) is close to 1 in the rest region. Based on the relation between ⃗𝐽(π‘˜) and the curve length of Ξ“AB, a necessary condition is πœƒ0β‰ˆπœ‹, or equivalently πœ€=πœƒ0/πœ‹βˆ’1β‰ˆ0. Otherwise, Ξ¦(βƒ—π‘Ÿ) and Ξ¨(βƒ—π‘Ÿ)will have a great difference. Unfortunately this formula was once thought as exact, which was essentially Theorem 3 in [2] and Theorem 4.3 in [3]. The details can be seen in [1].

Generally speaking, the criterion to identify exact formulas from approximate ones is to evaluate if every frequency component contributes exactly once to the final reconstruction, and the difference between exact formulas comes from different forms of weight functions.

4. Motion of the Frequency Plane and Condition for Sharp Points

In this section, we will study motion of the frequency plane, present the condition for sharp points of the locus of the rotation axis, and introduce the Euler angles. Let us simplify notations ̇⃗𝑒3=(𝑑/π‘‘πœƒ)⃗𝑒3,Μˆβƒ—π‘’3=(𝑑2/π‘‘πœƒ2)⃗𝑒3, and so on. In this section, we assume that Ξ“ is a three times differentiable curve on the unit ball, which may be not complete.

Obviously we have the following.

Proposition 1. For any point 𝑆(πœƒ)βˆˆΞ“, ⃗𝑒1,⃗𝑒2,⃗𝑒3 defined by (27) are orthogonal unit vectors.

Definition 1. For every point 𝑆(πœƒ)βˆˆΞ“, let us define 𝜎=⃗𝑒3⋅̇⃗𝑒3Γ—Μˆβƒ—π‘’3ξ€Έ,(38) which is the signed volume spanned by the three vectors ⃗𝑒3,̇⃗𝑒3,Μˆβƒ—π‘’3. If we consider the curve length πœƒ as time, ⃗𝑒3,̇⃗𝑒3,Μˆβƒ—π‘’3 are the position, velocity, and the acceleration of point S.

Proposition 2. The motion of ⃗𝑒1,⃗𝑒2,⃗𝑒3 can be described by the following differential equations π‘‘π‘‘πœƒβƒ—π‘’3=⃗𝑒1,π‘‘π‘‘πœƒβƒ—π‘’1=πœŽβƒ—π‘’2βˆ’βƒ—π‘’3,π‘‘π‘‘πœƒβƒ—π‘’2=βˆ’πœŽβƒ—π‘’1.(39)

Proof. The first equation is the definition. Now, let us prove the other two. Since ⃗𝑒1,⃗𝑒2,⃗𝑒3 are orthogonal unit vectors, we can write ⃗𝑒1⋅⃗𝑒1=1,⃗𝑒2⋅⃗𝑒2=1,⃗𝑒3⋅⃗𝑒3=1,⃗𝑒1⋅⃗𝑒2=0,⃗𝑒2⋅⃗𝑒3=0,⃗𝑒3⋅⃗𝑒1=0.(40) Applying 𝑑/π‘‘πœƒ on the above equations, one has ̇⃗𝑒1⋅⃗𝑒1Μ‡=0,⃗𝑒2⋅⃗𝑒2Μ‡=0,⃗𝑒3⋅⃗𝑒3Μ‡=0,⃗𝑒1⋅⃗𝑒2+⃗𝑒1⋅̇⃗𝑒2Μ‡=0,⃗𝑒2⋅⃗𝑒3=0,⃗𝑒3β‹…Μ‡βƒ—e1=βˆ’1.(41) Using the definition of 𝜎, we have ̇⃗𝑒1⋅⃗𝑒2=̇⃗𝑒1⋅⃗𝑒3×⃗𝑒1ξ€Έ=⃗𝑒3⋅̇⃗𝑒3Γ—Μˆβƒ—π‘’3ξ€Έ=𝜎,⃗𝑒1⋅̇⃗𝑒2=βˆ’πœŽ.(42) Using the orthogonal unit vectors ⃗𝑒1,⃗𝑒2,⃗𝑒3, ̇⃗𝑒1 and ̇⃗𝑒2 can be expressed as π‘‘π‘‘πœƒβƒ—π‘’1=πœŽβƒ—π‘’2βˆ’βƒ—π‘’3,π‘‘π‘‘πœƒβƒ—π‘’2=βˆ’πœŽβƒ—π‘’1,(43) which completes the proof.

From these differential equations, we recognize that the motion of ⃗𝑒1,⃗𝑒2,⃗𝑒3 is a rotation in the 3D space.

Now, we consider the frequency plane ξ‚”Ξ (𝑂𝑆) fixed with ⃗𝑒1,⃗𝑒2, ⃗𝑒3 as a rigid body. When 𝑆(πœƒ) moves along the curve Ξ“, the rigid body rotates in the 3D space.

Let us define βƒ—Ξ©=⃗𝑒2+πœŽβƒ—π‘’3.(44) We will see that βƒ—Ξ© is the angular velocity of the rigid body by the following proposition.

Proposition 3. The motion of the vectors ⃗𝑒1,⃗𝑒2,⃗𝑒3 can be expressed as π‘‘π‘‘πœƒβƒ—π‘’3=⃗Ω×⃗𝑒3,π‘‘π‘‘πœƒβƒ—π‘’1=⃗Ω×⃗𝑒1,π‘‘π‘‘πœƒβƒ—π‘’2=⃗Ω×⃗𝑒2.(45)

Proof. It is easy to verify the three equations one by one with Proposition 2 and the definition of βƒ—Ξ©.

Proposition 4. The angular acceleration of the rigid body is Μ‡βƒ—Ξ©=Μ‡πœŽβƒ—π‘’3(πœƒ).(46)

Proof. In fact, we have ̇⃗𝑑Ω=ξ€·π‘‘πœƒβƒ—π‘’2+πœŽβƒ—π‘’3ξ€Έ=̇⃗𝑒2+Μ‡πœŽβƒ—π‘’3Μ‡+πœŽβƒ—π‘’3=βˆ’πœŽβƒ—π‘’1+Μ‡πœŽβƒ—π‘’3+πœŽβƒ—π‘’1=Μ‡πœŽβƒ—π‘’3.(47)

Proposition 5. For any point 𝑆(πœƒ)βˆˆΞ“, we have ̈𝜎=0βŸΊβƒ—π‘’3=βˆ’βƒ—π‘’3,⃛(48)Μ‡πœŽ=0βŸΊβƒ—π‘’3ξ€·πœŽ=βˆ’2ξ€ΈΜ‡+1⃗𝑒3.(49)

Proof. Let us first prove (48).
If 𝜎=0, by Proposition 2 we have Μˆβƒ—π‘’3=βˆ’βƒ—π‘’3. On the other hand, if Μˆβƒ—π‘’3=βˆ’βƒ—π‘’3, we have 𝜎=⃗𝑒3Μ‡β‹…[⃗𝑒3Γ—Μˆβƒ—π‘’3]=0.
Now, let us prove (49).
By Proposition 2, we have ̇⃗𝑒1=πœŽβƒ—π‘’2βˆ’βƒ—π‘’3.(50) Applying 𝑑/π‘‘πœƒ to both sides of (50) and using Μ‡πœŽ=0, we have Μˆβƒ—π‘’1=Μ‡πœŽβƒ—π‘’2Μ‡+πœŽβƒ—π‘’2βˆ’Μ‡βƒ—π‘’3ξ€·=πœŽβˆ’πœŽβƒ—π‘’1ξ€Έβˆ’βƒ—π‘’1ξ€·πœŽ=βˆ’2ξ€Έ+1⃗𝑒1.(51) On the other hand, if ⃛⃗𝑒3=βˆ’(𝜎2Μ‡+1)⃗𝑒3, we have π‘‘Μ‡πœŽ=ξ€Ίπ‘‘πœƒβƒ—π‘’3⋅̇⃗𝑒3Γ—Μˆβƒ—π‘’3=̇⃗𝑒3⋅̇⃗𝑒3Γ—Μˆβƒ—π‘’3ξ€Έ+⃗𝑒3β‹…ξ€·Μˆβƒ—π‘’3Γ—Μˆβƒ—π‘’3ξ€Έ+⃗𝑒3⋅̇⃗𝑒3×⃛⃗𝑒3ξ€Έ=⃗𝑒3⋅̇⃗𝑒3×⃛⃗𝑒3ξ€Έ=0.(52) This completes the proof. Physically, ⃛⃗𝑒3 is the change rate of the acceleration Μˆβƒ—π‘’3.

Proposition 6. 𝜎=0 for every point 𝑆(πœƒ)βˆˆΞ“, if and only if Ξ“ is a great circle or a great arc.

Proof. Since 𝜎=0 for every point 𝑆(πœƒ)βˆˆΞ“, the motion equation of the three unit vectors, (39), is simplified as π‘‘π‘‘πœƒβƒ—π‘’3=⃗𝑒1,π‘‘π‘‘πœƒβƒ—π‘’1=βˆ’βƒ—π‘’3,π‘‘π‘‘πœƒβƒ—π‘’2=0.(53) Since (𝑑/π‘‘πœƒ)⃗𝑒2=0 for every point, ⃗𝑒2 is a fixed unit vector. Since (𝑑/π‘‘πœƒ)⃗𝑒3=⃗𝑒1=⃗𝑒2×⃗𝑒3, ⃗𝑒3 rotates about the fixed axis ⃗𝑒2. Since ⃗𝑒3βŸ‚βƒ—π‘’2, 𝑆, the endpoint of ⃗𝑒3, draws a great circle or a great arc; see Figure 4(a).
On the other hand, if Ξ“ is a great circle or a great arc, one can verify that Μˆβƒ—π‘’3=βˆ’βƒ—π‘’3. Therefore, we have 𝜎=⃗𝑒3Μ‡β‹…(⃗𝑒3Γ—Μˆβƒ—π‘’3)=0. This completes the proof.

Proposition 7. Μ‡πœŽ=0 for every point 𝑆(πœƒ)βˆˆΞ“, if and only if Ξ“ is a circle or a circular arc.

Proof. Since Μ‡πœŽ=0, by Proposition 4 we have Μ‡βƒ—Ξ©=Μ‡πœŽβƒ—π‘’3(πœƒ)=0. We know that βƒ—Ξ© is a fixed vector. By Proposition 3, we have π‘‘π‘‘πœƒβƒ—π‘’3=⃗Ω×⃗𝑒3.(54) Therefore, ⃗𝑒3 rotates about the fixed axis βƒ—Ξ©, and 𝑆, the endpoint of ⃗𝑒3, draws a circular arc. Furthermore, the radius of the circle is √1/1+𝜎2; see Figure 4(b).
On the other hand, if Ξ“ is a circular arc, the line connecting the circle center and the origin 𝑂 is the rotation axis. We know that Μˆβƒ—π‘’3 is pointing to the axis. By circular symmetry, 𝜎=⃗𝑒3Μ‡β‹…(⃗𝑒3Γ—Μˆβƒ—π‘’3) is constant at every point of the circular arc. Therefore, Μ‡πœŽ=0 on the circular arc. This completes the proof.

Clearly, Proposition 6 is a special case of Proposition 7; see Figure 4(a).

One important case is that 𝜎=⃗𝑒3Μ‡β‹…(⃗𝑒3Γ—Μˆβƒ—π‘’3)=0, for a certain point 𝑆(πœƒ)βˆˆΞ“. At such a moment, we have (𝑑/π‘‘πœƒ)⃗𝑒2=βˆ’πœŽβƒ—π‘’1=0, the rotation axis of the frequency plane does not move, the point 𝑆 moves along a great arc, and the Ξ“ coordinate system is more similar to a cylindrical coordinate system. We say that πœƒ is a stationary point of ⃗𝑒2 if 𝜎(πœƒ)=0. If 𝜎 has different signs before and after 𝜎=0, the stationary point will become a sharp point, since in this case ⃗𝑒2 goes forward, stops and goes back according to (𝑑/π‘‘πœƒ)⃗𝑒2=βˆ’πœŽβƒ—π‘’1. Therefore, we have the following straightforward proposition.

Proposition 8. A necessary condition for the curve ⃗𝑒2(πœƒ) to have a sharp point at πœƒ1∈(0,πœƒ0) is 𝜎(πœƒ1)=0; a sufficient condition for the curve ⃗𝑒2(πœƒ) to have a sharp point at πœƒ1∈(0,πœƒ0) is that 𝜎(πœƒ) has different signs for πœƒ<πœƒ1 and for πœƒ>πœƒ1.

Based on the definition of 𝜎 and the equation of the osculating plane, 𝜎(πœƒ1)=0 means that the osculating plane at 𝑆(πœƒ1)βˆˆΞ“ passes through the origin since (⃗𝑒3Μ‡βˆ’0)β‹…(⃗𝑒3Γ—Μˆβƒ—π‘’3)=0.𝜎(πœƒ) changes its sign before and after πœƒ1 means that the osculating plane sweeps the origin, or relatively speaking, the origin goes from one side of the osculating plane to the other. An example can be seen in the next section. Sharp points serve as the landmarks for defining weight functions.

Proposition 9. For any given point on the frequency plane βƒ—π‘˜=πœ”1⃗𝑒1+πœ”2⃗𝑒2, the motion equation is π‘‘βƒ—βƒ—βƒ—π‘˜π‘‘πœƒπ‘˜=Ω×=βˆ’πœŽπœ”2⃗𝑒1(πœƒ)+πœŽπœ”1⃗𝑒2(πœƒ)βˆ’πœ”1⃗𝑒3(πœƒ).(55)

Proof. It is straightforward from Propositions 2 or 3.
We call (55) the motion equation of the frequency plane.
There are two applications of Proposition 9. First, it can be used to compute the Jacobian of the Ξ“ coordinate system ||||ξ‚΅πœ•βƒ—π‘˜πœ•πœ”1Γ—πœ•βƒ—π‘˜πœ•πœ”2ξ‚Άβ‹…πœ•βƒ—π‘˜||||=||||πœ•πœƒβƒ—π‘’3β‹…πœ•βƒ—π‘˜||||=||πœ•πœƒβˆ’πœ”1||=||πœ”1||.(56) Second, since ⃗𝑒1β‹…πœ•βƒ—πœ•πœƒπ‘˜=βˆ’πœŽπœ”2,⃗𝑒2β‹…πœ•βƒ—πœ•πœƒπ‘˜=πœŽπœ”1,(57) we know that among the three vectors πœ•βƒ—π‘˜/πœ•πœ”1,πœ•βƒ—π‘˜/πœ•πœ”2,πœ•βƒ—π‘˜/πœ•πœƒ, the last one is not orthogonal to the first two. Hence, a general Ξ“ coordinate system is not orthogonal. Notice that πœ•βƒ—βƒ—π‘˜/πœ•πœƒ=πœ•π‘˜/πœ•πœƒ|πœ”1,πœ”2βƒ—=π‘‘π‘˜/π‘‘πœƒ
Note that the frequency plane rotates about the axis ⃗𝑒2 with the unit angular velocity while the vectors ⃗𝑒1 and ⃗𝑒2 themselves rotate about the axis ⃗𝑒3 with the angular velocity 𝜎; see Figure 5(a). We see that the Ξ“ coordinate system is different from the cylindrical coordinate system in the motion process. However, since the spin about the axis ⃗𝑒3 of the frequency plane does not have an effect on the scanned volume, their Jacobians have the same expression. Figure 5(b) shows the shape of a differential volume of the Ξ“ coordinate system.

Proposition 10. When S moves along a curve Ξ“ on the unit sphere, the scanned volume in the unit ball by the frequency plane is 𝑉=4πœƒ03=ξ€œ||βƒ—π‘˜||<1π½ξ‚€βƒ—π‘˜ξ‚π‘‘3βƒ—π‘˜,(58) where πœƒ0 is the length of Ξ“, and ⃗𝐽(π‘˜) is the number of times point βƒ—π‘˜ is scanned.

Proof. In fact, the scanned volume can be calculated in both the Ξ“ and the Cartesian coordinate systems as follows: ξ€œπ‘‰=πœƒ00ξ€œπ‘˜2=πœ”21+πœ”22<1||πœ”1||π‘‘πœ”1π‘‘πœ”2=ξ€œπ‘‘πœƒπœƒ00ξ€œ1π‘˜=0ξ€œ2πœ‹πœ™=0π‘˜2||||=4cosπœ™π‘‘π‘˜π‘‘πœ™π‘‘πœƒ3ξ€œπœƒ004π‘‘πœƒ=3πœƒ0=ξ€œπ‘˜21+π‘˜22+π‘˜23<1π½ξ‚€βƒ—π‘˜ξ‚π‘‘3βƒ—π‘˜.(59) This completes the proof.

This concise formula provides a basic relation between the curve length, the scanned volume by the frequency plane, and the J function.

When we study the rotation problem, Euler angles (πœ—,πœ™,𝛿) are helpful. As shown in Figure 6(a), (πœ—,πœ™) with πœ—βˆˆ(βˆ’πœ‹/2,πœ‹/2),πœ™βˆˆ(βˆ’βˆž,∞) are the spherical coordinate of 𝑆=𝑆(πœ—,πœ™),𝛿 is the angle between βƒ—π‘’πœ— and ⃗𝑒2. ξ‹π‘’πœ™, and βƒ—π‘’πœ— are the latitude and longitude unit vectors at point 𝑆 in the spherical coordinate system.

The Euler angles are functions of the curve length πœ—=πœ—(πœƒ),πœ™=πœ™(πœƒ),𝛿=𝛿(πœƒ).(60) The three unit vectors can be represented as follows: ⃗𝑒3=⃗𝑒3(πœ—,πœ™),⃗𝑒1=sinπ›Ώβƒ—π‘’πœ—+cosπ›Ώξ‹π‘’πœ™=Μ‡πœ—βƒ—π‘’πœ—+Μ‡πœ™cosπœ—ξ‹π‘’πœ™,⃗𝑒2=βˆ’cosπ›Ώβƒ—π‘’πœ—+sinπ›Ώξ‹π‘’πœ™Μ‡=βˆ’πœ™cosπœ—βƒ—π‘’πœ—+Μ‡πœ—ξ‹π‘’πœ™.(61) Furthermore, one can verify the motion equations of the three unit vectors in Proposition 2.

Based on the motion theory of rigid body, the angular velocity of a rigid body can be written as βƒ—Μ‡Ξ©=βˆ’πœ—ξ‹π‘’πœ™+Μ‡πœ™βƒ—Μ‡π‘˜+𝛿⃗𝑒3Μ‡=βˆ’πœ—ξ‹π‘’πœ™+Μ‡πœ™cosπœ—βƒ—π‘’πœ—+ξ€·Μ‡Μ‡π›Ώξ€Έπœ™sinπœ—+⃗𝑒3=⃗𝑒2+πœŽβƒ—π‘’3.(62) Therefore, we have the following relationships:⃗𝑒2Μ‡=βˆ’πœ—ξ‹π‘’πœ™+Μ‡πœ™cosπœ—βƒ—π‘’πœ—=βƒ—Ξ©0βˆ’ξ‚€βƒ—Ξ©0⋅⃗𝑒3⃗𝑒3,Μ‡Μ‡βƒ—Ξ©πœŽ=πœ™sinπœ—+𝛿=0⋅⃗𝑒3βˆ’Μ‡βƒ—Ξ©0⋅⃗𝑒1,̇̇⃗Ω𝛿=βˆ’0⋅⃗𝑒1=βƒ—Ξ©0⋅̇⃗𝑒1,(63) where we have used a new notation βƒ—Ξ©0Μ‡=βˆ’πœ—ξ‹π‘’πœ™+Μ‡πœ™βƒ—π‘˜,(64) which is the first two terms of the total angular velocity βƒ—Ξ© of (62).

In fact, βƒ—Ξ©0 is the angular velocity of the detector with the unit vectors ξ‹π‘’πœ™,βƒ—π‘’πœ— when 𝑆 moves along Ξ“. To be specific, we consider the detector for parallel-beam projection as a rectangle of finite width and length through the origin orthogonal to direction ⃗𝑒3=𝑂𝑆. During any motion, one side of the detector always keeps horizontal, that is, parallel to the direction ξ‹π‘’πœ™. The pose of the detector is completely determined by the point 𝑆=𝑆(πœ—,πœ™), or two Euler's angles (πœ—,πœ™). The relative motion of the frequency plane to the detector is described by the angle 𝛿; see Figure 6(b).

In summary, there are three sets of orthogonal vectors (beside βƒ—βƒ—βƒ—ki,j,, the unit vectors of the Cartesian coordinate system) and three important planes. The orthogonal unit vectors ξ‹π‘’πœ™,βƒ—π‘’πœ—, ⃗𝑒3 are related to the movement of the detector and ⃗𝑒3, ̇⃗𝑒3=⃗𝑒1, ⃗𝑒2=⃗𝑒3×⃗𝑒1 are related to the movement of the frequency plane, see Figure 6(b). It is easy to verify that βƒ—Ξ©=⃗𝑒2+πœŽβƒ—π‘’3=⃗𝑒1×̇⃗𝑒1 and hence the vectors ⃗𝑒1, ̇⃗𝑒1, βƒ—Ξ© are the third set of orthogonal vectors. The direction of βƒ—Ξ© is the normal of the osculating plane and its magnitude, √1+𝜎2, is the reciprocal of the radius of the osculating circle. Here the osculating plane is spanned by ⃗𝑒1, ̇⃗𝑒1 [15] and the osculating circle is referred to as the intersection between the osculating plane and the unit ball. When S moves along a curve Ξ“ on the unit ball, it can be viewed that the point S moves along an osculating circle, whose orientation and radius are described by βƒ—Ξ©. From βƒ—Ξ©=⃗𝑒2+πœŽβƒ—π‘’3, we can see that 𝜎 is an angular velocity of the spin of the frequency plane. The locus lengths of the endpoints of ⃗𝑒3,⃗𝑒1, ⃗𝑒2 are πœƒ0, βˆ«πœƒ00√1+𝜎2π‘‘πœƒ, βˆ«πœƒ00|𝜎|π‘‘πœƒ, respectively, according to Proposition 2.

5. Example of a Curve with ⃗𝐽(π‘˜)=1,3

Here, we discuss the reconstruction based on the simplified model of the 3D helical cone beam reconstruction, that is, 3D helical parallel-beam reconstruction.

The curve 𝐢 in 𝑅3 defined by 𝑆𝐢=ξ…žβˆˆπ‘…3βˆΆξ‚”π‘‚π‘†ξ…žξ‚€β„Ž=⃗𝑦(πœ‘)=π‘Ÿcosπœ‘,π‘Ÿsinπœ‘,πœ‘ξ‚,ξ‚ƒβˆ’πœ‹2πœ‹πœ‘βˆˆ2,πœ‹2(65) is a segment of a helix with radius π‘Ÿ and pitch β„Ž, which is shown in Figure 7 for π‘Ÿ=1,β„Ž=3.

In fact, it is formed by a combination of a uniform circular motion in the XOY plane and a uniform straight motion along the 𝑍 direction.

The projection of the curve 𝐢 on the unit sphere Ξ© is ξ‚»ξ‚”Ξ“=π‘†βˆˆΞ©βˆΆπ‘‚π‘†=⃗𝑒3(πœ‘)=⃗𝑦(πœ‘)||||ξ‚ƒβˆ’πœ‹βƒ—π‘¦(πœ‘),πœ‘βˆˆ2,πœ‹2ξ‚„ξ‚Ό,(66) or ξ‚»ξ‚”Ξ“=π‘†βˆˆΞ©βˆΆπ‘‚π‘†=⃗𝑒3(πœ—,πœ™),πœ—=tgβˆ’1ξ‚΅β„Žπœ™ξ‚Άξ‚ƒβˆ’πœ‹2πœ‹,πœ™βˆˆ2,πœ‹2ξ‚„ξ‚Ό.(67) The curve length πœƒ can be given by ||π‘‘πœƒ=𝑑⃗𝑒3||=||||𝑑(πœ‘)⃗𝑦(πœ‘)||||ξ‚Ά||||,ξ€œβƒ—π‘¦(πœ‘)πœƒ(πœ‘)=πœ‘βˆ’πœ‹/2π‘‘πœƒπ‘‘πœ‘ξ…žπ‘‘πœ‘ξ…ž=ξ€œπœ‘βˆ’πœ‹/2|||||π‘‘ξƒ©ξ€·πœ‘βƒ—π‘¦ξ…žξ€Έ||⃗𝑦(πœ‘ξ…ž)||ξƒͺ|||||.(68) For example, we can calculate πœƒ0=πœƒ(πœ‹/2)β‰ˆ3.177 and πœ€=0.011, when π‘Ÿ=1 and β„Ž=3, indicating that the overlapping degree in the Fourier domain is small.

For every πœ‘βˆˆ[βˆ’πœ‹/2,πœ‹/2], the three orthogonal unit vectors are⃗𝑒3=𝑂𝑆=⃗𝑦(πœ‘)||||,⃗𝑦(πœ‘)⃗𝑒1=(𝑑/π‘‘πœ‘)⃗𝑒3||(𝑑/π‘‘πœ‘)⃗𝑒3||,⃗𝑒2=⃗𝑒3×⃗𝑒1.(69) The inverse Fourier transform can be expressed with either the Ξ“ coordinates (πœ”1,πœ”2,πœƒ) or the new coordinates (πœ”1,πœ”2,πœ‘)Ξ¨ξ€·ξ€Έ=ξ€œβƒ—π‘Ÿπœƒ00ξ€œβˆžβˆ’βˆžξ€œβˆžβˆ’βˆžπ‘Šξ€·πœ”1,πœ”2ξ€Έξ‚€βƒ—ξ‚ξΞ¨ξ‚€βƒ—π‘˜ξ‚Γ—||πœ”,πœƒexp2πœ‹π‘–π‘˜β‹…βƒ—π‘Ÿ1||π‘‘πœ”1π‘‘πœ”2=ξ€œπ‘‘πœƒπœ‹/2βˆ’πœ‹/2ξ€œβˆžβˆ’βˆžξ€œβˆžβˆ’βˆžπ‘Šξ€·πœ”1,πœ”2ξ€Έξ‚€βƒ—ξ‚Γ—ξΞ¨ξ‚€βƒ—π‘˜ξ‚||πœ”,πœƒ(πœ‘)exp2πœ‹π‘–π‘˜β‹…βƒ—π‘Ÿ1||π‘‘πœƒπ‘‘πœ‘π‘‘πœ”1π‘‘πœ”2π‘‘πœ‘,(70) with βƒ—βƒ—π‘˜ξ€·πœ”π‘˜=1,πœ”2ξ€Έ,πœƒ(πœ‘)=πœ”1⃗𝑒1(πœƒ(πœ‘))+πœ”2⃗𝑒2(πœƒ(πœ‘)).(71) To visualize the motion of the frequency plane ∏(𝑂𝑆), its five positions are marked in Figure 7 at πœ‘1=βˆ’πœ‹/2, πœ‘2=βˆ’πœ‹/4, πœ‘3=0,πœ‘4=πœ‹/4 and πœ‘5=πœ‹/2, respectively. Note that the initial position Π(𝑂𝑆1) and the final position Π(𝑂𝑆5) coincide with each other.

In Figure 7, the two curved sides of the red triangle on the top of the unit sphere is the locus of ⃗𝑒3, the instantaneous rotation axis of the frequency plane. Similarly, βˆ’βƒ—π‘’2 and Π(𝑂𝑆1) form the counterpart triangle on the bottom of the unit sphere. The locus of ⃗𝑒2 is the division line of the regions with different ⃗𝐽(π‘˜) values. Based on (58), the total area of the two triangles is only 0.55% of the surface area of the unit sphere.

Note that a sharp point on the locus of ⃗𝑒2 appears at πœ‘=πœ‘3=0. The reason is as follows.

Let us define𝜎0𝑑(πœ‘)=⃗𝑦(πœ‘)β‹…π‘‘π‘‘πœ‘βƒ—π‘¦(πœ‘)Γ—2π‘‘πœ‘2⃗𝑦(πœ‘).(72) It can be verified that 𝜎0(0)=0, 𝜎0(πœ‘)<0 for πœ‘<0, and 𝜎0(πœ‘)>0 for πœ‘>0.

It is not difficult to prove 𝜎=⃗𝑒3⋅⃗𝑒1×̇⃗𝑒1𝑑=⃗𝑦(πœ‘)β‹…π‘‘π‘‘πœ‘βƒ—π‘¦(πœ‘)Γ—2π‘‘πœ‘2ξ‚Ά1⃗𝑦(πœ‘)||||Γ—1⃗𝑦(πœ‘)||(||𝑑/π‘‘πœ‘)⃗𝑦(πœ‘)2ξ‚΅π‘‘πœ‘ξ‚Άπ‘‘πœƒ3=𝜎0(πœ‘)𝐾(πœ‘),(73) where 𝐾(πœ‘) is positive. Therefore, 𝜎 and 𝜎0 have the same zeros and sign.

Based on the motion process of ξ‚”Ξ (𝑂𝑆), we have π½ξ‚€βƒ—π‘˜ξ‚=ξƒ―βƒ—βƒ—1,π‘˜outsidethetwotriangles,3,π‘˜insidethetwotriangles.(74) The overlapped region in the frequency space is due to the motion of the instantaneous rotation axis ⃗𝑒2 of the frequency plane, or equivalently the discrepancy between the curve Ξ“ and half a great circle.

As shown in Figure 8, the weight factor can be (1βˆ’1+1), (1/3+ 1/3+1/3), or (0+1+0), and so forth, for the overlapped region. By dividing the overlapped region into 2 symmetric parts, one can define the weight function as (βˆ’1+1+1) on the left side and (1+1βˆ’1) on the right side, or (1βˆ’1+1) on the left and (1+1βˆ’1) on the right, and so forth. Hence, one can design various parallel-beam reconstruction formulas and derive the associated cone-beam formulas using the approach described in [1].

Appendices

A. Weighted Ξ› Reconstruction Formula

𝑓(π‘₯,𝑦) is a function in the 2D space with the Fourier transform ξπ‘“ξ€·π‘˜1,π‘˜2ξ€Έ=ξ€œβˆžβˆ’βˆžξ€œβˆžβˆ’βˆžξ€·π‘“(π‘₯,𝑦)expβˆ’π‘–2πœ‹π‘˜1π‘₯ξ€Έξ€·expβˆ’π‘–2πœ‹π‘˜2𝑦𝑑π‘₯𝑑𝑦.(A.1) Its parallel-beam projection is denoted as π‘ƒπœ—ξ€œ(𝑙)=βˆžβˆ’βˆžπ‘“(𝑙cosπœ—βˆ’π‘£sinπœ—,𝑙sinπœ—+𝑣cosπœ—)𝑑𝑣(A.2) with π‘™βˆˆ(βˆ’βˆž,∞),πœ—βˆˆ[βˆ’πœ‹/2,πœ‹/2].

Let 𝑓(π‘₯,𝑦) pass through the filter |π‘˜1|, we get a new function =ξ€œΞ¦(π‘₯,𝑦)βˆžβˆ’βˆžξ€œβˆžβˆ’βˆž||π‘˜1||ξπ‘“ξ€·π‘˜1,π‘˜2ξ€Έξ€·exp𝑖2πœ‹π‘˜1π‘₯ξ€Έξ€·exp𝑖2πœ‹π‘˜2π‘¦ξ€Έπ‘‘π‘˜1π‘‘π‘˜21=𝑓(π‘₯,𝑦)βˆ—βˆ’2πœ‹2π‘₯2.(A.3) Defining the signed polar coordinates (π‘˜πœŒ,πœ—) with π‘˜πœŒβˆˆ(βˆ’βˆž,∞),πœ—βˆˆ(βˆ’πœ‹/2,πœ‹/2), which satisfy π‘˜1=π‘˜πœŒcosπœ—,π‘˜2=π‘˜πœŒsinπœ—, we haveξ€œΞ¦(π‘₯,𝑦)=πœ‹/2βˆ’πœ‹/2ξ€œβˆžβˆ’βˆž||π‘˜1||ξπ‘“ξ€·π‘˜πœŒξ€Έξ€Ί,πœ—Γ—expβˆ’π‘–2πœ‹π‘˜πœŒξ€»||π‘˜(π‘₯cosπœ—+𝑦sinπœ—)𝜌||π‘‘π‘˜πœŒ=ξ€œπ‘‘πœ—πœ‹/2βˆ’πœ‹/2ξ€œβˆžβˆ’βˆžξπ‘“ξ€·π‘˜πœŒξ€Έξ€Ί,πœ—expβˆ’π‘–2πœ‹π‘˜πœŒξ€»(π‘₯cosπœ—+𝑦sinπœ—)Γ—π‘˜2𝜌cosπœ—π‘‘π‘˜πœŒ=1π‘‘πœ—βˆ’4πœ‹2ξ€œπœ‹/2βˆ’πœ‹/2πœ•2πœ•π‘™2||||𝑙=π‘₯cosπœ—+𝑦sinπœ—ξ€œβˆžβˆ’βˆžξπ‘“ξ€·π‘˜πœŒξ€Έξ€·,πœ—Γ—expβˆ’π‘–2πœ‹π‘˜πœŒπ‘™ξ€Έπ‘‘π‘˜πœŒ=1cosπœ—π‘‘πœ—βˆ’4πœ‹2ξ€œπœ‹/2βˆ’πœ‹/2πœ•2πœ•π‘™2||||𝑙=π‘₯cosπœ—+𝑦sinπœ—π‘ƒπœ—(𝑙)cosπœ—π‘‘πœ—.(A.4) Comparing (A.3) and (A.4), we have1𝑓(π‘₯,𝑦)βˆ—βˆ’2πœ‹2π‘₯2=1βˆ’4πœ‹2ξ€œπœ‹/2βˆ’πœ‹/2πœ•2πœ•π‘™2||||𝑙=π‘₯cosπœ—+𝑦sinπœ—π‘ƒπœ—(𝑙)cosπœ—π‘‘πœ—.(A.5) That is, we have rediscovered the identity (20) in a general way. In contrast to the classical Ξ› reconstruction formula [1], a weight factor cosπœ— is involved in (A.5). Hence, we can call (A.5) a weighted Ξ› reconstruction formula.

B. A Weighted Radon Formula and the Orlov Formula

Let us define 𝑒=βƒ—π‘Ÿβ‹…βƒ—π‘’3,𝑑=βƒ—π‘Ÿβ‹…βƒ—π‘’1,𝑠=βƒ—π‘Ÿβ‹…βƒ—π‘’2,(B.1) and the parallel-beam projection is π‘ƒπœƒξ€œ(𝑑,𝑠)=βˆžβˆ’βˆžΞ¨ξ€·π‘’βƒ—π‘’3+𝑑⃗𝑒1+𝑠⃗𝑒2𝑑𝑒,(B.2) as shown in Figure 6(b). In the frequency plane ξ‚”Ξ (𝑂𝑆), we define a signed polar coordinate system (πœ”,𝛼) with πœ”βˆˆ(βˆ’βˆž,∞),π›Όβˆˆ[0,πœ‹] for a point (πœ”1,πœ”2) with πœ”1=πœ”cos𝛼,πœ”2=πœ”cos𝛼. Similarly, let us define the signed polar coordinates (𝑋,𝛽) with π‘‹βˆˆ(βˆ’βˆž,∞),π›½βˆˆ[βˆ’πœ‹/2,πœ‹/2] for a point (𝑑,𝑠) with 𝑑=𝑋cos𝛽,𝑠=𝑋sin𝛽.

Let the weight function in (32) be π‘Š(πœ”1,πœ”2,πœƒ)=π‘Š(𝛼,πœƒ), we have Ξ¨ξ€·ξ€Έ=ξ€œβƒ—π‘Ÿπœƒ00ξ€œβˆžβˆ’βˆžξ€œβˆžβˆ’βˆžπ‘Šξ€·πœ”1,πœ”2ξ€Έξ‚€βƒ—ξ‚ξΞ¨ξ‚€βƒ—π‘˜ξ‚Γ—||πœ”,𝛼exp2πœ‹π‘–π‘˜β‹…βƒ—π‘Ÿ1||π‘‘πœ”1π‘‘πœ”2=ξ€œπ‘‘πœƒπœƒ00ξ€œπœ‹0ξ€œβˆžβˆ’βˆžπ‘Šξ€·(𝛼,πœƒ)|cos𝛼|exp2πœ‹π‘–πœ”1𝑑+2πœ‹π‘–πœ”2π‘ ξ€ΈΓ—ξΞ¨ξ€·πœ”1,πœ”2ξ€Έπœ”,πœƒ2=ξ€œπ‘‘πœ”π‘‘π›Όπ‘‘πœƒπœƒ00ξ€œπœ‹0π‘Šξ€œ(𝛼,πœƒ)|cos𝛼|βˆžβˆ’βˆžπœ”2×=exp(2πœ‹π‘–πœ”(𝑑cos𝛼+𝑠sin𝛼))Ξ¨(πœ”,𝛼,πœƒ)π‘‘πœ”π‘‘π›Όπ‘‘πœƒβˆ’14πœ‹2ξ€œπœƒ00ξ€œπœ‹0πœ•π‘Š(𝛼,πœƒ)2πœ•π‘™2||||𝑙=βƒ—π‘Ÿβ‹…βƒ—π‘›=𝑑cos𝛼+𝑠sin𝛼=𝑅Ψ𝑙,⃗𝑛×|cos𝛼|π‘‘π›Όπ‘‘πœƒβˆ’14πœ‹2βˆ‡2ξ€œπœƒ00ξ€œπœ‹0ξ€·ξ€ΈΓ—π‘Š(𝛼,πœƒ)𝑅Ψ𝑑cos𝛼+𝑠sin𝛼,⃗𝑛|cos𝛼|π‘‘π›Όπ‘‘πœƒ,(B.3) with ⃗𝑛=⃗𝑒1cos𝛼+⃗𝑒2sin𝛼.

We call it the weighted Radon formula. When Ξ“AB is a great semicircle, we have π‘Š(𝛼,πœƒ)=1 and πœƒ0=πœ‹, and the above formula becomes the classic Radon formula [12, 13].

Moreover, it is easy to rederive Orlov's formula [14] from (32) Ξ¨ξ€·ξ€Έ=ξ€œβƒ—π‘Ÿπœƒ00ξ€œβˆžβˆ’βˆžξ€œβˆžβˆ’βˆžπ‘Šξ€·πœ”1,πœ”2ξ€Έξ‚€βƒ—ξ‚ξΞ¨ξ‚€βƒ—π‘˜ξ‚Γ—||πœ”,πœƒexp2πœ‹π‘–π‘˜β‹…βƒ—π‘Ÿ1||π‘‘πœ”1π‘‘πœ”2=ξ€œπ‘‘πœƒπœƒ00ξ€œβˆžβˆ’βˆžξ€œβˆžβˆ’βˆžξ€·π‘Š(𝛼,πœƒ)exp2πœ‹π‘–πœ”1𝑑+2πœ‹π‘–πœ”2π‘ ξ€ΈΓ—ξΞ¨ξ€·πœ”1,πœ”2ξ€Έ||πœ”,πœƒ1||π‘‘πœ”1π‘‘πœ”2=ξ€œπ‘‘πœƒπœƒ00ξ€œβˆžβˆ’βˆžξ€œβˆžβˆ’βˆžβˆ’14πœ‹2ξ‚΅πœ•2πœ•π‘‘2+πœ•2πœ•π‘ 2ξ‚Άξ€·π‘Š(𝛼,πœƒ)Γ—exp2πœ‹π‘–πœ”1𝑑+2πœ‹π‘–πœ”2π‘ ξ€ΈξΞ¨ξ€·πœ”1,πœ”2ξ€Έ||πœ”,πœƒ1||πœ”2π‘‘πœ”1π‘‘πœ”21π‘‘πœƒ=βˆ’4πœ‹2βˆ‡2ξ€œπœƒ00ξ€œ0π‘₯0200π‘‘βˆžβˆ’βˆžξ€œβˆžβˆ’βˆž|cos𝛼||ξΞ¨ξ€·πœ”πœ”|π‘Š(𝛼,πœƒ)β‹…1,πœ”2ξ€Έξ€·,πœƒΓ—exp2πœ‹π‘–πœ”1𝑑+2πœ‹π‘–πœ”2π‘ ξ€Έπ‘‘πœ”1π‘‘πœ”21π‘‘πœƒ=βˆ’4πœ‹2βˆ‡2ξ€œπœƒ00π‘πœƒ(𝑑,𝑠)βˆ—βˆ—π‘”πœƒ(𝑑,𝑠)π‘‘πœƒ,(B.4) where π‘”πœƒ=ξ€œ(𝑑,𝑠)βˆžβˆ’βˆžξ€œβˆžβˆ’βˆž|cos𝛼|ξ€·|πœ”|π‘Š(𝛼,πœƒ)exp2πœ‹π‘–πœ”1𝑑+2πœ‹π‘–πœ”2π‘ ξ€Έπ‘‘πœ”1π‘‘πœ”2=|cos𝛼|||𝑋||||π‘Š(𝛼,πœƒ)𝛼=πœ‹/2+𝛽.(B.5) The sign βˆ—βˆ— denotes the 2D convolution operation.

We recall that a frequency point has ⃗𝐽(k) sets of Ξ“ coordinates βƒ—βƒ—π‘˜ξ€·π‘˜=1πœ”1,1πœ”2,1πœƒξ€Έ=βƒ—π‘˜ξ€·2πœ”1,2πœ”2,2πœƒξ€Έβƒ—π‘˜ξ‚΅=β‹―=π½ξ‚€βƒ—π‘˜ξ‚πœ”1,π½ξ‚€βƒ—π‘˜ξ‚πœ”2,π½ξ‚€βƒ—π‘˜ξ‚πœƒξ‚Ά.(B.6) If a given triplet (πœ”1,πœ”2,πœƒ) is the 𝑗th set of Ξ“ coordinates of βƒ—π‘˜, ξ€·πœ”1,πœ”2ξ€Έ=ξ€·,πœƒπ‘—πœ”1,π‘—πœ”2,π‘—πœƒξ€Έ,(B.7) a new normalized weight function [1] for the triplet can be defined by π‘Šξ€·πœ”1,πœ”2ξ€Έ||βƒ—,πœƒ=π‘Š(𝛼,πœƒ)=1/π‘˜β‹…βƒ—π‘’1ξ€·π‘—πœƒξ€Έ||βˆ‘βƒ—π½(π‘˜)𝑗=1||βƒ—1/π‘˜β‹…βƒ—π‘’1ξ€·π‘—πœƒξ€Έ||=1/|cos𝛼|βˆ‘ξ‹π‘˜π½(0)𝑗=1||ξ‹π‘˜1/0⋅⃗𝑒1ξ€·π‘—πœƒξ€Έ||,(B.8) and in this case we have π‘”πœƒ(1𝑑,𝑠)=||𝑋||1βˆ‘ξ‹π‘˜π½(0)𝑗=1||ξ‹π‘˜1/0⋅⃗𝑒1ξ€·π‘—πœƒξ€Έ||||||||ξ‹π‘˜0=⃗𝑒3×𝑋0.(B.9) Here, ξ‹π‘˜0=βƒ—βƒ—π‘˜/|π‘˜|, and 𝑋0=(𝑑⃗𝑒1+𝑠⃗𝑒2√)/𝑑2+𝑠2=⃗𝑒1cos𝛽+⃗𝑒2sin𝛽 are unit vectors.

Equation (B.4) with (B.9) is the Orlov formula for a general curve on the unit sphere ((2a) with (2c) in [14]). Note that in [14] different letters were used such as ξ‹π‘˜π‘šβ†”0,πœβ†”βƒ—π‘’3,Μ‡πœβ†”βƒ—π‘’1,π‘₯ξ…žβ†”(𝑑⃗𝑒1+𝑠⃗𝑒2),|π‘₯ξ…žβˆ«|↔|𝑋|,πΊπ‘‘πΊπœβ†”βˆ«πœƒ00π‘‘πœƒ, and so forth, and Μ‡πœξ…žπ‘˜ should be printed as Μ‡πœπ‘˜.

To study the property of a curve Ξ“ on the unit sphere and the associated motion of the frequency plane, one can introduce another curve, Γ, defined as ξ€œβƒ—π‘Ÿ(πœƒ)=πœƒ0⃗𝑒3ξ€·πœƒξ…žξ€Έπ‘‘πœƒξ…ž+βƒ—π‘Ÿ(0),(C.1) where βƒ—π‘Ÿ(0) is an arbitrary vector in the 3D space, and πœƒ (as well as πœƒξ…ž) is an arc-length parameter of the curve Ξ“. Clearly, Μ‡π‘‘βƒ—π‘Ÿ=π‘‘πœƒβƒ—π‘Ÿ(πœƒ)=⃗𝑒3,Μˆβƒ—π‘Ÿ=⃗𝑒1,βƒ›Μ‡βƒ—π‘Ÿ=⃗𝑒1.(C.2) Suppose that the arc-length parameter of the curve Γ is 𝑠. The tangential direction of the curve Γ is ⃗𝑑T=π‘‘π‘ βƒ—π‘Ÿ=π‘‘πœƒπ‘‘π‘‘π‘ π‘‘πœƒβƒ—π‘Ÿ(πœƒ)=π‘‘πœƒπ‘‘π‘ βƒ—π‘’3(πœƒ).(C.3) Furthermore, we have |||𝑑|||=|||π‘‘π‘ βƒ—π‘Ÿπ‘‘πœƒ|||||𝑑𝑠⃗𝑒3||,|||(πœƒ)1=π‘‘πœƒ|||,𝑑𝑠(C.4) where |(𝑑/𝑑𝑠)βƒ—π‘Ÿ|=|π‘‘βƒ—π‘Ÿ/𝑑𝑠|=|π‘‘βƒ—π‘Ÿ|/𝑑𝑠=𝑑𝑠/𝑑𝑠=1 has been used.

When the positive directions are the same, we have π‘‘πœƒπ‘‘π‘ =1,πœƒ=𝑠+𝑠0,βƒ—T=⃗𝑒3(πœƒ).(C.5) Picking up the constant 𝑠0=0, we have πœƒ=𝑠. Therefore, the symbol πœƒ represents the arc length parameters of both Γ and Ξ“. Based on the tangential direction βƒ—T, the norm and binormal unit vectors [15] can be defined as βƒ—Μ‡βƒ—TN=|||Μ‡βƒ—T|||=⃗𝑒1||⃗𝑒1||=⃗𝑒1,βƒ—βƒ—βƒ—B=TΓ—N=⃗𝑒3×⃗𝑒1=⃗𝑒2.(C.6) The curvature πœ… and the torsion 𝜏 of the curve Γ are ||̈||=||πœ…=βƒ—π‘Ÿβƒ—π‘’1||Μ‡ξ€·Μˆβƒ›ξ€Έ=1,𝜏=βƒ—π‘Ÿβ‹…βƒ—π‘ŸΓ—βƒ—π‘Ÿ=⃗𝑒3⋅⃗𝑒1×̇⃗𝑒1ξ€Έ=𝜎,(C.7) and the Darboux vector is βƒ—βƒ—βƒ—Ξ©=𝜏T+πœ…B=⃗𝑒2+πœŽβƒ—π‘’3.(C.8) The Frenet-Serret Equations of the curve Γ̇⃗⃗̇⃗⃗⃗̇⃗⃗T=N,N=βˆ’T+𝜎B,B=βˆ’πœŽN,(C.9) and Μ‡βƒ—βƒ—βƒ—Μ‡βƒ—βƒ—βƒ—Μ‡βƒ—βƒ—βƒ—T=Ω×T,N=Ω×N,B=Ω×B,(C.10)

are (39) and (45). In this frame, βƒ—π‘Ÿ(πœƒ) is the original curve (function), and ⃗𝑒3(πœƒ) becomes its derived curve (function). In practice, it is difficult to find the physical meaning of the 3D curve βƒ—π‘Ÿ(πœƒ).

Example C. A circle on the unit sphere and its corresponding curve in the 3D space.

As shown in Figure 4(b), a circle on the unit sphere, marked as 𝐢, can be described as ⃗𝑒