H Boffy

M C Baietto

P Sainsot

A A Lubrecht

Detailed modelling of a moving heat source using multigrid methods

The actual contact between solid surfaces is generaly rough and time dependent. This roughness can only be correctly described using a 3d model. Moreover, severe thermo-mechanical loading induces high stress and temperature gradients within a thin subsurface layer. Even higher gradients are obtained in the case of surface coatings. Consequently, the problem is strongly multiscale: from contact to coating to roughness, the characteristic dimensions range from millimeter to nanometer.

A straightforward discretisation of this multiscale problem would exceed the memory and CPU capabilities of current computers. The aim of this work is to propose a more efficient numerical model able to deal with this multiscale problem: using 10 9 points and 10 3 timesteps. This paper traces the history of the numerical solutions of the heat equation from the pioneering work of Carslaw and Jaeger, to the current era. The proposed model is based on multigrid techniques within a finite difference frame work. Localised refinement is implemented to optimize memory and computing time costs. The numerical performance of the solver is presented through a comparison with analytical results using different types of boundary conditions. A multisource contact is studied as a first approximation to real asperity interaction.

Introduction

Stationary and moving heat source problems are frequently encountered in tribological applications [START_REF] Cowan | Frictional Heating Calculations[END_REF] (frictional heating) as well as in many manufacturing processes [START_REF] Komanduri | Thermal Modeling of the Metal Cuting Process, Part II -Temperature Rise Distribution Due to Frictional Heat Source at the Tool-Chip Interface[END_REF] (metal cuting, laser treatment). Moreover, these problems which are complex and multiscale (time dependent, moving component, heat source distribution [START_REF] Hou | General Solutions for Stationary/Moving Plane Heat Source Problems in Manufacturing and Tribology[END_REF], dimensions ranging over several decades and boundary conditions) are at the root of different types of mechanical failure [START_REF] Kennedy | Thermal and Thermomechanical Effects in Dry Sliding[END_REF].

Many researchers have studied frictional heating. Some have analyzed a semi-infinite thermally isolated solid subjected to a heat source [START_REF] Jaeger | Moving sources of heat and the temperature at sliding contacts[END_REF][6] [START_REF] Hocine | Thermal behaviour of a rotating disk subjected to discrete surface heat sources[END_REF], or to several sources [START_REF] Baïri | Effect of thickness and thermal properties of film on the thermal behavior of moving rough interfaces[END_REF]. More recently, a relativistic moving heat source has been studied [START_REF] Ali | Relativistic moving heat source[END_REF].

This paper proposes a multigrid solution method to deal with the transient thermal problem. The proposed method is fast and particularly flexible compared with other methods (FFT,FEM). The direct solution of the heat equation allows one to deal with any kind of source displacement and speed. Multigrid methods enable general source geometries as their efficiency permits thousands of points in the contact area. Moreover, these methods allow efficient solutions, especially when including a local refinement strategy. Good convergence rates are obtained using single-grid calculations when considering Dirichlet boundary conditions. Local mesh refinement can be used [START_REF] Bai | Local Mesh refinement MultiLevel Techniques[END_REF]. However, very fine grids are often needed to capture thin boundary layer properties and to detect heat transfer spots. A local refinement strategy allows one to reduce the calculation time and use an adaptive mesh size where gradients are largest. The influence of the Peclet number has also been analyzed [START_REF] Rabi | The effects of Peclet number and cycling strategy on multigrid numerical solutions of convective-conductive problems[END_REF]. The differential equation of heat conduction for a stationary, homogeneous, isotropic semi-infinite solid (cf fig1) with constant thermal properties and no heat generation reads (Fourier conduction equation):

∂ 2 T ∂x 2 + ∂ 2 T ∂y 2 + ∂ 2 T ∂z 2 = 1 D ∂T ∂t ( 1
)
where D is the thermal diffusivity.

Using the following dimensionless parameters:

X = x a ; Y = y a ; Z = z a ; t = t τ ; T = T ∆T ; τ = a 2 D ( 2
)
Where a is the characteristic dimension in x, y and z. τ is the characteristic time. ∆T corresponds to the temperature increase compared with room temperature. The reduced form of equation (1) becomes:

∂ 2 T ∂X 2 + ∂ 2 T ∂Y 2 + ∂ 2 T ∂Z 2 = ∂T ∂ t (3)
The thermal boundary conditions (BC) are either: a prescribed temperature derivative (Neumann BC) or an imposed temperature value (Dirichlet BC). On the surface boundary (z = 0) the heat flux is prescribed:

q(x, y, t) = -k ∂T ∂z | z=0 (4)
Using the following dimensionless parameters:

Q = q q max = - k q max ∆T a ∂ T ∂ Z (5) ∆T = aq max k (6)
One obtains:

