Microlocal analysis of a spindle transform

An analysis of the stability of the spindle transform, introduced in ("Three dimensional Compton scattering tomography"arXiv:1704.03378 [math.FA]), is presented. We do this via a microlocal approach and show that the normal operator for the spindle transform is a type of paired Lagrangian operator with"blowdown--blowdown"singularities analogous to that of a limited data synthetic aperture radar (SAR) problem studied by Felea et. al. ("Microlocal analysis of SAR imaging of a dynamic reflectivity function"SIAM 2013). We find that the normal operator for the spindle transform belongs to a class of distibutions $I^{p,l}(\Delta\cup\widetilde{\Delta},\Lambda)$ studied by Felea and Marhuenda ("Microlocal analysis of SAR imaging of a dynamic reflectivity function"SIAM 2013 and"Microlocal analysis of some isospectral deformations"Trans. Amer. Math.), where $\widetilde{\Delta}$ is reflection through the origin, and $\Lambda$ is associated to a rotation artefact. Later, we derive a filter to reduce the strength of the image artefact and show that it is of convolution type. We also provide simulated reconstructions to show the artefacts produced by $\Lambda$ and show how the filter we derived can be applied to reduce the strength of the artefact.


Introduction
Here we present a microlocal analysis of the spindle transform, first introduced by the authors in [1], which describes the Compton scattering tomography problem in three dimensions for a monochromatic source and energy sensitive detector pair. Compton scattering is the process in which a photon interacts in an inelastic collision with a charged particle. As the collision is inelastic, the photon undergoes a loss in energy, described by the equation where E s is the energy of the scattered photon which had an initial energy E λ , ω is the scattering angle and E 0 ≈ 511keV is the electron rest energy. For E s and E λ fixed (i.e if the source is monochromatic and we can measure E s ), the scattering angle ω remains fixed and, in three dimensions, the surface of scatterers is the surface of revolution of a circular arc [1]. The surface of revolution of a circular arc is a spindle torus. We define to be the spindle torus, radially symmetric about the x 3 axis, with tube centre offset r ≥ 0 and tube radius √ 1 + r 2 . See figure 1 which displays a rotation of T r . In [12] Norton considered the problem of reconstructing a density supported in a quadrant of the plane from the Compton scattered intensity measured at a single point detector moved laterally along the axis away from a point source at the origin. Here the curve of scatterers is a circle. He considers the circle transform where F (ρ, θ) = f (ρ cos θ, ρ sin θ) is the polar form of f : . Then Af (r, φ) = Rf ( 1 r , φ), where R denotes the polar form of the straight line Radon transform. So A is equivalent to R via the diffeomorphism x → x |x| 2 and from this we can derive stability estimates from known theory on the Radon transform [13].
In [10,11], Nguyen and Truong consider an acquisition geometry of a point source and detector which remain opposite one another and are rotated on S 1 , and aim to reconstruct a density supported on the unit disc. Here the curve of scatterers is a circular arc. They define the circular arc transform Lettingf (x) = √ |x| 2 +1−|x| 1−|x| 2 f |x| 2 + 1 − |x| · x |x| , we have that r √ 1+r 2 Bf (r, φ) = Af (r, φ), and hence B is equivalent to A via the diffeomorphism x → |x| 2 + 1 − |x| x |x| . So, in two dimensions, the inverse problem for Compton scattering tomography is injective and mildly ill posed (the solution is bounded in some Sobolev space). In [10], Palamodov also derives stability estimates for a more general class of Minkowski-Funk transforms.
In [1] the authors consider a three dimensional acquisition geometry, where a single source and detector are rotated opposite one another on S 2 and a density supported on a hollow ball is to be recovered. Here the surface of scatterers is a spindle torus. They define the spindle transform Sf (r, θ) = 2π 0 π 0 ρ 2 sin ϕ 1 + r 2 1 + r 2 sin 2 ϕ (h · F ) (ρ, ψ, ϕ) | ρ= √ r 2 sin 2 ϕ+1−r sin ϕ dϕdψ, (5) where F (ρ, ψ, ϕ) = f (ρ cos ψ sin ϕ, ρ sin ψ sin ϕ, ρ cos ϕ) is the spherical polar form of f : R 3 → R and h ∈ SO(3) describes the rotation of the north pole to θ, where h defines a group action on real-valued functions in the natural way (h · f )(x) = f (hx). They show that a left inverse to S exists through the explicit inversion of a set of onedimensional Volterra integral operators, and show that the null space of S consists of those functions whose even harmonic components are zero (odd functions). However the stability of the spindle transform was not considered. We aim to address this here from a microlocal perspective. In [2], various acquisition geometries are considered for synthetic aperture radar imaging of moving objects. In each case the microlocal properties of the forward operator in question and its normal operator are analysed. In Figure 1: A spindle torus with axis of rotation θ, tube centre offset r and tube radius √ 1 + r 2 . The distance between the origin and either of the points where the torus self intersects is 1. three of the four cases considered the Schwartz kernel of the normal operator was shown to belong to a class of distributions associated to two cleanly intersecting Lagrangians I p,l (∆, Λ) (that is, the wavefront set of the kernel of the normal operator is contained in ∆ ∪ Λ). We show a similar result for the spindle transform S * S, although the diagonal ∆ is replaced by the disjoint union ∆ ∪ ∆ where ∆ is reflection through the origin. We also determine the associated Lagrangian Λ. In [2] they suggest a way to reduce the size of the image artefact microlocally by applying an appropriate pseudodifferential operator as a filter before applying the backprojection operator. Similarly we derive a suitable filter for the spindle transform and show how it can be applied using the spherical harmonics of the data.
In section 2.1 we show that S is equivalent to a weighted cylinder transform C, which gives the weighted integrals over cylinders with an axis of revolution through the origin. After this we prove that C is a Fourier integral operator and determine its canonical relation. Later in section 2.2 we present our main theorem (Theorem 3), where we show that C belongs to a class of distributions I p,l (∆ ∪ ∆, Λ), where ∆ is a reflection and the Lagrangian Λ is associated to a rotation artefact.
In section 3, we adopt the ideas of Felea et al in [2] and derive a suitable pseudodifferential operator Q which, when applied as a filter before applying the backprojection operator of the cylinder transform, reduces the artefact intensity in the image. We show that Q can be applied by multiplying the harmonic components of the data by a factor c l , which depends on the degree l of the component, and show how this translates to a spherical convolution of the data with a distribution on the sphere h.
Simulated reconstructions from spindle transform data are presented in section 4. We reconstruct a small bead of constant density by unfiltered backprojection and show the artefacts produced by Λ in our reconstruction. We then reconstruct the same density by filtered backprojection, applying the filter Q as an intermediate step, and show how the size of the artefacts are reduced in the image. Later we provide reconstructions of densities of oscillating layers using the conjugate gradient least squares (CGLS) method and Landweber iteration. We arrange the layers as spherical shells centred at the origin and as planes and compare our results. We also investigate the effects of applying the filter Q as a pre-conditioner, prior to implementing CGLS and the Landweber method.