Q = - ∂ T ∂ Z (7)
On all other boundaries the temperature is imposed.

Discretization

The partial differential equation PDE is discretised using finite difference techniques. This discrete approximation introduces an error with respect to the solution of the PDE. A discretised first order accurate approximation for the Poisson equation and the time derivate is:

Ti-1,j,k -2 Ti,j,k + Ti+1,j,k h 2 x + Ti,j-1,k -2 Ti,j,k + Ti,j+1,k h 2 y + Ti,j,k-1 -2 Ti,j,k + Ti,j,k+1 h 2 z = Ti,j,k -T O i,j,k
h t (8) h x , h y and h z are the mesh-size in x, y and z-direction, respectively, Ti,j,k is the discrete approximation of the unknown T (x i , y j , z k ) to be solved and T O i,j,k is the discrete representation of T (x i , y j , z k ) at time t-1. To simplify the notation a stencil is introduced. Using this stencil L equation ( 8) is written as:

L < T > i,j,k = -T O i,j,k h t ( 9
)
with : where the gridsize is assumed to be equal in each direction (h x = h y = h z = h). At the surface, both Dirichlet and Neumann conditions can be implemented. In the case of Neumann BC, one uses a second-order discretization leading to: Ti,j,k-1 = Ti,j,k+1 -2hQ [START_REF] Bai | Local Mesh refinement MultiLevel Techniques[END_REF] The unknown Ti,j,k-1 is then substituted in the discretized heat equation at the surface (z = 0, k = 0).

Local MultiGrid solver

The heat problem is solved as a linear system of the form:

L T = f ( 11
)
Where L is a differential operator, T the exact solution of the problem and f the right hand side of the equation. Multigrid methods, employed to solve this problem, have been pioneered by Brandt [START_REF] Brandt | Multigrid techniques : 1984 guide with applications to fluid dynamics[END_REF]. Further details concerning the implementation can be found in Venner-Lubrecht [START_REF] Venner | Multilevel methods in lubrication[END_REF].

Boundary conditions

Solving the thermal problem considering Neumann BC adds several difficulties. In order to obtain a solver with the same efficiency as the one employing Dirichlet BC, the BC have to fulfill two constraints: consistency on various grids, no local error introduction.

Venner et al. [START_REF] Venner | Multilevel methods in lubrication[END_REF] relax the interior problem on the surface using the Neumann BC to fulfill both criteria. fig. [START_REF] Hou | General Solutions for Stationary/Moving Plane Heat Source Problems in Manufacturing and Tribology[END_REF] shows the results obtained: the convergence speed is not affected by the type of boundary condition but the residual is higher for the Neumann case because the solver does not start with the solution on the surface. Some local prerelaxations are necessary to decrease the residual near the boundaries before solving the problem on the entire grid.

Mesh refinement strategy

The solution of a fast moving heat source requires a multiscale approach, both in space and time. Dimensions ranging from microns to meters are concerned: a near subsurface layer (µm), a loaded zone (mm), specific volume dimensions (x, y, z) (cm,m). Furthermore, depending on the Peclet number, the time step can range from microseconds to seconds. A very fine grid is required around the heat source to capture the large gradients. Elsewhere such level of refinement is not required, especially when considering a semiinfinite solid. MG methods are very convenient for local mesh refinement. The finest grids can be restricted to smaller and smaller subdomains, whereas the coarse grids cover the entire domain (fig. (4b)). The fine grid residuals and the coarse grid correction are only carried out in the local part where the finest grid exists. The unit disk flow problem is solved using 6 levels, involving a coarsest grid with 16*16*8 and a finest grid of 512*512*256 points. The computational domain is

X ∈ [-8, 8], Y ∈ [-8, 8], Z ∈ [-8, 0
]. Local refinement does not affect the convergence speed as shown in Fig. 4(a). The convergence rate in both cases is approximately the same.

Validation

Green's function

To validate the numerical model, cases with analytical solutions have been computed. The general solution for the heat equation with a moving source in X direction (P e = U s τ /a) is given by Jaeger [START_REF] Jaeger | Moving sources of heat and temperature at sliding contacts[END_REF] :

T (X, Y, Z, t) = t 0 2 [4π( t -t′ )] 3/2 Q(X ′ , Y ′ , t′ )exp - (X -X ′ -[ t t′ P e(τ )dτ ]) 2 + (Y -Y ′ ) 2 + Z 2 4( t′ -t′ ) dX ′ dY ′ d t′ (12)
For a stationary flux one obtains :

T (X, Y, Z, t) = Q(X, Y ) * t 0 G(X, Y, Z, t -t′ )d t′ (13)
Liu, Lannou et al. [START_REF] Liu | Solutions for Temperature Rise in Stationnary/Moving Bodies Caused by Surface Heating With Surface Convection[END_REF] give an algebraic solution for the time integrated Green's function for a constant speed (P e = constant) :

t 0 G(X, Y, Z, t-t′ ) dt ′ = G s (X, Y, Z, t) = e -P e(R-X)/2 4πR 2 -erf c P e √ t 2 - R 2 √ t + e (R * P e) erf c P e √ t 2 + R 2 √ t (14) With R 2 = X 2 + Y 2 + Z 2

Stationary case

Analytical solution

For a stationary source and established solution, i.e. P e = 0 and t = ∞, one obtains:

∞ 0 G(X, Y, Z, t -t′ ) dt ′ = G P e (X, Y, Z) = 1 2πR (15)
The same Green's function is found for a normal displacement due to a normal strain in solid mechanics. Four analytical cases are presented: constant flux on a square surface, constant flux on a disk surface, inverse elliptical flux on a disk surface, hertzian flux on a disk surface.

• For a constant flux applied on the surface of a unit square, i.e. Q(X, Y, Z = 0, t) = 1, one obtains:

T (X, Y, Z = 0) = 1 2π X i * log(Y i + X i 2 + Y i 2 ) + Y i * log(X i + X i 2 + Y i 2 ) X+dx/2 X-dx/2 Y +dy/2 Y -dy/2 (16)
T (0, 0, 0) = 2 * ln(1 + √ 2)/π This solution is shown in fig. 5.(a).

• For Q(X, Y, Z = 0) = 1 on a disk with a radius R ay centered on O. One obtains :

if

R < R ay T (X, Y, Z = 0) = R ay * 2 π * E R R ay (17) R ay = 1; T (0, 0, 0) = 1; T (R = R ay , 0) = 2 π if R > R ay T (X, Y, Z = 0) = R ay * 2 π * R E R R ay -1 - R 2 ay R 2 * K R ay R ( 18
)
With E the Young modulus and

K = E/2(1 + ν), ν the poisson coefficient. The solution is shown in fig.5.(b). • For Q(X,Y) = 1/ 1 -[R/R ay ] 2 on a disk with a radius R ay centered on O. Calculations give: if R < R ay T (X, Y, Z = 0) = R ay * π 2 T (R = 0, 0) = R ay * π 2 ; T (R = R ay , 0) = R ay * π 2 (19) if R > R ay T (X, Y, Z = 0) = R ay * arcsin( R ay R ) ( 20
)
The solution is shown in fig. 5.(c).

• For Q(X,Y) = 1 -[R/R ay ] 2 on a disk with a radius R ay centered on O. Calculations give: if R < R ay T (X, Y, Z = 0) = R ay * π 8 (2 - R 2 R 2 ay ) (21)
T (R = 0, 0) = R ay * π 4 ; T (R = R ay , 0) = R ay * π 8 if R > R ay T (X, Y, Z = 0) = 1 4 R ay 2 - R 2 R 2 ay * arcsin R ay R + R R ay 1 - R 2 ay R 2 ( 22
)
The solution is shown in fig. (5d). Calculation have been done on a cubic volume which is not a semiinfinite one. Consequently, Dirichlet boundary condition are used on the borders of the volume as given in eq [START_REF] Liu | Solutions for Temperature Rise in Stationnary/Moving Bodies Caused by Surface Heating With Surface Convection[END_REF] excepted on the surface where a flux is applied. To approximate at best the steady-state results, convergence has been pursued until residuals were less than 10 -8 . The temperature distributions of equations ( 16) to [START_REF] Archard | The temperature of Rubbing Surfaces[END_REF] are plotted in fig. [START_REF] Jaeger | Moving sources of heat and the temperature at sliding contacts[END_REF]. These thermal problems are numerically solved using the MG technique presented in this paper. In further sections the numerical solutions will be compared with the analytical solutions. The dimensions of the domain with respect to the heat source, the mesh size and consequently the numerical cost are analyzed through a detailed comparison. The discretization error is studied through a convergence analysis. The efficiency of the local refinement strategy is demonstrated through the evolution of the CPU time.

Numerical solution

Effect of boundary condition

The problem analysed is a unit disk flow applied on the surface [START_REF] Gao | An FFT-based transient flash temperature model for general three-dimensional rough surface contact[END_REF]. The coarsest grid is discretized using 16*16*8 points. In order to approximate the semi-infinite conditions, boundary conditions have been implemented on the borders according to the far field solution [START_REF] Liu | Solutions for Temperature Rise in Stationnary/Moving Bodies Caused by Surface Heating With Surface Convection[END_REF]. Fig. [START_REF] Laraqi | An exact explicit analytical solution for steady-state temperature of a half space subjected to a circular moving heat source[END_REF] shows the influence of the domain size on the central temperature for a constant mesh size. The error in the temperature field is caused by the finite dimensions of the volume studied. The order of magnitude of this error has to be evaluated to ensure it does not dominate the overall solution. The next subsection shows how one can increase accuracy by increasing the domain size.