The microlocal properties of S and S * S
Here we investigate the microlocal properties of the spindle transform and its normal operator. We start by showing the equivalence of S to a cylinder transform C, and how we can write S and C as Fourier integral operators. Then we determine the canonical relations associated with C and from these we discover that C * C is a paired Lagrangian operator with blowdown-blowdown singularities. First we give some preliminaries.
Let B n 1 , 2 = {x ∈ R n : 0 < 1 < |x| < 2 < 1} denote the set of points on a hollow ball with inner radius 1 and outer radius 2 . Let Z n = R×S n−1 denote the n-cylinder and, for X ⊂ R n an open set, let D (X) denote the vector space of distributions on X, and let E (X) denote the vector space of distributions with compact support contained in X.

Definition 1.
For a function f in the Schwarz space S(R n ) we define the Fourier transform and its inverse in terms of angular frequency as Definition 2. Let m, ρ, δ ∈ R with 0 ≤ ρ ≤ 1 and δ = 1 − ρ. Then we define S m ρ (X × R n ) to be the set of a ∈ C ∞ (X × R n ) such that for every compact set K ⊂ X and all multi-indices α, β the bound holds for some constant C α,β,K . The elements of S m ρ are called symbols of order m, type ρ.
where φ is a phase function, and a ∈ S m ρ ((X × Y ) × R N ) is a symbol. Definition 5. The canonical relation of an FIO with phase function φ is defined as If Y and X are manifolds without boundary, then an operator A : is a Fourier integral operator if its Schwartz kernel can be represented locally in coordinates by oscillatory integrals of the form (8), and the canonical relations of the phase functions for the local representations all lie within a single immersed Lagrangian submanifold of T * Y × T * X. For much more detail on Fourier integral operators and their definition see [4].