Effect of mesh size and CPU time

In this subsection the domain size and the number of grid points in the source are analysed. Two aspects need to be studied, the first one concerns the decrease of the numerical error with the mesh size. The 7) as a function of the source discretization. One can remark a similar error reduction as in fig. [START_REF] Laraqi | An exact explicit analytical solution for steady-state temperature of a half space subjected to a circular moving heat source[END_REF]. Furthermore, the error stagnates for a certain source discretization for each domain size. Increasing the domain allows a further error reduction. This figure also shows the efficiency of the local refinement strategy, which does not significantly affect the precision. Fig. [START_REF] Baïri | Effect of thickness and thermal properties of film on the thermal behavior of moving rough interfaces[END_REF] gives the CPU time using global and local calculations, with a residual less than 10 -12 . One can conclude that the problem can be solved with almost the same accuracy but a 100 times smaller computational effort using local grid refinement.

Moving heat source

Analytical solution

A moving heat source problem can be characterized using the Peclet number, defined by U s = τ V /a. This number is directly linked to the mesh size (a is the characteristic dimension in x and y direction). A non-dimensional unit velocity corresponds to a single mesh size displacement per time step. A well-known 2D solution describes a moving heat band on the surface of a semi-infinite volume.

The solution has been calculated by Jaeger [START_REF] Johnson | Contact Mechanics[END_REF]. For large values of Pe (>10), the heat will diffuse only a short distance into the solid in the time taken for the surface to move through the heated zone. The heat flow will then be approximately perpendicular to the surface at all points. The temperature of a surface point is given by Carslaw and Jaeger [START_REF] Jaeger | Moving sources of heat and temperature at sliding contacts[END_REF]. (1 + X)

1/2 (23)
T √ P e = 4 π (1 + X) 1/2 (24) With -1 ≤ X ≤ 1
The mean temperature for a band source is approximately:

Tmean ≈ 4Q 3π 1/2 P e -1/2
(25)

Numerical solution

In this subsection, the problem of the infinitely wide moving source is studied. The source moves in the x direction. The comparison between analytical ( 23) and numerical solution ( T * √ P e) is plotted in fig. [START_REF] Ali | Relativistic moving heat source[END_REF]. The calculations use 3 global and 4 local levels in y and z direction. The band source is discretized using 129 points in the x-direction.

Fig. (9.a) shows that for increasing Peclet number, the numerical solution tends to the analytical one. Fig. (9.b) shows the evolution of the error between both solutions. This evolution is approximately first order as can be expected from [START_REF] Baïri | Effect of thickness and thermal properties of film on the thermal behavior of moving rough interfaces[END_REF]. At very low speeds (Pe<0,1) the temperature distribution becomes symmetrical and similar to that of a stationary source.

General sources

In this section, a moving multisource problem is studied to approximate a rough thermal contact. Several circular fluxes are applied on the surface in a circular area, as shown in fig [START_REF] Bai | Local Mesh refinement MultiLevel Techniques[END_REF]. The other boundaries are kept at a prescribed temperature. Gao et al. [START_REF] Gao | An FFT-based transient flash temperature model for general three-dimensional rough surface contact[END_REF] and Yu et al. [START_REF] Yu | FFT Thermoelastic solutions for moving heat source[END_REF] studied this kind of problem using the FFT method. The 7 heat spots move with a velocity defined by the Peclet number. In a first subsection the influence of the size of the spots is studied, keeping the Peclet number constant. In a second subsection the effect of speed is studied, looking at several Peclet numbers and keeping the spot size constant.

Speed effect

The volume studied has dimensions of 64*32*16, the coarsest grid is discretized using 64*32*16 points. The local refinement strategy is used: the calculations have been performed using 6 or 7 levels (for high Peclet numbers (> 10)): a global one and 5 or 6 local ones. Each small disk is discretized using approximately 30 points on the finest grid. Several Peclet numbers have been used: 1, 2, 10 and 100. In order to visualize the influence of the spot size, 3 diameters have been studied: 1/3, 1/4 and 1/8. When using a radius of 1/3, the disks are in contact. The overall intensity in each unitary disk is kept constant (π). The fig. [START_REF] Rabi | The effects of Peclet number and cycling strategy on multigrid numerical solutions of convective-conductive problems[END_REF] shows that the global temperature profile (outside the source) is not affected by the source micro-geometry. The following study will consequently look after the temperature distribution inside the source. The surface temperature is plotted in the x direction for each Peclet number. These curves are compared with the smooth unit disk result. The results are shown in fig. [START_REF] Johnson | Contact Mechanics[END_REF]:

One finds the usual results of a moving heat source. The curves show that each disk flow is little affected by the other ones for a radius less than 1/4. The maximum flash temperature rise [START_REF] Blok | The flash temperature concept[END_REF] is located towards the trailing regions of the contact depending of the Peclet number. In addition, the maximum flash temperature decreases when the Peclet number increases. For high Peclet numbers (> 10), the temperature profile in each disk becomes more and more pointed. Flash temperatures are expressed in terms of the heat supply over the contact area, the velocity, and the thermal properties of the material. These have been developed by Blok and Jaeger for linear and square contacts [START_REF] Blok | Theoretical Study of the Temperature Rise of Surfaces of Actual Contact Under Oiliness Lubricating Conditions, General Discussion on Lubrication[END_REF][START_REF] Jaeger | Moving sources of heat and temperature at sliding contacts[END_REF] and by Archard [START_REF] Archard | The temperature of Rubbing Surfaces[END_REF] for circular contacts.

Size effect

In this subsection, the effect of the radius variation is analyzed. The same curves as in fig. 12. are presented, showing the same radius in each sub-figure for different Peclet numbers, fig. 13.

One remarks that the curves with a radius of 1/3 look like the ones obtained for the semi infinite band case. When the radius decreases, the flash temperatures are better distinguished. The profiles obtained can be considered as a first approximation of a rough contact between two solids.

Sub-layer analysis

In this subsection, the temperature field below the heat source is analyzed. The objective is to see how velocity affects the temperature in the sublayer. Fig. 141516. shows the temperature fields for the same Peclet numbers and radii as in fig. 11 and12.

The temperature field for a Peclet number less than 0.2 (not shown) closely resembles the stationary case [START_REF] Johnson | Contact Mechanics[END_REF]. For medium Peclet numbers (0.2 to 10), the temperature distribution remains largely symmetrical. The temperature influence occurs closer to the surface when speed increases [START_REF] Laraqi | Thermal Constriction Phenomenon in Sliding Contact[END_REF]. the macroscopic aspect was demonstrated through the influence of the domain size. The borders of the volume have to be distant enough from the source to respect the semi-infinite conditions. Concerning the microscopic aspects, it has been shown that the refinement strategy allows a good accuracy and permits one to observe microscopic phenomena close to the source. Using the present code, one can obtain an accuracy better than 10 -2 as shown by the comparison between computed and analytical solutions. The algorithm reaches this accuracy because the hypothesis of a semi-infinite volume has been closely respected using a large enough volume. Moreover, calculations can be performed using 10 9 points on the finest grid. Furthermore stationary and moving sources, single and multiple sources have been studied. The increase in velocity creates a reduction in the temperature amplitude and a displacement of the temperature maximum in the direction of the movement. Our major objective is to extend this work towards 3D thermoelastic calculations.

Conclusion and further research

1 2Figure 1 :

 11 Figure 1: Description of the physical domain and coordinate system.

Figure 2 :

 2 Figure 2: Heat transfer stencil.

Figure 3 :

 3 Figure 3: Residual reduction on the finest level as a function of the number of V-cycles for Dirichlet and Neumann BC.

Figure 4 :

 4 Figure 4: (a) Global and Local residual reduction obtained for the case of a unit disk flow (b). Local 2-D grids

Figure 5 :

 5 Figure 5: Numerical temperature distribution for a unit square flux, unit disk flux, disk flux with Q(X,Y) = 1/ 1 -[R/R ay ] 2 and disk flux with Q(X,Y) = 1 -[R/R ay ] 2 left to right, top to bottom

Figure 6 :

 6 Figure 6: Effect of domain size on T (0, 0, 0) for a unitary disk flow, for a constant discretization error: 65 points in the source (finest grid)

Figure 7 :

 7 Figure 7: Relative error between analytical and numerical solutions as a function of the source discretization, circular unit flux, 16*16*8 points on the coarse grid. Volume size increase by a factor of 2 from top (4*4*2) to bottom (64*64*32). Global and local refinement results

Figure 8 :

 8 Figure 8: CPU time (in seconds).

Figure 9 :Figure 10 :

 910 Figure 9: (a)Comparison between analytical and numerical solutions for increasing Peclet number. (b)Error between analytical and numerical solution for Pe=1, 5, 10, 25, 50 and 100

Figure 11 :

 11 Figure 11: T(X,0,0) for a Peclet number of 1. Solutions for r=1/8,1/4,1/3 and a unit disk are plotted top to bottom.

Figure 12 :

 12 Figure 12: T(X,0,0) for different Peclet numbers with r=1/8,1/4,1/3 and unit disk solution top to bottom.