The spindle transform as an FIO
Recall the author's acquisition geometry in [1] (displayed in figure 1). We have the implicit equation (r + |x × θ|) 2 + (x · θ) 2 = 1 + r 2 for the set of points on a spindle with tube centre offset r and axis of revolution given by θ ∈ S 2 . With this in mind we define and then we can write the spindle transform S : where s = 1/r 2 and δ is the Dirac-delta function. Note that is smooth, bounded above, and does not vanish on on (0, 1) × B 3 1 , 2 × S 2 . We define the backprojection operator S * : where dΩ is the surface measure on S 2 . Proposition 1. The backprojection operator S * is the adjoint operator to S.
. Then (in the third step note that ∇ x h does not actually depend on s) which completes the proof. Let wheref We define the weighted cylinder transform C : and its backprojection operator C * : As in Proposition 1, we can show that C * is the formal adjoint to C.
The above is to say that the spindle transform is equivalent, via the diffeomorphism to the transform C which defines the weighted integrals over cylinders with radius √ s and axis of rotation through the origin with direction θ. With this in mind we consider the microlocal properties of the cylinder transform C and its normal operator C * C for the remainder of this section.
First, we characterise C as a Fourier integral operator in the next theorem.
Proof. The delta function may be written as the oscillatory integral Thus, by equation (18) we have and from this we see that C is an FIO with phase function and amplitude where the single phase variable is σ. Indeed, as we noted above and is evident from the formula (13), |∇ v h(s, v, θ)| is smooth, bounded from above, and bounded from below larger than zero when Therefore, since a also does not depend on the phase variable σ, a ∈ S 0 1 ((B α 1 ,α 2 × (0, 1) × S 2 ) × R). Hence the order of C is 0 + 1 2 − 1 4 (3 + 3) = −1. Now suppose that θ ∈ S 2 is parametrized by α and β ∈ R (for example using standard spherical coordinates). Then and the derivatives of φ are From Definition 5, it follows that the canonical relation of C is: which completes the proof.

C * C as a paired Lagrangian operator
If we analyse the canonical relation C given in Theorem 1, we can see that it is noninjective as the points on the ring {x ∈ R 3 : |x| 2 −s = 0, x·θ = 0} map to ((s, θ), (σ, 0)) if we fix s and σ. Let C * be the canonical relation of C * and let ∆ denote the diagonal. Then, given the non-injectivity of C, C * • C ∆ and C * C is not a pseudodifferential operator, or even an FIO. In this section we show that the Schwarz kernel of C * C instead belongs to a class of distributions I p,l (∆ ∪ ∆, Λ) studied in [2,3]. First we recall some definitions and theorems from [2].
Recall the definitions of the left and right projections of a canonical relation.
Definition 8. Let C be the canonical relation associated to the FIO A : E (X) → D (Y ). Then we denote π L and π R to be the left and right projections of C, π L : C → T * Y \0 and π R : C → T * X\0.
We have the following result from [4].
. Then at any point in C: (i) if one of π L or π R is a local diffeomorphism, then C is a local canonical graph; (ii) if one of the projections π R or π L is singular, then so is the other. The type of the singularity may be different (e.g. fold or blowdown [6]) but both projections drop rank on the same set Now we have the definition of a blowdown singularity and the definitions of a nonradial and involutive submanifold: Let M and N be manifolds of dimension n and let f : N → M be a smooth function. f is said to have a blowdown singularity of order k ∈ N along a smooth hypersurface Σ ⊂ M if f is a local diffeomorphism away from Σ, df drops rank by k at Σ, ker(df ) ⊂ T (Σ), and the determinant of the Jacobian matrix vanishes to order k at Σ.
involutive if the differentials dp i , i = 1, . . . , k, are linearly independant and the Poisson brackets satisfy {p i , p j } = 0, i = j.
From [5], we have the definition of the flowout.
Then the flowout of Γ is given by {(x, ξ; y, η) ∈ T * X × T * X : (x, ξ) ∈ Γ, (y, where H p i is the Hamiltonian vector field of p i . We now state a result of [3, Theorem 1.2] concerning the composition of FIO's with blowdown-blowdown singularities. ) be a canonical relation which satifies the following: (i) away from a hypersurface Σ ⊂ C, the left and right projections π L and π R are diffeomorphisms; (ii) at Σ, both π L and π R have blowdown singularities; (iii) π L (Σ) and π R (Σ) are nonradial and involutive.
If A ∈ I m (C) and B ∈ I m (C t ), then BA ∈ I m+m + k−1 where ∆ is the diagonal and Λ π R (Σ) is the flowout of of π R (Σ).
We now have our main Theorem.
Theorem 3. Let C be the canonical relation of the cylinder transform C. Then the left and right projections of C have singularities along a dimension 1 submanifold Σ, π L (Σ) and π R (Σ) are involutive and nonradial, and and Λ is the flowout of π R (Σ).
We can now calculate the determinant of Dπ L as follows: This is zero when x · θ = 0 or (x · θ α ) 2 + (x · θ β ) 2 cos 2 β = 0. The latter case corresponds to when x and θ are parallel, which we do not consider (x and θ are parallel only when the cylinder is degenerate, i.e. when s = 0), and so π L is a local diffeomorphism away from the manifold Σ = {x · θ = 0}.
Finally we show that the singularities of π L on Σ are blowdown of order 1. Indeed, on Σ we have and the kernel of Dπ L on Σ is So the left projection π L drops rank by 1 on Σ and its critical points on Σ are blowdown type singularities. Furthermore, which is involutive and nonradial (hereα andβ are the dual variables of α and β).
Using the same parameterization of C as above, the right projection is given by and its differential is .
(40) The determinant can now be calculated as Hence, on Σ So π R drops rank by 1 on Σ = {x · θ = 0} and its singularities are blowdown type as the kernel of Dπ R on Σ is Moreover, we have where The Hamiltonian vector fields of the p i are given by (45) Then, as x = tξ for some t ∈ R, we can see that ρ / ∈ span{H p 1 , H p 1 , H p 3 }, so π R (Σ) is nonradial.
To check that π R (Σ) is involutive, we first check that the Poisson brackets satisfy {p i , p j } = 0, i = j: Furthermore, if we work locally in a neighbourhood away from x 1 = 0, then p 1 , p 2 = 0 =⇒ p 3 = 0, so we need only consider the dependance of the differentials of p 1 and p 2 . dp 1 and dp 2 are linearly independant if and only if the Hamiltonian vector fields of p 1 and p 2 are linearly independant. But span{H p 1 , H p 2 } has dimension 2, so π R (Σ) is involutive. So the conditions of Theorem 2 are satisfied except for the fact that π L is two-to-one away from Σ. However, we can remedy this by working locally in neighbourhoods of any given point x 0 within B 3 α 1 ,α 2 small enough so that −x 0 is not in the same neighbourhood. When we compose the operators restricted to neighbourhoods of x 0 and −x 0 , we compose with the operator giving reflection in the origin so that even in that case we may apply Theorem 2.
So, applying Theorem 2 we have C * C ∈ I 2m+ k−1 The flow of H z is thus given by the system of linear ODE's with initial conditions (x(0), ξ(0)) = (x 0 , ξ 0 ). Then the solution to (51) at time t = 1 is and hence the flow of H z can be computed as the exponential of the matrix H. Now, if we parameterize z 1 = t cos ω and z 2 = t sin ω in terms of standard polar coordinates, then H = tG = t(ba T − ab T ), where a = (1, 0, 0) T and b = (0, cos ω, sin ω) T . Let P = −G 2 = aa T + bb T . Then P 2 = P (P is idempotent) and P G = GP = G.

From this it follows that
We can write e H = I 3×3 +G sin t+G 2 (1−cos t) =   cos t − sin t cos ω − sin t sin ω sin t cos ω sin 2 ω + cos t cos 2 ω − cos t cos ω sin ω sin t sin ω − cos t cos ω sin ω cos 2 ω + cos t sin 2 ω   , The right hand side of equation (55) is the standard parameterization of S 2 in terms of spherical coordinates, where t is the polar angle from the x axis pole and ω is the angle of rotation in the yz plane (azimuth angle). So e H defines a full set of rotations on S 2 and hence, given a vector-conormal vector pair (x, ξ), the Lagrangian Λ includes the rotation of x and ξ over the whole sphere. This completes the proof.
The above results tell us that the wavefront set of the kernel of the normal operator C * C is contained in ∆ ∪ ∆ ∪ Λ, where ∆ is the diagonal, ∆ the diagonal composed with reflection through the origin, and Λ is the flowout from π R (Σ), which is a rotation by Corollary 1. Also microlocally C * C ∈ I −2 (∆\Λ) and C * C ∈ I −2 (Λ\∆), which implies that the strength of the artefacts represented by Λ are the same as the image intensity on Σ = {x · θ = 0}. We will give examples of the artefacts implied by Λ later in our simulations in section 4, but in the next section we shall show how to reduce the strength of this artefact microlocally.

Reducing the strength of the image artifact
Here we derive a filter Q, which we show can be applied to reduce the strength of the image artefact Λ for the cylinder tranform C. We further show how Q can be applied as a spherical convolution with a distribution h on the sphere, which we will determine.
Using the ideas of [2], our aim is to apply a filtering operator Q : E ((0, 1) × S 2 ) → E ((0, 1) × S 2 ), whose principal symbol vanishes to some order s on π L (Σ), to C before applying the backprojection operator C * . From [2], we have the following theorem.
Theorem 4. Let A ∈ I m (C) be such that both the left projections of A, π L and π R are diffeomorphisms except on a set Σ where they drop rank by k, and let π L (Σ) and π R (Σ) be involutive and nonradial. Let Q be a pseudodifferential operator of order 0 whose principal symbol vanishes to order s on π L (Σ). Then A * QA ∈ I 2m+ k−1 where ∆ is the diagonal and Λ is the flowout from π R (Σ).
Let ∆ S 2 denote the Laplacian on S 2 , and I the identity operator. Then we will take whose symbol vanishes to order 2 on There are two technical issues with the application of Theorem 4 in our case. One is the fact, which we already mentioned in the proof of Theorem 3 that π L is two-to-one away from Σ. We can deal with this in the same way we dealt with in the proof of Theorem 4 by restricting C to small neighbourhoods of each point. The other issue is that this operator Q, defined by (56), is not a pseudodifferential operator on Y = (0, 1) × S 2 since differentiation of its symbol in the dual angular variables does not increase the decay in theŝ direction. However, this objection can be overcome by noting that Both Ψ 1 and Ψ 2 are then pseudodifferential operators, and Ψ 1 is elliptic and of order zero. Thus, if Ψ −1 1 is a pseudodifferential parametrix for Ψ 1 , we have Q = Ψ 2 Ψ −1 1 + R where R is an operator with smooth kernel and Ψ 2 Ψ −1 1 is a pseudodifferential operator satisfying the hypotheses of Theorem 4. From this the results of Theorem 4 hold when Q is given by (56).
We thus have, upon applying the filter Q to C before applying C * that C * QC ∈ L −4,2 (∆, Λ). So C * QC ∈ L −2 (∆\Λ) and C * QC ∈ L −4 (Λ\∆) and the strength of the artefact is reduced and is now less than the strength of the image.
We now show how the filter Q can be applied as a convolution with a distibution h on the sphere. First we give some definitions and theorems on spherical harmonic expansions. For integers l ≥ 0, |m| ≤ l, we define the spherical harmonics Y m l as Y m l (α, β) = (−1) m (2l + 1)(l − m)! 4π(l + m)! P m l (cos β)e imα , where and are Legendre polynomials of degree l. The spherical harmonics Y m l are the eigenfunctions of the Laplacian on S 2 , with corresponding eigenvalues c l = −l(l + 1). So ∆ S 2 Y m l = c l Y m l . From [8] we have the following theorem.
where dΩ is the surface measure on the sphere. Then the series converges uniformly absolutely on compact subsets of Z 3 to F .
So after writing the cylinder tranform C in terms of its spherical harmonic expansion, we can apply the filter Q as follows: where C lm = S 2 CȲ m l dΩ. For a function f on the sphere and h, a distribution on the sphere, we define the spherical convolution [9] where ω is the north pole, and from [9] we have the next theorem.
Theorem 6. For functions f, h ∈ L 2 (S 2 ), the harmonic components of the convolution is a pointwise product of the harmonic components of the transforms: For our case, this gives the following.
Proof. From equation (63) we have Defining h lm = h l for all l ∈ N, |m| ≤ l, we have by Theorem 6 where , which completes the proof.

Simulations
Given the equivalence of the spindle transform S and the cylinder transform C, and given also that the diffeomorphism defining their equivalence v(x) = 1 + 1 |x| 2 − 1 |x| · x |x| depends only on |x|, the artefacts described by the Lagrangian Λ apply also to the normal operator of the spindle transform S * S. Here we simulate the image artefacts produced by Λ in image reconstructions from spindle transform data and show how the filter we derived in section 3 can be used to reduce these artefacts. We also provide simulated reconstructions of densities which we should find difficult to reconstruct from a microlocal perspective (i.e. densities whose wavefront set is in directions normal to the surface of a sphere centred at the origin), and investigate the effects of applying the filter Q as a pre-conditioner, prior to implementing some discrete solver, in our reconstruction.
To conduct our simulations we consider the discrete form of the spindle transform as in [1], and solve the linear system of equations where A is the discrete operator of the spindle transform, x is the vector of pixel values and b = Ax is the vector of spindle transform values (simulated as an inverse crime).
To apply the filter Q derived in section 3, we decompose b into its first L spherical harmonic components and then multiply each component by the filter components h l /a l = l(l+1) l(l+1)−1 , 0 ≤ l ≤ L before recomposing the series.
Consider the small bead of constant density pictured in figure 2. In figure 3 we present a reconstruction of the small bead by unfiltered backprojection (represented as an MIP image to highlight the bead). Here we see artefacts described by the Lagrangian Λ as, in the reconstruction, the bead is smeared out over the sphere. If we apply the filter Q to the spherical components of the data and sum over the first L = 25 components before backprojecting, then we see a significant reduction in the strength of the artefact, the density is more concentrated around the small bead and the image is sharper. See figure 4. If we simply truncate the harmonic series of our data before backprojecting without a filter, this has a regularising effect and the level of blurring around the sphere is reduced. However we still see the artefacts due to Λ. See figure 5. The line profiles in figures 2-5 have been normalised. We note that in the reconstructions presented, the object is reflected through the origin in the reconstruction. This is as predicted by the Lagrangian ∆. The reflection artefact seems intuitive given the symmetries involved in our geometry. It was shown in [1] that the null space of S consists of odd functions (i.e. functions whose even harmonic components are zero), and so what we see in the reconstruction is the projection of the density onto its even components. We also see this effect in the reconstructions presented in [1]. Now let us consider the layered spherical shell segment phantom (the layers have values oscillating between 1 and 2) centred at the origin, shown in figure 6. We reconstruct the phantom by applying CGLS implicitly to the normal equations (i.e. we avoid a direct application of A T A) with 1% added Gaussian noise and regularise our solution using Tikhonov regularisation. See figure 8. Here the image quality is not clear and the layers seem to blur into one, and the jump discontinuities in the image are not reconstructed adequately. However if we arrange the layers as sections of planes and perform the same reconstruction (see figures 7 and 10), then the image quality is significantly improved, and the jump discontinuities between the oscillating layers are clear. This is as expected, as the wavefront set of the spherical density is contained in π R (Λ), so we see artefacts in the reconstruction. When the layers are arranged as planes this is not the case and we see an improvement in the reconstruction. In figure 9 we have investigated the effects of applying the filter Q as a pre-conditioner prior to a CGLS implementation. To obtain the reconstruction, we solved the system of equations Q 1 2 Ax = Q 1 2 b using CGLS with 1% added Gaussian noise. The filter has the effect of smoothing the radial singularities in the reconstruction. Here we see that the outer shell is better distinguished than before but the inner shells fail to reconstruct and overall the image quality is not good.
In figures 11, 12 and 13 we have presented reconstructions of the layered spherical shell and layered plane phantoms by Landweber iteration, with 1% Gaussian noise. Here the jump discontinuities in the spherical shell reconstruction are clearer. However as the Landweber method applies the normal operator (A T A) at each iteration, we see the artefacts predicted by Λ in the reconstruction and the spherical segment is blurred out over spheres centred at the origin. The artefacts are less prevalent in the plane phantom reconstruction. Although we do see some blurring at the plane edges. In figure 13 we have applied Q 1 2 as a pre-conditioner to a Landweber iteration. Here it is not clear that we see a reduction in the spherical artefact and there is a loss in clarity due to the level of smoothing.

Conclusions and further work
We have presented a microlocal analysis of the spindle transform introduced in [1]. An equivalence to a cylinder transform C was proven and the microlocal properties of C were studied. We showed that C was an FIO whose normal operator belonged to a class of distributions I p,l (∆, Λ), where Λ is the flowout from the right projection of C, which we calculated explicitly. In section 3, we showed how to reduce the size of the rotation artefact associated to Λ microlocally, through the application of an operator Q, and showed that Q could be applied as a spherical convolution with a distribution h on the sphere, or using spherical harmonics. We provided simulated reconstructions to show the artefacts produced by Λ, and showed how applying Q reduced the artefacts in the reconstruction. Reconstructions of densities of oscillating layers were provided using CGLS and a Landweber iteration. We also gave reconstructions of a spherical layered shell centred at the origin, using Q 1 2 as a pre-conditioner, prior to a CGLS and Landweber implementation and compared our results.
In future work we aim to derive an inversion formula of either a filtered backprojection or backprojection filter type. That is, we aim to determine whether there exists an operator A such that either S * A or AS * is a left inverse for S. After which we could see how the filter Q derived here may be involved in the inversion process. We also aim to assess if Sobolev space estimates can be derived for the spindle transform to gain a further understanding of its stability.