Figure 13 :

 13 Figure 13: T (X,0,0)* √ P e for different source radii with Pe=100, 10, 2, 1 top to bottom

Figure 14 :

 14 Figure 14: T (X, Y = 0, Z) for a Peclet number of 1, r=1/3, 1/4, 1/8

A 3 -Figure 15 : 8 Figure 16 :

 315816 Figure 15: T (X, Y = 0, Z) for a Peclet number of 10, r=1/3 , 1/4 , 1/8


Rémi Régnier

Mai K Nguyen

A novel coupled transmission-reflection tomography and the V-line Radon transform

Keywords: X-ray imaging, tomography, reflection, mirrors, Radon transform

The Radon transform (RT ) on straight lines deals as mathematical foundation for conventional tomographic imaging (e.g. X-ray scanner, Single Photon Emission Computed Tomography (SPECT), Positron Emission Tomography (PET)) which use only one physical phenomenon, i.e. either transmission or emission of radiation. An imaging concept exploiting jointly two phenomena leads to a Radon transform defined on a geometrical support different from straight lines. In this paper, we propose a new two-dimensional X-ray imaging based on the coupling between transmission and reflection. Its modeling leads to a Radon transform defined on a pair of half-lines forming a vertical letter V, called V-line RT (V -RT ). Moreover we establish the analytic inverse formula of this new transform, which forms the mathematical basis for image reconstruction. Through simulations, image formation and reconstruction results show the feasibility of this new imaging. The main advantage is the use of an one-dimensional detector which does not rotate for two-dimensional image reconstruction.

I. INTRODUCTION

While the classical Radon transform defined on straight rays models the imaging systems (X-ray scanner, SPECT) which need a relative rotation between the object and the detector ( [START_REF] Radon | Über die Bestimmung von Funktionnen durch ihre Integralwerte längs gewisser Mannigfaltikeiten[END_REF], [START_REF] Natterer | The Mathematics of Computerized Tomography[END_REF], [START_REF] Barrett | The Radon Transform and its Applications[END_REF]), there exists a considerable interest for systems which do not require such a relative rotation and a two-dimensional detector for two-dimensional image reconstruction. The modeling of the last systems leads to the RT defined on other geometries than the straight-lines.

In this paper, we propose a new tomographic modality which operate with X-ray. But instead of simply analyzing the absorption coefficients on a straight line as in X-ray scanner, we place a special mirror capable of reflecting radiation and so we project the absorption function on a pair of half-lines (and no longer on a straight line). This variant of the tomography would allow to replace the rotation of source-detector couple by setting a series of point-sources and a series of point-detectors. So using a one dimensional non rotating detector for two-dimensional image reconstruction turns out to be an attractive alternative to conventional tomography.

Mathematically, the modeling of this coupled transmissionreflection tomography leads to a Radon transform defined Rémi Régnier and Mai K. Nguyen are with the Laboratory Equipe Traitement de l'Information et Systèmes, CNRS UMR 8051/ENSEA/Université de Cergy-Pontoise, mail: (see regnierremi@yahoo.fr, mai.nguyen-verger@ucergy.fr). on a pair of half-lines forming a V-letter, called V-line RT (V -RT ). This transform is a member of the family of the so-called T V -transform ( [START_REF] Nguyen | Imagerie par rayonnement gamma diffus[END_REF], [START_REF] Morvidone | On the V-line Radon transform and its imaging applications[END_REF], [START_REF] Nguyen | Compton Scattering Emission Imaging Based on the V line Radon transform and its medical imaging applications[END_REF], [START_REF] Régnier | Compton Scattering Emission Imaging Based on the V line Radon transform and its applications[END_REF], [START_REF] Nguyen | Inversion of a new V-line Radon transform and its numerical analysis[END_REF]). The simple T V -transform was first proposed by Basko [START_REF] Basko | Analytical reconstruction formula for the one-dimensional Compton camera[END_REF] in 1997 as a model for image formation in a two-dimensional Compton camera. However this Basko transform is in fact a V-line Radon transform with swinging axis around a detection site whereas the one considered here has a fixed axis direction [START_REF] Nguyen | Imagerie par rayonnement gamma diffus[END_REF]. The V-line Radon transform is also considered as a two-dimensional version of the Conical Radon Transform (CRT ) which have been introduced some years ago ( [START_REF] Morvidone | On the V-line Radon transform and its imaging applications[END_REF], [START_REF] Nguyen | On an integral transform and its inverse in nuclear imaging[END_REF], [START_REF] Nguyen | A novel inverse problem in γ-rays emission imaging[END_REF], [START_REF] Truong | The mathematical foundation of 3D Compton scatter emission imaging[END_REF], [START_REF] Nguyen | On a novel approach to Compton scattered emission imaging[END_REF]). These generalized Radon transforms (T V and CRT ) have numerous applications in imaging science, in particular in new Compton scattering emission imaging ( [START_REF] Nguyen | Imagerie par rayonnement gamma diffus[END_REF], [START_REF] Morvidone | On the V-line Radon transform and its imaging applications[END_REF], [START_REF] Nguyen | Compton Scattering Emission Imaging Based on the V line Radon transform and its medical imaging applications[END_REF], [START_REF] Régnier | Compton Scattering Emission Imaging Based on the V line Radon transform and its applications[END_REF], [START_REF] Nguyen | Inversion of a new V-line Radon transform and its numerical analysis[END_REF], [START_REF] Nguyen | On an integral transform and its inverse in nuclear imaging[END_REF], [START_REF] Nguyen | A novel inverse problem in γ-rays emission imaging[END_REF], [START_REF] Truong | The mathematical foundation of 3D Compton scatter emission imaging[END_REF], [START_REF] Nguyen | On a novel approach to Compton scattered emission imaging[END_REF]).

In this paper a new application of the T V -transform is developed in primary radiation transmission imaging, more precicely as a coupled transmission-reflection tomography.

In this new tomography, X-radiation emitted from the source falls on the mirror under an incidence angle and registered later by a detector (see Fig. 1). In this way when the object is scanned by all possible V-lines, then the standard rotational motion can be avoided. In section 2 we present the modeling of the new transmission-reflection imaging concept which leads to the V-line Radon transform. Then we establish its analytic inverse and derive a corresponding filtered back-projection form. This last form has the advantage of reconstructing the image by fast algorithms.

In section 3, we present numerical simulation results on image formation and reconstruction to support the feasibility of this new imaging. The paper ends with a short conclusion on the obtained results and opens some future research perspectives.

II. MODELING OF THE TRANSMISSION-REFLEXION IMAGING AND THE V-LINE RADON TRANSFORM

A. Image formation modeling and the V-line Radon transform (V -RT )

In tomographic transmission imaging, the attenuation map of a two-dimensional object is represented by a non negative continuous function f (x, y) with bounded support. A linear detector collects the transmitted radiation after a reflection on a mirror. We may assume that the reflection coefficient of the mirror is of the form e -c(θ) , where θ is the angle of incidence. From now on, we consider the ideal case where c(θ) = 0. This does not affect the validity of the modeling.

Each measurement on the detector represents the sum of the radiation attenuation along a trajectory in the shape of letter V, since photons travel from the point of emission at the source up to the point of detection after a reflection on the mirror. Thus f (x, y) is integrated on a discontinuous line having the form of the V-line with a symmetry axis parallel to the Oydirection (see Fig. 2). The projection data g(x R , θ) on the detector can be written as

Mirror

g(x R , θ) = ∞ 0 f (x R + r sin(θ), r cos(θ))dr + ∞ 0 f (x R -r sin(θ), r cos(θ))dr. (1)
Equation (1) gives g(x R , θ), as the V-line Radon transform of the unknown attenuation distribution f (x, y). This is the image formation equation in the transmission-reflection imaging procedure.

B. The inverse transformation (V -RT ) -1 and image reconstruction

For convenience of notations, let us define g(x R , t) = g(x R , θ), with t = tan θ. The inverse transform (V -RT ) -1 can be worked out using Fourier transforms f (q, y) (resp. g(q, t)) with respect to the variable x (resp.

(x R )) of f (x, y) (resp. g(x R , t)), g(x R , t) = ∞ -∞ g(q, t) e 2iπqx R dq (2) f (x, y) = ∞ -∞ f (q, y) e 2iπqx dq.
Equation ( 1) becomes now with the change of variable z = r cos(θ) g(q, t)

√ 1 + t 2 = ∞ 0 f (q, z) 2 cos(2πqzt) dz for θ ∈ [0, π/2] (3) g(q, -t) √ 1 + t 2 = ∞ 0 f (q, -z) 2 cos(2πqzt) dz for θ ∈ [-π/2, 0] . (4)
We multiply equations (4) on both sides by ∞ 0 2 cos(2πzt) d (qt) and using the identity

∞ 0 4 cos(2πηx) cos(2πη x)dx = δ(η + η ) + δ(η -η ), (5) we obtain f (q, z) = |q| ∞ 0 g(q, t) √ 1 + t 2 2 cos(2πqtz)dt, (6)
f (q, -z) = |q| ∞ 0 g(q, -t) √ 1 + t 2 2 cos(2πqtz)dt. (7)
So for the full angular range, we can write

f (q, z) = |q| ∞ 0 2 cos(2πqtz).U (z)g(q, t) dt √ 1 + t 2 + |q| ∞ 0 2 cos(2πqtz).U (-z)g(q, -t) dt √ 1 + t 2 , ( 8
)
where U is the Heaviside unit step function. By inverse Fourier transform we recover

f (x, z) = ∞ -∞ dq e 2πiqx |q| ∞ 0 2 cos(2πqtz).U (z)g(q, t) dt √ 1 + t 2 + ∞ 0 2 cos(2πqtz).U (-z)g(q, -t) dt √ 1 + t 2 . ( 9
)
This formula is also a so-called filtered back-projection method for image reconstruction [START_REF] Morvidone | On the V-line Radon transform and its imaging applications[END_REF]. It has the advantage of offering fast algorithms. Technically the back-projection operation consists in assigning the value g(x R , t) = g(x R , θ) to every point on the "projection" V-line, which has given rise to this value, and then to sum over all contributions for every V-line "projections". Therefore we shall follow this procedure in our numerical simulations in the next section.

III. NUMERICAL SIMULATIONS In this section, we simulate image formation and reconstruction using the V -RT and its inverse for a point object and for a Shepp-Logan phantom (Fig. 3,6)). The original images are of size 128×128 of length units and the angular sampling rate is dθ = 0.005 rad. Fig. 4 shows the V -RT transform of the point object. It is also called the point spread function (PSF) of the transform. We note that the shape of the V-line Radon transform (or the shape of the projections) in the (θ, x S )-representation (x S is the source position) is consistent with the arctan-curve [START_REF] Régnier | Compton Scattering Emission Imaging Based on the V line Radon transform and its applications[END_REF]. The reconstruction of the point object is presented in Fig. 5.

The point object is well reconstructed but there are some artifacts. The use of back-projection on V-lines generates more artifacts than back-projection on straight lines in classical Radon transform. This is due to the existence of more spurious line intersections. In order to reduce these artifacts, the Hann filter H is used. It is defined on the Fourier domain by its action on the first variable of a function f as follows :

Hf (q, y) = |q| 2 (1 + cos(2πq)) f (q, y), (10)
where the Fourier transform is taken on the first variable. Moreover for a better reconstruction, we would have to scan the angular range of reflection between 0 and π/2. This is not possible because we would have to take away the source from the detector over huge distances. We are thus limited to a restricted angular range between 0 and π/4 in our simulations (314 projections for every position of the source). This causes a loss of information resulting in incomplete data.

However this lack of information can be recovered by a second mirror opposite to the first one and two sets of data are used for the reconstruction (see Fig. 7, Fig. 8). We can observe an improvement on the reconstruction of the Shepp-Logan phantom from two sets of data (Fig. 9). There are less artifacts and the small structures in the object are clearly reconstructed.

These results illustrate undoubtedly the feasibility of the new imaging modality.

IV. CONCLUSION AND PERSPECTIVES

In this paper, we propose a novel coupled transmissionreflection tomography based on the V-line Radon transform, which is a generalization of the classical Radon transform and can be used in biomedical imaging, in material nondestructive testing or in homeland security applications. The corresponding back-projection inversion formula generalizes the one of the ordinary Radon transform. In this new imaging the rotational motion between object and detector is replaced by the translational motion which is less stringent in some working conditions, for example for long objects or a series of object on a conveyor belt. The numerical simulation results support the feasibility and the performance of this new imaging.

The V-line Radon transform with its family of associated TV-transforms are instrumental in new possible tomographic processes involving scattering, emission and coupled transmission-reflection phenomena. More general Radon transforms on piecewise continuous curves consisting of connected pieces of arcs, e.g. arcs of circle [START_REF] Cormack | The Radon transform on a family of curves in the plane[END_REF], [START_REF] Quinto | Radon transforms on curves in the plane[END_REF], as well as Radon transforms on swinging V-line around a fixed point are interesting to be considered since their inversion would open the way to new advantageous imaging processes. These topics are planned for future research.

Fig. 1 .

 1 Fig.1. The X-ray source can be translated (see red rays) and emits according to some angle (see the blue rays) which generates enough data on the detector to reconstruct the attenuation map.

Fig. 2 .

 2 Fig. 2. Parameters of the V-line Radon transform

Fig. 3 .Fig. 5 .

 35 Fig. 3. Original point object

Fig. 6 .

 6 Fig. 6. Original Shepp-Logan phantom

Fig. 8 .

 8 Fig. 8. The V-line Radon transform of the Shepp-Logan phantom shown in Fig. 6 with dθ =0.005 rad for the second set of data.

Fig. 9 .

 9 Fig. 9. Reconstruction of the Shepp-Logan phantom from the two sets of data. The M SE = 1.15 * 10 -2 .