Frédéric Brigui
email: frederic.brigui@isae.fr.

Guillaume Ginolhac
email: guillaume.ginolhac@u-paris10.fr

Laetitia Thirion-Lefevre
email: laetitia.thirion@supelec.fr.

Philippe Forster
email: philippe.forster@u-paris10.fr.

New SAR Target Imaging Algorithm based on Oblique Projection for Clutter Reduction

Keywords: SAR algorithms, subspace methods, oblique projection, physical-based models, polarimetry, statistical performance, FoPen detection. I

We have developed a new Synthetic Aperture Radar (SAR) algorithm based on physical models for the detection of a Man-Made Target (MMT) embedded in strong clutter (trunks in a forest). The physical models for the MMT and the clutter are represented by low-rank subspaces and are based on scattering and polarimetric properties. Our SAR algorithm applies the oblique projection of the received signal along the clutter subspace onto the target subspace. We compute its statistical performance in terms of probabilities of detection and false alarms. The performances of the proposed SAR algorithm are improved compared to those obtained with existing SAR algorithms: the MMT detection is greatly improved and the clutter is rejected. We also studied the robustness of our new SAR algorithm to interference modeling errors. Results on real FoPen (Foliage Penetration) data showed the usefulness of this approach.

change detection (for instance, see [START_REF] Ulander | Change detection of vehicle-sized targets in forest concealment using VHF-and UHF-band SAR[END_REF]) that consists in detecting temporal changes by comparing two SAR images acquired over the same area, with the same frequency support. The clutter is assumed stationary and therefore is expected to be suppressed. The performance may be improved but they highly depend on the incidence angle and on the choice of test statistics. In addition, it is assumed that the target has moved between two acquisitions. Other techniques rely on the variation of the scattering properties of the target to detect with either this incidence angle (also named range angle) or the azimuthal angle. These methods, generally called subaperture decompositions provide subaperture datasets from SAR images using deconvolution. The distinction between targets and clutter relying on the aspect dependency of their scattering signatures has been presented in [START_REF] Runkle | Multi-aspect target detection for SAR imagery using hidden Markov models[END_REF]. It is to note that subaperture technique applied to target detection has been extended to polarimetry (for instance see [START_REF] Souyris | On the use of complex SAR image spectral analysis for target detection: Assessment of polarimetry[END_REF], [START_REF] Ferro-Famil | Scene characterization using subaperture polarimetric SAR data[END_REF]). However, several limitations have been listed in [START_REF] Runkle | Multi-aspect target detection for SAR imagery using hidden Markov models[END_REF], [START_REF] Souyris | On the use of complex SAR image spectral analysis for target detection: Assessment of polarimetry[END_REF] regarding in particular the robustness to variations in the target response.

Polarimetric interferometry has also been used to derive filter to suppress the foliage clutter, but not the target to detect [START_REF] Cloude | Target detection beneath foliage using polarimetric SAR interferometry[END_REF]. This technique has been tested in [START_REF] Colin | FOPEN with polarimetric interferometry: Validations with experimental data at P-band[END_REF] on real full polar data acquired by RAMSES ONERA sensor in P-band and in repeat-pass mode. The forest clutter is widely suppressed, but a lot of false alarms remain in this case that prevent from an efficient detection. All these techniques analyze data that have been already synthesized using classical SAR tools that assume that any target can be seen as a set of isotropic scatterers. So, even if these postprocessing methods try to compensate this assumption, a lot of information is missing. We proposed in previous works [START_REF] Durand | New SAR processor based on matched subspace detector[END_REF], [START_REF] Ginolhac | SAR processors based on signal or interference subspace detector[END_REF] new SAR algorithms that take into account MMTs and interference scattering properties and do not assume any more the isotropy of the scatterers. Moreover, we assume that the forest disturbance is a sum of deterministic clutter (interferences) caused by the trunks and random noise caused by the foliage, the branches of the trees ... Some solutions have been proposed to use the radar signature of the targets, and all of them postulate that MMTs are a set of canonical scatterers (as a set of plates) [START_REF] Tiberi | Electromagnetic scattering from large faceted conducting bodies by using analytically derived characteristic basis function[END_REF], [START_REF] Laviada | Solution of electrically large problems with multilevel characteristic basis function[END_REF]. In this case, the orientation of the canonical scatterers is important to determine, in addition to their location. Considering the directivity of a MMT scattering has consequently led to an increased number of unknowns in the detection problem. A simple approach is to develop a filter adapted to the orientations and sizes of the target to detect [START_REF] Allen | FOPEN-SAR detection by direct use of simple scattering physics[END_REF] (a dihedral corner reflector in the above quoted paper). This processing has several limitations: it is not robust to configurations not exactly described by the sampling, the computation time may be prohibitive and the interpretation may be dubious (as a target contains elements with various orientations). A similar approach has been developed in [START_REF] Sharma | Matched subspace detectors for discrimination of targets from trees in SAR imagery[END_REF], where the authors use a subspace generated from the signals scattered by a dihedral corner reflector. However, the dimension of this subspace is large, and therefore does not lead to a significant improvement of the detection performances. In addition, the generation of such a subspace [START_REF] Forster | Generalized rectification of cross spectral matrices for array of arbitrary geometry[END_REF] DRAFT requires a significant computation time, which prevents from any application to real data. Recently, methods based on sparse representation using structured dictionaries have been proposed to taking into account scattering features of targets in the SAR image formation process. Dictionary which takes into account anisotropic properties of the target is proposed in [START_REF] Varshney | Sparse representation in structured dictionaries with application to synthetic aperture radar[END_REF]; shape-based and wavelet-dictionaries are studied in [START_REF] Samadi | Sparse representation-based synthetic aperture radar imaging[END_REF]. Although these methods give good results for SAR image quality enhancement, they have not been applied to deterministic noise rejection for FoPen detection.

We have previously proposed to reconsider the SAR algorithms by including prior knowledge based on simple physical models of the MMT and the interferences. In the first algorithm [START_REF] Durand | New SAR processor based on matched subspace detector[END_REF], the Signal Subspace

Detector SAR (SSDSAR), a MMT is assumed to be a set of plates whose scattering belongs to a low-rank subspace, called target subspace. The SSDSAR consists in projecting orthogonally the SAR received signal into the target subspace. Compared to classical SAR processors, the SSDSAR algorithm has shown a 5dB detection gain for a single polarization [START_REF] Durand | New SAR processor based on matched subspace detector[END_REF]. To reduce false alarms due to the interferences, another SAR algorithm (SISDSAR for Signal and Interference Subspace Detector SAR) has been proposed in [START_REF] Ginolhac | SAR processors based on signal or interference subspace detector[END_REF]. An interference is assumed to be a dielectric cylinder whose scattering also belongs to a low-rank subspace, called interference subspace. The SISDSAR consists of the difference between the orthogonal projection of the SAR received signal into the target subspace and the orthogonal projection of the SAR received signal into the interference subspace. We observed that the SISDSAR performs generally poorly as the interference responses are widely decreased but the response of the target too.

In this paper, we decide to use a new interference model and to develop a new subspace SAR algorithm which rejects the interferences without reducing the response of the target. This new subspace SAR algorithm is based on oblique projection [START_REF] Behrens | Signal processing applications of oblique projection operators[END_REF] and consists in projecting the SAR received signal along the interference subspace onto the target one. Compared to the SSDSAR, the oblique projection reduces random noise and the interferences at the same time; and unlike the SISDSAR, the MMT response is not reduced. In this paper, we also show the importance of the polarimetric information to obtain good performances. Actually, some SAR algorithms use the polarimetric information in the pre-processing stage but they only focus on polarimetric properties of the random noise [START_REF] Novak | Studies of target detection algorithms that use polarimetric radar data[END_REF], [START_REF] Pastina | Adaptive polarimetric target detection with coherent radar part I[END_REF]. We propose to use the polarimetric information of the MMT and the deterministic interferences. Compared to classical SAR processors, we already showed for the SSDSAR a 8dB detection gain when considering two polarization channels (HH and VV in [START_REF] Brigui | New polarimetric signal subspace detectors for SAR processors[END_REF]). Therefore in this paper, the OBSAR algorithm (for Oblique SAR) is developed by incorporating the polarimetric information in the target and the interference subspaces.

For the polarimetric target subspace, the method presented in [START_REF] Brigui | New polarimetric signal subspace detectors for SAR processors[END_REF] is used. The polarimetric interference subspace generation is given in this paper. We compute the statistical performances of the OBSAR and March 25, 2013 DRAFT compare them to those of the SSDSAR. These results show the importance of the polarimetry. We also study its robustness to interference modeling errors. Results on real FoPen data show the interest and the robustness of this new SAR algorithm compared to classical SAR and SSDSAR [START_REF] Durand | New SAR processor based on matched subspace detector[END_REF].

The paper is organized as follows: section II first describes the SAR configuration, the MMT and interference modeling as well as the construction of the target and interference subspaces. Section III presents our new SAR algorithm, OBSAR, its statistical performances and a study on its robustness to interference modeling errors. Section IV gives the different simulation results, in particular different

Receiving Operating Characteristics (ROC) and finally section V shows the results on real FoPen data.

The following convention is adopted: italic indicates a scalar quantity, lower case boldface indicates a vector quantity and upper case boldface a matrix. T denotes the transpose operator and † the transpose conjugate. CN (a, b) is the complex normal distribution of mean a and variance b and Id d is the identity matrix with dimension d × d.

II. PROBLEM STATEMENT AND MODELING

A. SAR configuration and notation

In SAR (Synthetic Aperture Radar) systems, an antenna on an airborne platform moves along a linear or non-linear trajectory (e.g. circular SAR). In this paper we assume a linear trajectory for clarity of presentation but SAR algorithms presented in this paper can be applied directly for non-linear trajectories as well. As shown in Fig. 1, at each position u i , i ∈ J1, N K, the antenna transmits a signal e and receives the response z i from the scene under observation (we make the stop and go assumption which means that the antenna is not moving when transmitting and receiving). The distance between adjacent positions is constant and equal to δu. For more details on SAR configuration, see [START_REF] Soumekh | Synthetic Aperture Radar Signal Processing[END_REF].

The transmitted signal is a chirp in polarization H and V, with bandwidth B, center frequency f 0 and pulsewidth T e . After demodulation and sampling, we denote the transmitted signal e as follows:

e = e 0 . . . e FsTe-1 0 FsTr-FsTe T

where F s = 1 Ts is the sampling frequency and T r is the integration time. From an electromagnetic point of view, the response of a scatterer in the scene under observation depends on several physical parameters. Firstly, the signal e contains several frequencies; secondly, at each position u i of the antenna, e is transmitted at different incident angles (θ i , ϕ i ) as illustrated in Fig. 2. Finally, e is transmitted for different polarimetric channels; we consider only the HH and VV channels. p = V ), the received signal at the position u i is denoted by z p i ∈ C K , with K = F s T r the number of time samples. Finally, the total received signal z p for one polarization channel is the concatenation of the N vectors z p i :

z p ∈ C N K , z p = z p 1 T z p 2 T . . . z p N T T . (2)
The total polarimetric received signal z is then the concatenation of z H and z V :

z ∈ C 2N K , z = z H T z V T T . (3)
The principal parameters used in the SAR algorithms are listed in Table I.

B. Classical modeling

1) Received signal: In classical SAR algorithms, no prior knowledge is considered on the scatterers.

Thus all the deterministic scatterers are modeled using the white and isotropic point [START_REF] Soumekh | Synthetic Aperture Radar Signal Processing[END_REF]. Without considering random noise, the SAR received signals z p ∈ C N K×1 for single polarization p and z ∈ C 2N K×1 for dual polarization, of any deterministic scatterer at the position (x, y) can then be written as: 2) Polarimetric modeling: As proposed in [START_REF] Brigui | New polarimetric signal subspace detectors for SAR processors[END_REF], polarimetric SAR models r xy ∈ C 2N K×1 including contributions of H and V channels can be generated. To include polarimetric properties to r xy , we write it as follows:

z p = d p
r xy = P   r p xy r p xy   (5)
where P ∈ C 2N K×2N K is a matrix which characterizes the polarimetric properties of the model. Several matrices P can be defined depending of the application. We choose to characterize the polarimetric properties by using Pauli decompositions [START_REF] Cloude | An entropy based classification scheme for land applications of polarimetric SAR[END_REF] which consider canonical scatterers to describe the polarimetric scattering mechanisms. Considering H and V channels, we propose to model two polarimetric mechanisms: the single bounce (or odd bounce) scattering associated to a trihedral corner and the double bounce (even bounce) scattering associated to a dihedral. Thus we can generate two models r + xy ∈ C 2N K×1 and r - xy ∈ C 2N K×1 with different polarimetric properties by using the white and isotropic point:

r + xy = P +   r p xy r p xy   =   1 [N K×N K] 0 [N K×N K] 0 [N K×N K] 1 [N K×N K]    
  =   1 [N K×N K] 0 [N K×N K] 0 [N K×N K] -1 [N K×N K]     r p xy r p xy   =   r p xy -r p xy   (7)
The choice of one of these polarimetric SAR models allows us to select the target to detect according to its polarimetric properties. Using these models to generate SAR images, we can easily shown that these images are equivalent to the polarimetric SAR images classically generated using Pauli decompositions with H and V channels.

C. Man-made target and interference subspace modeling

In FoPen detection, a man-made target is located in a forest environment. The forest represents a strong disturbance with random and deterministic scatterers. We assume that at the frequencies used (P-band)

March 25, 2013 DRAFT only the tree trunk scattering has a deterministic behavior. The other elements of the forest such as the branches or the foliage are assumed to be random scatterers. Thus, we consider two kinds of scatterers of interest with a deterministic scattering behaviour:

• the target to detect which is assumed to be man-made (man-made target or MMT),

• the interferences which are the trunks of the trees.

Rather than the white and isotropic point model used in classical SAR processors [START_REF] Soumekh | Synthetic Aperture Radar Signal Processing[END_REF], we propose to include prior knowledge of the scatterers to improve the SAR algorithm performance. The scattering properties of the target and the interferences are taken into account by considering them as a set of canonical elements which are more complex than the isotropic point. For the sake of clarity, the models are presented for dual polarizations although they are derived in the same way for single polarization.

1) MMT modeling:

To include prior knowledge on the MMT, we assume that the target can be considered as a set of perfectly conducting (PC) plates with several orientations in free space. As the plate possesses an axis of symmetry, only two angles denoted (α, β) are needed to describe its orientation.

These angles are shown in Fig 3.

As the MMT orientation is unknown, one single PC plate with a fixed orientation is not enough to describe the scattering properties of the target. We assume that any response of a PC plate y xy (α, β)

with an orientation (α, β) located at the position (x, y) belongs to a subspace which is assumed to be of low rank:

∀(α, β) ∈ [α min , α max ] × [β min , β max ] y xy (α, β) ∈ H xy (8)
March 25, 2013 DRAFT where H xy is the target subspace of rank D H .

As a consequence, the SAR response of a MMT located at the position (x, y) can be written as:

z = H xy λ xy (9)
where H xy ∈ C 2N K×DH is an orthonormal basis of the target subspace H xy and λ xy ∈ C DH ×1 is an unknown complex coordinate vector which represents the target complex amplitude decomposition in

H xy .
2) Interference modeling: Similarly to the MMT modeling, we include prior knowledge on the interferences by assuming that they can be seen as a set of canonical elements. In the FoPen detection framework, the interferences are the trunks of the trees [START_REF] Thirion | On the use of a coherent scattering model to determine the origin of artificial signatures of a target hidden in a forest[END_REF]. We assume that the interferences can be modeled by dielectric cylinders standing on a PC ground. As the cylinder has an axis of symmetry, only two angles (γ, δ) as shown in Fig. 4 are needed to describe all the possible orientations of the cylinder.

As the orientations (γ, δ) of the trunks are unknown, we assume that the response of a dielectric cylinder x and ẑ are rotated in the plane (x Oẑ ) leading to a new local coordinate system (O, x , ŷ , ẑ ). After the two rotations, the main axis of the cylinder is now along ẑ .

on ground i xy (γ, δ) at (x, y) with an orientation (γ, δ) belongs to a low rank subspace:

∀(γ, δ) ∈ [γ min , γ max ] × [δ min , δ max ] i xy (γ, δ) ∈ J xy
where J xy is the interference subspace of rank D J .

Therefore, the SAR response of an interference located at the position (x, y) can be expressed as:

z = J xy µ xy (10)
March 25, 2013 DRAFT where J xy ∈ C 2N K×DJ is an orthonormal basis of the interference subspace J xy and µ xy ∈ C DJ ×1 is an unknown complex coordinate vector which represents the trunk complex amplitude decomposition in J xy .

D. Construction of the basis H xy and J xy

We explain briefly in this section how to generate the target and interference subspaces and to compute the corresponding basis. For more details on the target subspace, see [START_REF] Durand | New SAR processor based on matched subspace detector[END_REF] for single polarization and [START_REF] Brigui | New polarimetric signal subspace detectors for SAR processors[END_REF] for dual polarization. The generation of the interference subspace in single polarization is presented in [START_REF] Ginolhac | SAR processors based on signal or interference subspace detector[END_REF] and the corresponding generation in dual polarization is given in this paper.

1) Subspace generation:

We first present the generation of the target and interference subspaces for a single polarization p. There are two steps:

• we compute the SAR response y p xy (α, β) using Physical Optics (PO) [START_REF] Ulaby | Radar Polarimetry for Geoscience Application[END_REF] (or i p xy (γ, δ) using the approximation of the truncated infinite cylinder [START_REF] Ulaby | Radar Polarimetry for Geoscience Application[END_REF]) of the canonical element located at the position (x, y), for all the orientations (α, β)

1 covering J0 • , 180 • K × J0 • , 180 • K (or J0 • , 10 • K × J0 • , 360 • K for (γ, δ)) 2 . We assume that (Q, R) (or (Q , R
)) samples are available to cover all the orientations (α, β) (or (γ, δ)).

• using the responses of the PC plate, we generate the target matrix Y p xy ∈ C N K×QR :

Y p xy = y p xy (α 1 , β 1 ) . . . y p xy (α Q , β R ) (11) =      y p 1xy (α 1 , β 1 ) . . . y p 1xy (α Q , β R ) . . . . . . . . . y p N xy (α 1 , β 1 ) . . . y p N xy (α Q , β R )      . ( 12
)
As the responses of a PC plate are the same in HH and VV using PO [START_REF] Ulaby | Radar Polarimetry for Geoscience Application[END_REF], we have

Y H xy = Y V xy ≡ Y xy .
Two interference matrices I H xy and I V xy with the same structure are also generated from the responses of the dielectric cylinder standing over the ground. As its responses are not equal in HH and VV [START_REF] Ulaby | Radar Polarimetry for Geoscience Application[END_REF], we have

I H xy = I V xy .
Subspaces generated from shape-based scattering such as target and interference subspaces have also been studied as dictionaries in [START_REF] Mcclure | Matching pursuits with a wave-based dictionary[END_REF] and for SAR algorithms based on sparse representation in [START_REF] Varshney | Sparse representation in structured dictionaries with application to synthetic aperture radar[END_REF], [START_REF] Samadi | Sparse representation-based synthetic aperture radar imaging[END_REF].

We now present the generation of the subspaces for dual polarization.

• Polarimetric target subspaces: the generation of the polarimetric subspace has previously been studied in [START_REF] Brigui | New polarimetric signal subspace detectors for SAR processors[END_REF]. As a plate in free space describes only one polarimetric mechanism (trihedral type), we can apply the same method used for the white and isotropic point model, to characterize the polarimetric properties of the target. Thus, two polarimetric target matrices have to be generated to cover the two main polarimetric mechanisms associated to the MMT (trihedral-and dihedral-type mechanisms) [START_REF] Cameron | Feature motivated polarization scattering matrix decomposition[END_REF]. The first one Y + xy ∈ C 2N K×DH is obtained using the matrix P + , that corresponds to concatenating the target matrices in the HH and VV polarizations:

Y + xy = P +   Y xy Y xy   =   Y xy Y xy   . ( 13
)
We have shown in [START_REF] Brigui | New polarimetric signal subspace detectors for SAR processors[END_REF] that the polarimetric subspace H + xy generated from the signal matrix Y + xy describes a target with trihedral polarimetric mechanism [START_REF] Lee | Polarimetric radar imaging from basics to applications[END_REF].

A second polarimetric signal matrix Y - xy ∈ C 2N K×DH can be obtained using the matrix P -, that corresponds to concatenating the target matrices in HH and VV polarizations in opposite phase:

Y - xy = P -   Y xy Y xy   =   Y xy -Y xy   . (14)
We have shown in [START_REF] Brigui | New polarimetric signal subspace detectors for SAR processors[END_REF] that the polarimetric subspace H - xy generated from the signal matrix Y - xy describes a target with dihedral polarimetric mechanism [START_REF] Lee | Polarimetric radar imaging from basics to applications[END_REF].

• Polarimetric interference subspace: the cylinder over the ground already describes all the scattering mechanisms: trihedral and dihedral types. Consequently, only one polarimetric interference matrix

I xy ∈ C 2N K×Q R
is needed to fully describe the interferences in dual polarization:

I xy =   I H xy I V xy   . ( 15
)
2) Computation of the basis: The computation of the basis is the same for target and interference subspaces and for single and dual polarizations. For this reason, we only present the computation of the basis H p xy of the target subspace H p xy . The computation of H p xy is done as follows:

• the required orthonormal basis H p xy must minimize the criterion [START_REF] Forster | Generalized rectification of cross spectral matrices for array of arbitrary geometry[END_REF]:

C(H p xy ) = i,j y p xy (α i , β j ) -H p xy † y p xy (α i , β j ) 2 . ( 16
)
March 25, 2013 DRAFT

• the criterion of Eq. ( 16) is satisfied by using the Singular Value Decomposition (SVD) of the target matrix Y p xy :

Y p xy = U p Y xy Σ p Y xy V p † Y xy (17)
where U p Y xy ∈ C N K×min (N K,QR) and V p Y xy ∈ C min (N K,QR)×QR are the left and right singular matrices, and Σ p Y xy ∈ R min (N K,QR),min (N K,QR) is a diagonal matrix containing the singular values. The high computational cost of the basis generation makes it almost impossible to be used on real data.

Nevertheless, techniques to reduce the computation time of the basis generation have been proposed and discussed in [START_REF] Durand | New SAR processor based on matched subspace detector[END_REF]. We sum up the important points of the implementation of the subspace basis:

• choose large step size when sampling orientation angles. A criterion as a function of the rankreduction error is used for this choice.

• reduction of the size of the target and interference matrices by suppressing null elements.

• computation of only one single basis for a position reference (x 0 , y 0 ).

More details on processing complexity and computational times of subspace SAR algorithms can be found in [START_REF] Ginolhac | SAR processors based on signal or interference subspace detector[END_REF], [START_REF] Durand | Back projection version of subspace detector SAR processors[END_REF]. Moreover, accuracy of the subspace models with respect to angular sampling and rank reduction is studied in [START_REF] Durand | New SAR processor based on matched subspace detector[END_REF].

E. Discussion on the subspace models

We propose to discuss here on the subspace models. It is evident that the canonical elements have to present similarities with the scatterer to detect. Nevertheless, generating exact models is neither possible, as the target and the interferences are not a priori exactly known, nor desired. We look for:

• Simplicity of the derivation and fast computation. Indeed, the computation of the subspace basis requires SVD as discussed previously. Moreover, our objective is to develop simple SAR imaging processing close to classical SAR algorithms to be able to apply it on real data.

• Robustness. We aim at detecting a large range of MMT types and assume limited prior-knowledge on the target. In the same way, we assume limited information on the forest.

• Subspace model. Different parameters like size and orientation influence the scattering properties of the target and the interference. As these parameters are unknown, a single canonical element is not able to properly capture scattering features of the target and the interference. In order to take For all these reasons, we decided to favour approximated models, both for the target and the tree scattering features.

1) Target subspace: The MMT scattering depends of the several parameters. We discuss below the approximations of our model.

• Shape. A MMT can easily be assumed as a faceted object standing over the ground. Three scattering mechanisms are involved in the MMT response: trihedral type from the faces, dihedral type from the interaction of the target and the ground and scattering from the edges. We assume that the target subspace generated from plates can describe the two first mechanisms. Indeed, we have shown good detection performances in [START_REF] Durand | New SAR processor based on matched subspace detector[END_REF] using this subspace model and the edge scattering is assumed to be negligible. The choice of the plate model is motivated for its simple derivation, its fast computation and its robustness to MMT scattering. Of course, others shapes can be used as dihedral corner reflector as proposed in [START_REF] Sharma | Matched subspace detectors for discrimination of targets from trees in SAR imagery[END_REF] and the model can be refined by taking into account the edge effect; however one can generate more accurate target subspace but there is a loss in robustness regarding the type of MMT.

• Permittivity. A MMT is generally metallic and then it is reasonable to use the PC approximation.

For known target permittivity, it is possible to include it in the target model.

• Size. It is difficult to precisely know the size of any MMT. Moreover, the size of the canonical element influences the resolution of the SAR image. We propose to use one fixed size which corresponds to the size of a pixel for the resolution obtained with classical SAR processing. We could consider different sizes (lower than the classical SAR image resolution) of the plate to generate the target subspace but it would increase the subspace rank and then increase the noise intensity.

• Orientation. The orientation of the target is taken into account in our model.

The way we include polarimetric information is based on existing methods. Several polarimetric decompositions have been proposed such the H/α decomposition [START_REF] Cloude | An entropy based classification scheme for land applications of polarimetric SAR[END_REF] which can be interesting for FoPen detection. However the H/α decomposition and many of others are difficult to include in our subspace models and SAR image processing scheme so we only consider simple polarimetric decompositions. We include polarimetric information by following the Pauli decomposition approach [START_REF] Lee | Polarimetric radar imaging from basics to applications[END_REF]. We have previously shown for classical SAR imaging that the Pauli decomposition can be seen as using polarimetric models with the matrices P ± defined in Equations ( 6) and [START_REF] Durand | New SAR processor based on matched subspace detector[END_REF]. We use the same polarimetric matrices and apply March 25, 2013 DRAFT them to the target subspace models in HH and VV. In that way, we do not need to change our canonical element model for single and dual polarizations. Moreover, the use of the matrices P ± allows us to select polarimetric mechanism we want to detect. We have previously shown in [START_REF] Brigui | New polarimetric signal subspace detectors for SAR processors[END_REF] the effectiveness of this polarimetric model for MMT detection in white Gaussian noise.

2) Interference subspace: For FoPen application, the interference is the forest. In our study, we only consider the trees of the forest as they are the main cause of false alarms. We only take into account the direct scattering from the trunk and the scattering from the interaction of the trunk with the ground. The interference model is an approximation of the true trunk scattering. Considering true tree scattering, we have to consider several parameters. We then discuss on the approximations of our interference model.

• Size and permittivity. The size of the trees of a forest can be different from tree to tree. The size mainly changes the magnitude and the polarimetric properties of the tree scattering. We propose to use a fixed size and permittivity of the cylinder according to the type of forest. Average height of the trees can be estimated using interferometry to have a more accurate model. Permittivity of the ground influences the scattering interaction of the trunk with the ground and can be highly variable with the meteorological conditions (for example scattering interaction with dry ground is negligible).

As we want to include this scattering interaction feature in our model, permittivity of the ground is chosen to be PC.

• Branches and foliage. We do not consider in our interference model the branches and foliage which is reasonable assumption at low frequencies. Their scattering is assumed to be a zero mean complex Gaussian noise. Nevertheless, for highly dense forest this approximation may not be valid any more and induces some loss of interference rejection. It would be then interested to assume that the branches and foliage scattering is a heterogeneous noise with a given distribution and to adaptively estimate it. This method has been proposed in [START_REF] Bandiera | Adaptive radar detection without secondary data: Exploiting a diversity idea[END_REF] for sea clutter using Doppler diversity and we can apply in a future work it to SSDSAR algorithm using directivity and polarimetry diversity.

• Orientation. Our interference model takes into account the orientation of the trees.

• Multiple reflections between trunks. The interactions between close objects have been shown to be negligible when the distance between these objects is superior to 2λ 0 (wavelength of the emitted signal) [START_REF] De Badereau Nd | Radar remote sensing of forests at low frequencies: a 3D electromagnetic scattering model[END_REF]. Then our model does not take into account multiple reflections between trunks.

The interference subspaces in HH and VV already include all the polarimetric mechanisms of a trunk over the ground. The polarimetric subspace model for the interference is then just the concatenation of the susbspace models in HH and VV. xy and J V xy for single polarizations; for dual polarization, only one interference subspace is generated and its orthonormal basis is denoted by J xy . To evaluate the accuracy of the subspace model, we propose in Section IV-C to compute the ratio between the target (or the interference) signal energy projected in the target (or interference) subspace model and the total target (or interference) signal energy.

III. SAR ALGORITHMS

This section is dedicated to the development of our new SAR algorithms by using both the target and the interference subspaces. As the derivations of these algorithms are the same for single or dual polarization, we consider only the dual polarization case.

To form SAR images, conventional techniques are based on the Fourier transform; nevertheless these techniques do not allow to easily integrate SAR subspace models in the image processing. To overcome this issue, SAR images in this paper are generated using estimation techniques. More details on the use of estimation methods for SAR image formation can be found in [START_REF] Degraaf | SAR imaging via modern 2-D spectral estimation methods[END_REF], [START_REF] Wu | SAR image formation via semiparametric spectral estimation[END_REF].

A. Classical SAR (CSAR)

In classical SAR algorithms such as the TDCA (Time Domain Correlation Algorithm) or the RMA (Range Migration Algorithm) [START_REF] Soumekh | Synthetic Aperture Radar Signal Processing[END_REF], any scatterer is modeled by the isotropic point presented in Eq. ( 4).

We consider a target located at the position (x, y) whose scattering is corrupted by n = CN (0, σ 2 Id 2N K ) which is a zero mean complex Gaussian noise with known variance σ 2 (the derivation of subspace SAR algorithms in the case of unknown variance is shown in [START_REF] Durand | New SAR processor based on matched subspace detector[END_REF]). The SAR received signal z is written as:

z = d xy r xy + n ( 18
)
where d xy is an unknown complex amplitude and r xy is the dual polarization SAR response of an isotropic point located at (x, y) defined in Eq (5).

The unknown complex amplitude d xy can be estimated using the least square method [START_REF] Trees | Estimation and Modulation Theory[END_REF]. This estimation is done for each position (x, y). Denoting the estimate of the unknown complex amplitude for the position (x, y) as dxy , the intensity of the CSAR image is then defined as follows:

I CS (x, y) = dxy 2 σ 2 = z † P rxy z σ 2 = r † xy z 2 σ 2 (19)
March 25, 2013 DRAFT

where P rxy = r xy r † xy is the orthogonal projector into the subspace of rank 1 spanned by r xy . From

Eq. ( 19), we see that the classical approach is the matched filter of the received signal z. It has been

shown in [START_REF] Durand | New SAR processor based on matched subspace detector[END_REF] that the CSAR image is totally equivalent to those obtained by conventional SAR algorithms (TDCA, BackProjection, RMA). Besides, we can easily show that the CSAR images for dual polarizations are equivalent to the polarimetric images classically generated using the Pauli decompositions with H and V channels.

B. Signal Subspace Detector SAR (SSDSAR)

The Signal Subspace Detector SAR (SSDSAR) has been developed previously in [START_REF] Durand | New SAR processor based on matched subspace detector[END_REF] for single polarization and in [START_REF] Brigui | New polarimetric signal subspace detectors for SAR processors[END_REF] for dual polarization. The SSDSAR image is generated by including prior knowledge on the scattering of the MMT. We consider that a target whose scattering is modeled as in Eq. ( 9) is located at the position (x, y) and its signal is corrupted by n. The SAR received signal z is written as:

z = H xy λ SSD xy + n (20)
where λ SSD xy is an unknown coordinate vector and H xy is equal to H + xy or H - xy . As in the case of CSAR, the unknown coordinate vector λ SSD xy is estimated using least square method [START_REF] Trees | Estimation and Modulation Theory[END_REF], [START_REF] Scharf | Matched subspace detectors[END_REF]. For each position (x, y), we compute the estimate λSSD xy of the unknown coordinate vector.

The intensity of the SSDSAR image is then defined as follows:

I SSD (x, y) = λSSD xy 2 σ 2 = z † P Hxy z σ 2 = H † xy z 2 σ 2 (21)
where P Hxy = H xy H † xy is the orthogonal projector into the subspace H xy . We clearly see that the intensity is the square norm of the projection of the received signal z along the direction orthogonal to

H xy onto H xy .
We study the property of the SSDSAR according to a scatterer located in the pixel (x, y):

• target: as P Hxy H xy λ xy = H xy λ xy , the signal from a MMT whose scattering belongs to H xy is then unchanged by the projection P Hxy .

• interference: as P Hxy J xy µ xy = H xy H † xy J xy µ xy , the orthogonal projection of an interference whose scattering belongs to J xy actually depends on the product H † xy J xy . As it is almost certain that the interference and the target subspaces are not orthogonal, the response of an interference by the SSDSAR will not be null and will be the principal cause of false alarms.

March 25, 2013 DRAFT

• random noise: as H xy is a low rank subspace, the projection of the random noise is low. Thereby the SSDSAR increases significantly the detection of MMT in the presence of white Gaussian noise as shown in [START_REF] Durand | New SAR processor based on matched subspace detector[END_REF], [START_REF] Brigui | New polarimetric signal subspace detectors for SAR processors[END_REF].

C. Oblique SAR (OBSAR)

The Oblique SAR (OBSAR) image is generated by including prior knowledge on the scattering of the MMT and the interferences. The scattering properties of the target are taken into account by using the target subspace and those of the interferences by using the interference subspace. We consider that a target whose scattering is modeled as in Eq. ( 9), is located at the position (x, y) and its signal is corrupted by an interference modeled as in Eq. ( 10) and random noise. The SAR received signal z is written as:

z = H xy λ OB xy + J xy µ xy + n (22)
where λ OB xy and µ xy are unknown coordinate vectors. Once again, we estimate λ OB xy using the least squares [START_REF] Behrens | Signal processing applications of oblique projection operators[END_REF]. The estimate λOB xy is derived by resolving the following equation:

Finally, the estimate λOB xy is written as:

λOB xy = H † xy (Id 2N K -J xy (J † xy P ⊥ Hxy J xy ) -1 J † xy P ⊥ Hxy )z = H † xy E HxyJxy z (25)
where P ⊥ Hxy = Id 2N K -P Hxy and E HxyJxy is the oblique projector along J xy onto H xy [START_REF] Behrens | Signal processing applications of oblique projection operators[END_REF].

The estimation of λOB xy is done for all positions (x, y). The intensity of the OBSAR image for the position (x, y) is then defined as:

I OB (x, y) = λOB xy 2 σ 2 = H † xy E HxyJxy z 2 σ 2 . ( 26
)
It is clear from Eq. ( 26) that the OBSAR image is the square norm of the oblique projection of z along the interference subspace J xy onto the target subspace H xy .

We study the property of the OBSAR according to a scatterer located in the pixel (x, y):

• target: as E HxyJxy H xy λ xy = H xy λ xy , the signal from a MMT whose scattering belongs to H xy is then unchanged by the oblique projection E HxyJxy . In terms of the response of the MMT, we obtain the same result as for the SSDSAR.

• interference: since E HxyJxy J xy µ xy = 0 2N K×1 , the signal from an interference whose scattering belongs to J xy is totally suppressed.

• random noise: as H xy and J xy are low-rank subspaces, the oblique projection of the random noise is low.

The main advantage of OBSAR is its ability to reduce false alarms causing by both random noise and deterministic interference without affecting the detection of MMT.

D. Discussion on robustness of the OBSAR

The target detection performance and false alarm reduction of the OBSAR depend greatly on the modeling accuracy of the MMT and of the interferences.

In [START_REF] Durand | New SAR processor based on matched subspace detector[END_REF], the SSDSAR using the target subspace shows good performance and robustness with respect to the size of the target model both for simulated data and real data, where faceted targets were considered.

Therefore, the MMT modeling errors are assumed to be negligible in this paper.

However, this assumption does not hold anymore for the interference modeling. In reality, the interference environment is not exactly known (height and radius of trunks, permittivity, ground not PC, . . .).

To study the robustness of the OBSAR to interference modeling errors, we consider two cases: the ideal case in which the interference scattering belongs totally to the interference subspace and the realistic case. For the latter, we suppose that a part δ J xy of the interference scattering does not belong to J xy . As shown in Fig. 5, the OBSAR intensity of the pixel (x, y) is no longer null compared to the ideal case:

I OB (x, y) = H † xy E HxyJxy δ J xy 2 σ 2 . ( 27
)
However, we notice that the intensity for the OBSAR is still lower than that of the SSDSAR. In subsections IV-B and IV-D, we investigate the impact of the interference modeling errors on the OBSAR performances.

March 25, 2013 DRAFT Another possible error in OBSAR results could be caused by some elements of the forest which are not taken into account in the interference subspace or the random noise (the branch and leave scattering, foliage attenuation, . . .). OBSAR is applied to real data in Section V to study its robustness to this kind of errors.

E. Statistical performances

The probabilities of detection and false alarms associated to the SSDSAR and the OBSAR are derived in this section. To evaluate the robustness of the algorithms to the interference modeling, we consider the ideal and realistic cases.

1) Probability of detection P d : We suppose that a MMT is located at the position (x 1 , y 1 ). We also assume that the MMT scattering lies totally in the target subspace H xy . The SAR received signal and its distribution are written as:

z 1 = H x1y1 λ x1y1 + n, z 1 ∼ CN (H x1y1 λ x1y1 , σ 2 Id 2N K ). (28)
The probability of detection is the probability that the intensity of the pixel (x 1 , y 1 ) is higher than a certain threshold η. From the definition of the intensity for SSDSAR and OBSAR, we can derive the distribution as follows:

• P SSD d = P (I SSD (x 1 , y 1 ) > η)
From Eq. ( 21), the distribution of I SSD (x 1 , y 1 ) is non-central chi-square χ 2 [START_REF] Scharf | Statistical Signal Processing: Detection, Estimation and Time Series Analysis[END_REF]:

I SSD (x 1 , y 1 ) ∼ 1 2 χ 2 2D H , 2 λSSD x1y1 2 σ 2 (29)
where χ 2 (a, b) denotes the chi-square distribution with degrees of freedom a and non-centrality parameter b.

• P OB d = P (I OB (x 1 , y 1 ) > η)
By using the Eq. ( 26) of the OBSAR intensity, the distribution of I OB (x 1 , y 1 ) is also non-central χ 2 [START_REF] Behrens | Signal processing applications of oblique projection operators[END_REF]:

I OB (x 1 , y 1 ) ∼ 1 2 χ 2 2D H , 2 λOB x1y1 2 σ 2 . ( 30
)
For the signal z 1 , the intensities of the SSDSAR and OBSAR are equal. Hence the detection performances for a target whose scattering belongs totally to the target subspace, are the same for the SSDSAR and the OBSAR.

March 25, 2013 DRAFT

2) Probability of false alarm P f a :

• Ideal Case: we suppose that an interference whose scattering belongs totally to the interference subspace J xy is located at the position (x 0 , y 0 ). The SAR received signal and its distribution are written as:

z 0 = J x0y0 µ x0y0 + n, z 0 ∼ CN (J x0y0 µ x0y0 , σ 2 Id 2N K ). ( 31
)
The probability of false alarm is the probability that the intensity of the pixel (x 0 , y 0 ) is higher than the threshold η.

-

P SSD f a = P (I SSD (x 0 , y 0 ) > η)
The distribution of I SSD (x 0 , y 0 ) is non-central χ 2 [START_REF] Scharf | Statistical Signal Processing: Detection, Estimation and Time Series Analysis[END_REF]:

I SSD (x 0 , y 0 ) ∼ 1 2 χ 2 2D J , 2 H † x0y0 J x0y0 µ x0y0 2 σ 2 . ( 32
)
-P OB f a = P (I OB (x 0 , y 0 ) > η) For an interference whose scattering lies totally in J x0y0 , the OBSAR gives a null intensity.

Hence, the distribution of I OB (x 0 , y 0 ) is central χ 2 :

I OB (x 0 , y 0 ) ∼ 1 2 χ 2 (2D J , 0). (33)
• Realistic Case: we consider an interference such that part of its scattering δ x0y0 does not belong to

J xy : z = J x0y0 µ x0y0 + δ J x0y0 + n, z ∼ CN (J x0y0 µ x0y0 + δ J x0y0 , σ 2 Id 2N K ) . ( 34
)
We have then the following probabilities of false alarm:

-P SSD f a = P (I SSD (x 0 , y 0 ) > η)
The distribution of I SSD (x 0 , y 0 ) is non-central χ 2 [START_REF] Scharf | Statistical Signal Processing: Detection, Estimation and Time Series Analysis[END_REF]:

I SSD (x 0 , y 0 ) ∼ 1 2 χ 2 2D J , 2 H † x0y0 (J x0y0 µ x0y0 + δ J x0y0 ) 2 σ 2 . ( 35
)
-P OB f a = P (I OB (x 0 , y 0 ) > η) The distribution of I OB (x 0 , y 0 ) is also non-central χ 2 :

I OB (x 0 , y 0 ) ∼ 1 2 χ 2 2D J , 2
H † x0y0 E Hx 0 y 0 Jx 0 y 0 δ J x0y0 2 σ 2 . ( 36
)
IV. SIMULATED DATA This section is dedicated to the application of the SAR algorithms to simulated data. We first present the statistical performances with ROC (Receiving Operating Characteristic) curves. Then, we show and comment on the images given by the SAR algorithms. 2) Subspaces: The target subspaces for single or dual polarization are generated using 2m × 1m PC plates; their scattering are computed with PO [START_REF] Ulaby | Radar Polarimetry for Geoscience Application[END_REF]. The sampling step for the orientation angles is 9 • and

(α, β) ∈ [0, 180 • ] × [0, 180 • ] [7]
. Moreover, the subspaces spanned by the PC plates need to be low-rank. Fig. 7 shows the curves of the singular values of the signal matrix Y xy defined in Eq. ( 12) for a single polarization and Y ± xy defined in Eq. ( 13) and in Eq. ( 14) for dual polarization 3 . As the singular values do not exhibit clear cut-off, the choice of the subspace ranks is a trade-off between the description of MMT scattering with unknown orientation and the rejection of random noise. In addition, it has been shown in [START_REF] Durand | New SAR processor based on matched subspace detector[END_REF] that choosing a too low rank yields robustness loss. On the contrary, a too high rank degrades detection performances. From these curves, the ranks of the target subspaces are chosen to be equal to D H H = D V H = 10 for single polarization [START_REF] Durand | New SAR processor based on matched subspace detector[END_REF] and D + H = D - H = 10 for dual polarization [START_REF] Brigui | New polarimetric signal subspace detectors for SAR processors[END_REF]. The interference subspaces for single and dual polarization are generated using a dielectric cylinder with a permittivity of (22.96,-11.7), a height of 11m and a radius of 20cm standing on a PC ground; their scattering is computed using the approximation of the truncated infinite cylinder [START_REF] Ulaby | Radar Polarimetry for Geoscience Application[END_REF]. The sampling step for the orientation angles is 2

• and (γ, δ) ∈ [0, 10 • ] × [0, 360 • ] [8]
. As for the target subspaces, we plot the curves of the singular values of the interference matrices to determine their ranks. Fig. 8 show the singular values of the interference matrices I HH xy and I V V xy for single polarization. The singular values of the interference matrix I xy for dual polarization have also been plotted in Fig. 8. The ranks of the interference subspaces are chosen to be equal to D H J = D V J = 10 for single polarization and D J = 10 for dual polarization.

3) Target and interferences:

The MMT is a PC box with a size of (2m × 1.5m × 1m) over a flat PC ground located approximately at the center of the scene and its scattering is simulated using Feko [START_REF] Feko | User's Manual, EM Software and System[END_REF] which employs the method of moment (MoM) to compute the electromagnetic scattering. The target is oriented such that its larger dimension is parallel to the flight path. The target is placed in a simulated 3 The matrices Y + xy and Y - xy have the same singular values [START_REF] Brigui | New polarimetric signal subspace detectors for SAR processors[END_REF].

March 25, 2013 DRAFT forest (see Fig. 6). Since the principal cause of false alarms at the frequencies used in our simulations (f 0 = 400MHz, B = 100MHz) is the scattering of the trunks [START_REF] Thirion | On the use of a coherent scattering model to determine the origin of artificial signatures of a target hidden in a forest[END_REF], we simulate mainly the scattering of the trunks.

For the ideal case, trunk scattering belongs totally to the interference subspaces: their scattering is computed with the same dielectric cylinder lying over a PC ground used to generate the interference subspaces.

For the realistic case, the computation of the forest scattering is done using COSMO [START_REF] Thirion | Application of a coherent model in simulating the backscattering coefficient of mangrove forest[END_REF], [START_REF] Thirion | Capabilities of a forest coherent scattering model applied to radiometry, interferometry and polarimetry at P-and L-bands[END_REF], a software dedicated to the study of electromagnetic scattering by forests and based on asymptotic methods.

More information on this modeling tool and additional references to other existing tools to simulate forest scattering are given in [START_REF] Thirion | Application of a coherent model in simulating the backscattering coefficient of mangrove forest[END_REF], [START_REF] Thirion | Capabilities of a forest coherent scattering model applied to radiometry, interferometry and polarimetry at P-and L-bands[END_REF]. We choose the trunks to be dielectric cylinders with a height of 11m and a radius of 20cm. Their orientations (γ, δ) are randomly distributed with γ ∈ [0, 10 • ] and δ ∈ [0, 360 • ].

Compared to the ideal case, the ground is not PC but dielectric with a permittivity of (43.55,-0.3) and the trunk scattering is attenuated by the canopy.

B. Receiver Operating Characteristic (ROC)

To evaluate the performances of the SSDSAR and the OBSAR, we plot the ROC, P d against P f a .

We consider the interferences for ideal and realistic cases described in Section IV-A.3. In both cases, a high Signal to Noise Ratio (SNR) of 35dB is chosen in order to better study the contributions of the interferences to the performances and the robustness of both algorithms. In both ideal and realistic cases, Signal to Interference Ratios (SIR) are equal to -8dB in HH polarization, 3dB in VV polarization and -6dB in dual polarization.

1) Single polarization:

The ROC curves are plotted in Fig. 9 for the ideal case and in Fig. 10 for the realistic case. In both cases and for both polarizations (HH or VV), we obtain the same performances for the SSDSAR and the OBSAR. The target and the interference subspaces are too similar to differentiate the target from the interferences. In single polarization, the reduction of false alarm due to interference is therefore not possible. We notice that in both cases the SAR algorithms in VV polarization outperform the SAR algorithms in HH polarization; this result is due to the polarimetric properties of the the tree trunk which are different in HH and VV polarizations.

2) Dual polarization (dihedral type): We use the polarimetric target subspace H - xy as the box over the ground has a dihedral type scattering [START_REF] Brigui | New polarimetric signal subspace detectors for SAR processors[END_REF]. The ROC curves are plotted in Fig. 11 for the ideal case and in Fig. 12 for the realistic case. For the ideal case, we obtain P d higher than 0.9 for P f a higher than 2.10 -4 with the OBSAR and for P f a higher than 0.8 with the SSDSAR. This result shows the greatly improved performances of the OBSAR compared to the SSDSAR. Therefore, the oblique projection of the OBSAR algorithm is needed for a complete removal of the interferences. For the realistic case, we obtain P d higher than 0.9 for P f a higher than 8.10 -2 with the OBSAR and for P f a higher than 0.1 with the SSDSAR. The performances of both algorithms are degraded compared to those obtained in the ideal case. This may be explained by the fact that foliage attenuation and the non PC ground are not taken into account in the interference subspace. However, we notice that the OBSAR still outperforms the SSDSAR.

C. Accuracy of subspace models

We propose to evaluate the accuracy of the subspace models. We compute the ratio between the energy of the projected signal of MMT onto the target subspace and the total energy of the MMT signal energy.

Table II and III show these ratios for single and dual polarizations using white and isotropic point models and target subspaces. We see first that the white and isotropic point model just captures a quarter of the energy of the MMT signal. In the other hand, the target subspaces capture almost the entire energy of the MMT signal for both HH and VV channels. The loss of MMT signal energy can be explained by the truncation of the target subspaces and by the approximation of the target modeling. We draw the same comments on the accuracy of the target subspace for dual polarizations compared to the white and isotropic point model. As the main scattering mechanism for the box over the ground is dihedral type, the MMT scattering is mainly describe by the target subspace H - xy . The same comparison can be done with the trunks and the interference subspaces. Around 92% of trunk signals are described on average [START_REF] Forster | Generalized rectification of cross spectral matrices for array of arbitrary geometry[END_REF] DRAFT by the interference subspaces for single polarizations and around 80% for dual polarization.

D. Images

In this subsection, we only consider the realistic case for the interferences. We propose to compare the SAR images obtained using the CSAR, the SSDSAR and the OBSAR algorithms. In order to quantify and compare the performances of these SAR algorithms, we compute for each image the ratio ρ between the intensity of the target and the maximum intensity of the interferences defined by: ρ = 10 log 10 I(x t , y t )

I(x i , y i ) (37)
where I(x, y) is the intensity of the pixel (x, y) on the image processed by a given SAR processor.

(x t , y t ) is the pixel containing the target and (x i , y i ) is the pixel containing the interference with the maximum intensity.

We specify that all the SAR images are shown in intensity (square modulus) normalized by the maximum value; the scale of the intensity values is linear.

1) Single polarization:

The images in the HH polarization are presented in Fig. 13(a) for the CSAR, in Fig. 13(b) for the SSDSAR and in Fig. 13(c) for the OBSAR. First of all, we clearly see that the CSAR image does not allow any detection of the target because its response is very low and there are a lot of false alarms due to the trunks. The SSDSAR image shows that the target response is increased significantly;

nevertheless, the tree trunks responses remain high and still lead to a lot of false alarms. Unfortunately, this problem is not solved in the OBSAR image: the target and the trunks responses are unchanged from those in the SSDSAR. The ratio ρ CSAR is equal to -3.7dB and we have ρ SSDSAR = ρ OBSAR = -0.5dB.

The images in the VV polarization are presented in Fig. 14(a) for the CSAR, in Fig. 14(b) for the SSDSAR and in Fig. 14(c) for the OBSAR. We can draw the same conclusions as for the HH polarization, except that the forest has a lower response in general while the target response does not vary much: therefore, it is easier to distinguish the target from the interferences in the SSDSAR and the OBSAR images. The ratio ρ CSAR is equal to -1.8dB and we have ρ SSDSAR = ρ OBSAR = 2.8dB.

These SAR images clearly illustrate that the rejection of interferences is not possible using single polarization. The scattering of MMT and interferences cannot be discriminated as the two subspaces Therefore, these results confirm conclusions obtained in Section IV-B: the reduction of the false alarms with the OBSAR algorithm is only possible using dual polarization. In this case, the signal and the interference subspaces are sufficiently far apart.

V. REAL DATA

A. Configuration

The real SAR data presented in this section has been acquired during the PYLA'2004 campaign over the Nezer forest in the Landes (France), using the SAR system RAMSES from ONERA (the French Aerospace Lab) at frequencies between 400MHz and 470MHz and with an incidence angle of 59.8 • . For more details on the Nezer forest, see [START_REF] Thirion | Capabilities of a forest coherent scattering model applied to radiometry, interferometry and polarimetry at P-and L-bands[END_REF]. Two targets, a truck and a trihedral corner reflector, with an orientation parallel to the flight path were placed in the forest of pines as shown in Fig. 16.

We use 4m × 2m PC plates to generate the signal subspaces and PO to compute their scattering.

The sampling step for the orientation angles which cover ([0 and Fig. 19(c), respectively. We clearly distinguish the truck and the trihedral in the SSDSAR image compared to the CSAR. As in the case of simulated data, the OBSAR image does not show significant improvement compared to the SSDSAR image. For the trihedral, the ratios are ρ SSDSAR = ρ OBSAR = 1.5dB. For the truck, the ratios are ρ SSDSAR = ρ OBSAR = 0.8dB.

• 180 • ] × [0 • 180 • ]) is 9 • . The
2) Dual polarization: We present the SAR images for dual polarization. We first consider the polarimetric target subspace H + xy to detect target whose scattering is of trihedral type. The CSAR, SSDSAR and OBSAR images are presented in Fig. 20(a), Fig. 20(b) and Fig. 20(c), respectively. Only the trihedral corner reflector is detected as the truck scattering is of dihedral type. We clearly distinguish the corner reflector in the SSDSAR and the OBSAR images compared to the CSAR one. Moreover, the interferences are slightly reduced. The ratio ρ SSDSAR is equal to 1.5dB while for the OBSAR, ρ OBSAR = 2dB.

We present now the SAR images for dual polarization using the target subspace H - xy for the dihedraltype target. The CSAR, SSDSAR and OBSAR images are shown in Fig. 21(a), Fig. 21(b) and Fig. 21(c), respectively. Contrary to the previous results, only the truck is detected. Once again, we clearly distinguish the truck from the environment in the SSDSAR and OBSAR images compared to that for CSAR. Finally, the interferences are slightly reduced with ρ SSDSAR = 1.7dB for the SSDSAR and ρ OBSAR = 2.3dB for the OBSAR.

As we have limited information about the Nezer forest to build an accurate interference subspace, interference rejection is not optimal. Better performance could be obtained with a better modeling of the interferences. But, this accuracy is possible only if information about the forest is available. Nevertheless in dual polarization, these results show the robustness of the OBSAR algorithm which always outperforms the SSDSAR algorithm in this case.

VI. CONCLUSION

We have proposed a new algorithm OBSAR to form SAR image for FoPen application. The aim of the OBSAR is to improve target detection and interference rejection compared to conventional SAR processors or to the SSDSAR. OBSAR intensity for each pixel of the SAR image to generate is computed by solving an estimation problem with subspace models for the target and the interferences. The estimated intensity of the target in each pixel of the OBSAR image is then the oblique projection along the interference subspace onto the target subspace. We computed the statistical performances of the OBSAR and demonstrated the importance of the polarimetric information for interference rejection. Different results showed that the OBSAR outperforms the SSDSAR. In order to use the OBSAR in realistic cases, March 25, 2013 DRAFT we also studied its robustness to interference modeling errors. Even if the performances are then degraded, false alarms due to the interferences are still reduced. Finally, we validated these results using real data.

For future work, interference modeling needs to be improved to obtain better false alarm reduction. For example, principal branch responses and foliage attenuation effects could be included in our interference model [START_REF] Thirion | Application of a coherent model in simulating the backscattering coefficient of mangrove forest[END_REF], [START_REF] Thirion | Capabilities of a forest coherent scattering model applied to radiometry, interferometry and polarimetry at P-and L-bands[END_REF]. It would also be interesting to use cross-polarized channels to generate target and interference subspaces that are further apart, leading improved target detection and false alarm reduction.

Concerning the case of the MMT, it has to be noted that the choice of the model for the subspace generation may have an impact on the detection performances. Indeed a more accurate representation leads to a more realistic simulation and probably to a better detection when the MMT is perfectly known. However, it decreases the robustness in most of cases, for unknown MMT. Distributed background scattering and correlation between pixels of the forest area [START_REF] Jackson | A model for generating synthetic VHF SAR forest clutter images[END_REF] could be also another solution to improve interference modeling. At last the use of classical polarimetric post-processing tools may be investigated.

Fig. 1 .

 1 Fig. 1. SAR geometric configuration.

xy r p xy , z = d xy r xy ( 4 )Fig. 2 .

 42 Fig. 2. Definition of the incidence angles (θi, ϕi).

Fig. 3 .

 3 Fig. 3. Definition of orientation angles (α, β) for the PC plate. (a) The PC plate lies originally on the plane (xOy). (b) When applying a counter-clockwise rotation α about the axis x, the axes ŷ and ẑ are rotated in the plane (yOz). The PC plate is rotated in the same way. (c) A second rotation is applied to the new axis ŷ . The axes x and ẑ are then rotated by β in the plane (x0y ). We finally obtain a new coordinate system (O, x , ŷ , ẑ ) and there again the PC plate has been rotated likewise.

Fig. 4 .

 4 Fig. 4. Definition of orientation angles (γ, δ) for the dielectric cylinder. (a) The main axis of the cylinder is along ẑ. (b) A first rotation δ about the ẑ-axis leads to a new local coordinate system (O, x , ŷ , ẑ). (c) A second rotation γ is applied on ŷ .

  [START_REF] Forster | Generalized rectification of cross spectral matrices for array of arbitrary geometry[END_REF] DRAFT into account these variations, subspace models are generated from canonical element responses with different parameters.

  xy λ OB xy -J xy µ xy 2 = H † xy z -H † xy J xy μxy (23) where μxy is obtained by solving the following equation: μxy = arg min µxy z -H xy λ OB xy -J xy µ xy 2 = J † xy z -J † xy H xy λOB xy .

  Geometry: We consider a flight between the first position u 1 = -50m and the last position u 200 = 50m with δ u = 0.5m between each position and an altitude of 100m. The polarimetric transmitted signal is a chirp with a bandwidth B = 100MHz and a center frequency f 0 = 400MHz; its pulsewidth is T e = 2.10 -7 s. Finally, the radar scene is a 50m × 45m rectangle as illustrated in Fig 6.The antenna transmits and receives the signal in H and V polarizations. For this study, we consider single HH or VV polarization and dual polarization (HH and VV).

  H p xy and J p xy are too close. Nevertheless, the SSDSAR and the OBSAR algorithms enhance the detection of the target compared to the CSAR algorithm. Finally, the OBSAR does not degrade the response of the MMT even if the target and interference subspaces are close. March 25, 2013 DRAFT 2) Dual polarization (dihedral type): We present now the images of the CSAR in Fig. 15(a), of the SSDSAR in Fig. 15(b) and of the OBSAR in Fig. 15(c) for dual polarization. Similarly to the single polarization case, the target intensity in the CSAR image is low and false alarms due to the trunks are numerous. The ratio ρ CSAR is equal to -3.5dB. The SSDSAR image in dual polarization is similar to that obtained in single polarization: the target appears clearly with a ratio ρ SSDSAR of 1.8dB but the responses of the trunks are still high. In the OBSAR image, the responses of the trunks are reduced compared to the SSDSAR image with a ratio ρ OBSAR of 3.6dB. Compared to the SSDSAR, the interference intensities are greatly reduced while the target intensity is unchanged.

  ranks of the subspaces are determined from the singular values of the signal matrices shown in Fig.17for single and for dual polarizations. The ranks are chosen to be D H H = D V H = 10 and D + H = D - H = 10. The interference subspaces are generated using dielectric cylinders with a height of 11m and a radius of 20cm standing on a PC ground. The sampling step for the orientation angles is 2 • and (γ, δ) ∈ [0, 10 • ] × [0, 360 • ]. The singular values of the interference matrices are shown in Fig. 18 for HH and VV polarizations as well ad for dual polarization. We choose D H J = D V J = 10 and D J = 10. B. Images 1) Single polarization: As the response of the forest is too high in HH polarization making it impossible to detect the target, results of the three SAR algorithms for HH polarization are not presented in this paper. March 25, 2013 DRAFT We present the images of the CSAR, SSDSAR and OBSAR for VV polarization in Fig. 19(a), Fig. 19(b)

( a )

 a First step of oblique projection for ideal and realistic cases. The signal is projected in the subspace orthogonal to the random noise direction. (b) Ideal Case. The orthogonal projection of the interference response on the target subspace is not null while the oblique projection of the interference response on the target subspace is null.(c) Realistic Case. Compared to the ideal case, the oblique projection of the interference response on the target subspace is not null but is still lower than the orthogonal projection of the interference response on the target subspace.

Fig. 5 .Fig. 6 .

 56 Fig. 5. Interference projections for the ideal and the realistic cases. The orthogonal projection P Hxy projects the interference response along the direction orthogonal to Hxy onto the target subspace Hxy . The oblique projection E Hxy Jxy projects the interference response along the direction parallel to Jxy onto Hxy . March 25, 2013 DRAFT

Fig. 7 .

 7 Fig. 7. Singular values of the target matrices for single polarization (HH or VV) and dual polarization for simulated data. As the target matrices are the same for HH polarization and VV polarization, their singular values are the same for single and dual polarization. The ranks of the target subspaces are chosen to be D H H = D V H = 10 for single polarization and D + H = D - H = 10 for dual polarization.

Fig. 8 .

 8 Fig. 8. Singular values of the interference matrices for single polarization (HH and VV) and dual polarization for simulated data. The ranks of the interference subspaces are chosen to be D H J = D V J = 10 for single polarization and DJ = 10 for dual polarization.

Fig. 9 .Fig. 10 .

 910 Fig. 9. ROC of SSDSAR and OBSAR for the ideal case in single polarization (HH and VV) and for simulated data. The SSDSAR and the OBSAR give the same performance for single polarization.

Fig. 11 .Fig. 12 .

 1112 Fig. 11. ROC of SSDSAR and OBSAR for the ideal case in dual polarization and for simulated data. The OBSAR algorithm gives good performances compared to the SSDSAR algorithm.

Fig. 13 .

 13 Fig. 13. SAR images for HH polarization using simulated FoPen data with a MMT located at position (108, -1)m. (a) The detection of the MMT in the CSAR image is not possible (ρCSAR = -3.7dB) and the false alarms due to the tree trunks are numerous. (b) Compared to CSAR, the MMT response in the SSDSAR image is increased but false alarms due to tree trunks are numerous and have high intensities (ρSSDSAR = -0.5dB). (c) The OBSAR gives the same image as the SSDSAR (ρOBSAR = -0.5dB).

Fig. 14 .

 14 Fig. 14. SAR images for VV polarization using simulated FoPen data with a MMT located at position (108, -1)m. (a) The detection of the MMT in the CSAR image is not possible (ρCSAR = -1.8dB) and the false alarms due to the tree trunks are numerous. (b) Compared to CSAR, the MMT response in the SSDSAR image is increased. The false alarms due to tree trunks remain numerous but have intensities lower than the MMT (ρSSDSAR = 2.8dB). (c) The OBSAR gives the same image as the SSDSAR (ρOBSAR = 2.8dB).

Fig. 15 .

 15 Fig. 15. SAR images for dual polarization (dihedral) using simulated FoPen data with a MMT located at position (108, -1)m.

( a )Fig. 16 .

 a16 Fig. 16. Radar scene of the real data provided by ONERA. A truck located at (5520, 150)m and a trihedral corner reflector located at (5584, 126) are placed in the Nezer forest.

Fig. 17 .

 17 Fig. 17. Singular values for the target matrices in single and dual polarizations (real data). The ranks are chosen to be D H H = D V H = 10 and D + H = D - H = 10.

Fig. 18 .

 18 Fig. 18. Singular values for the interference matrices in single and dual polarizations (real data). The ranks are chosen to be D H J = D V J = 10 and DJ = 10.

Fig. 19 .

 19 Fig. 19. SAR images for VV polarization using real FoPen data. (a) The detection of the two targets in the CSAR image is difficult because there are many false alarms due to the forest. (b) In the SSDSAR image, the two targets are easily detected but false alarms still remain (for the truck ρSSDSAR = 0.8dB and for the trihedral corner reflector ρSSDSAR = 1.5dB). (c) The OBSAR gives the same image as the SSDSAR (for the truck ρOBSAR = 0.8dB and for the trihedral corner reflector ρOBSAR = 1.5dB).

Fig. 20 .

 20 Fig. 20. SAR images for dual polarization (trihedral) using real FoPen data. As we use the polarimetric target subspace H + xy , only the trihedral corner reflector is detected. (a) In the CSAR image, the detection of the target is difficult because there are many false alarms due to the forest. (b) The trihedral corner reflector is easily detected in the SSDSAR image but false alarms due to the forest still remain (ρSSDSAR = 1.5dB). (c) The trihedral corner reflector response is unchanged compared to the SSDSAR and false alarms due to the forest are slightly reduced (ρOBSAR = 2dB).

Fig. 21 .

 21 Fig. 21. SAR images for dual polarization (dihedral) real FoPen data. As we use the polarimetric target subspace H - xy , only the truck is detected. (a) The detection of the target in the CSAR image is difficult because there are many false alarms due to the forest. (b) In the SSDSAR image, the truck is easily detected but false alarms due to the forest still remain (ρSSDSAR = 1.7dB).

  (c) In the OBSAR image, the truck response is unchanged compared to the SSDSAR and false alarms due to the forest are reduced (ρOBSAR = 2.3dB). March 25, 2013 DRAFT

TABLE I

 I List of principal parameters used in the SAR algorithms.

  To summarize, the orthonormal bases of the target subspaces are denoted by H H xy and H V xy for single polarizations and by H + xy and H - xy for dual polarization. The orthonormal bases of the interference subspaces are denoted by J H

	3) Subspace model summary:
	March 25, 2013	DRAFT

TABLE II

 II Accuracy of the target subspace model for single polarization.

	Model	r + xy	r -xy	H + xy	H -xy
	Energy ratio 1% 46%	2%	92%

TABLE III

 III Accuracy of the target subspace model for dual polarization.

We have to choose a sample step for the orientation angle. The value of the sample step is a trade-off between the accuracy of the subspace and the computational cost of its generation. For more details, see[START_REF] Durand | New SAR processor based on matched subspace detector[END_REF].

For the plate, we consider all the orientations. For the interferences, we assume that the trunk tilt γ can not be more than 10 • . March 25, 2013 DRAFT

[START_REF] Forster | Generalized rectification of cross spectral matrices for array of arbitrary geometry[END_REF] DRAFT

ACKNOWLEDGEMENT

The authors thank DGA (french General Delegation for Armaments) for funding this project and ONERA (French Aerospace Lab) for providing us with real data.


G Vitale
email: guido.vitale@ujf-grenoble.fr

L Laforgue

R Michel

J Étienne

V M Laurent

A Duperray

C Verdier

Inverse Problems -from Theory to Applications (IPTA2014) 3D migration of cells solving an inverse problem

Keywords: Traction Force Microscopy, Inverse Problem, Three-dimensional Cell Migration

Traction Force Microscopy (TFM) is an inverse method that allows to obtain the stress field applied by a living cell on the environment on the basis of a pointwise knowledge of the displacement produced by the cell itself during its migration. This biophysical problem, usually addressed in terms of Green functions, can be alternatively tackled in a variational framework. In such a case, a suitable penalty functional has to be minimized. The resulting Euler-Lagrange equations include both the direct problem based on the linear elasticity operator as well as an equation built on its adjoint. Results from a two-dimensional model, i.e. where living cancer cells are migrating on a plane substrate, are briefly presented. While the mathematics is well established also in the three-dimensional case, i.e. where cells are completely embedded in the gel matrix, the experimental data needed are more difficult to obtain than the two-dimensional counterpart. First steps towards the complete three-dimensional traction reconstruction are reported.

Introduction

Many living cells have the ability to migrate, both in physiological and pathological conditions; examples include wound healing, embryonic morphogenesis and the formation of new vessels in tumours. The motility of a cell is driven by the reorganization of its inner structure, the cytoskeleton, according to a complex machinery. The net effect of this process is that a cell is able to apply a stress on the environment, pulling the surrounding material to produce its own movement. Therefore, the determination of forces on the basis of measured displacement, has become a popular problem in the biophysical community.

The early idea to study the force applied by cells in their migration (so called Traction Force Microscopy, TFM) on flat substrate as an inverse problem dates back to the work of Harris and coworkers in the 1980s [START_REF] Harris | Silicone rubber substrata: a new wrinkle in the study of cell locomotion[END_REF]. Afterwards, Dembo and Wang [START_REF] Dembo | Stresses at the cell-to-substrate interface during locomotion of fibroblast[END_REF] came up with a new methodology: the living cell is again put on an elastic substrate and the displacement of fluorescent beads dispersed in such material is evaluated from different images. Finally, they solve the problem of elasticity in the substrate in terms of Green functions and then minimize the error between the measured and the calculated displacement under Tichonov regularization.

An alternative approach can be stated in a continuous variational framework [START_REF] Ambrosi | Cellular traction as an inverse problem[END_REF]. Again, the starting point is a Tichonov functional defined as the displacement error norm plus a regularization. If a variation of such a functional is operated, the definition of an adjoint problem for the unknown force naturally arises. By doing this, two coupled elliptic partial differential equations are obtained and their solution can be approximated, for instance, using finite element.

Content from this work may be used under the terms of the Creative Commons Attribution 3.0 licence. Any further distribution of this work must maintain attribution to the author(s) and the title of the work, book title and publisher. Published under licence by IOP Publishing Ltd.

Although that approach is less common than the seminal method based on Green functions, it has some attractive features that make it worth investigating further: it can be generalized easily to more complex geometry and material behaviour including non-linearities.

In this work we aim to recall the mathematical structure of the TFM problem. We then present some results of the well established two-dimensional theory for a living cancer cell migrating on a plane substrate. Finally, some early data corresponding to an istance of the transmigration of a living cell in a three-dimensional environment are shown.

Mathematical Setting

The Direct Problem of Linear Elasticity

The stationary force balance equations on a linear elastic continuum body Ω ⊂ R 3 writes:

     -div (A∇u) = 0, in Ω, (A∇u) n = T , on Γ N , u = 0, on Γ D . ( 1
)
where u is the displacement field. The domain Ω represents the gel. The boundaries

Γ D , Γ N ⊂ ∂Ω are such that Γ N ∩ Γ D = ∅, Γ N ∪ Γ D = ∂Ω
and n is the unit outward normal of the boundary ∂Ω. On the Neumann boundary Γ N a load is applied per unit surface (here referred as 'traction' and denoted with T ) is applied. The region where the cell and the gel are in contact, say Γ c , is a subset of Γ N . The set Γ c is indeed the support of the traction field T . Last, A is the fourth order Hooke elasticity positive tensor.

For notational convenience, we define the solution operator S, as the map that, for a given control T , assigns the displacement field u that solves the elasticity problem:

u = ST if and only if u solves (1).
Two-dimensional approximations of the elasticity equation have already been succesfully employed in the field of TFM [START_REF] Ambrosi | Cellular traction as an inverse problem[END_REF]. They rely on the existence of suitable averaging operators [START_REF] Michel | Mathematical framework for traction force microscopy[END_REF] whose action is denoted by an overbar (¯). Let us consider Ω being a cylinder with Γ the upper surface, i.e. the plane of cell migration. In [START_REF] Ambrosi | Cellular traction as an inverse problem[END_REF] it is shown that the system (1) can be averaged in the direction orthogonal to Γ N and that gives:

-div Ā∇ū = T , in Γ N , ū = 0, on ∂Γ N . ( 2
)
where Ā denotes the depth-averaged elasticity tensor. Observe that on the Γ N boundary, Dirichlet conditions apply.

Available Experimental Data

Standard imaging techniques are avaible in order to detect the cell, thus giving the region Γ c ⊆ Γ N .

The experimental devices are able to give information on the actual gel displacement field. The position of small fluorescent beads embedded into the gel can be, indeed, tracked for the case of a plane substrate. When dealing with a three-dimensional environment, the position of the intersection of the gel fibers can be similarly obtained. This latter situation represents a novelty from the point of view of data acquisition in TFM and it is currently under development. In view of the smallness of both the beads and the fiber intersection area, we are basically dealing with pointwise measurements of the displacement field. That remark justifies the representation of the observation operator as a list of Dirac mass centered at bead positions. Then, the observation operator writes O := (δ x1 , . . . , δ xN ), where x k denotes the position of the bead labeled with the number k and δ x k is the Dirac mass located at x k .

Biomechanical Constraints

When dealing with applications, the need to include constraints into the model appears. It seems that this was not explicitly stated in the previous literature on TFM, except [START_REF] Vitale | Force traction microscopy: an inverse problem with pointwise observation[END_REF][START_REF] Vitale | A numerical method for the inverse problem of cell traction in 3d[END_REF][START_REF] Michel | Mathematical framework for traction force microscopy[END_REF]. We note, indeed, that Content from this work may be used under the terms of the Creative Commons Attribution 3.0 licence. Any further distribution of this work must maintain attribution to the author(s) and the title of the work, book title and publisher. Published under licence by IOP Publishing Ltd.

• since some reasonably accurate images of the cell are avalaible, the support of the traction field T can be determined, leading to the constraint supp T ⊆ Γ c ;

• the mechanical balance of the cell requires ´Γc T = 0 and ´Γc x × T = 0;

In the following, the space of tractions satisfying the above constraints will be called admissible and denoted with

F adm := {T : Γ N → R d | ˆΓc T = 0, ˆΓc x × T = 0, T = 0 a.e. on Γ N \ Γ c }. (3)
Where d = 2 or 3 depending on the specific application. Note that F adm is a linear subspace of the space of tractions [START_REF] Vitale | Force traction microscopy: an inverse problem with pointwise observation[END_REF]. The projection on the space F adm is denoted by P .

The Inverse Problem

The information provided experimentally to solve the inverse problem, i.e. the pointwise measurements of the displacement u, are not sufficient to yield a unique traction field T . The problem is therefore underdetermined and, as customary, we have to enforce a suitable minimization problem to circumvent this drawback.

Tichonov Functional

The Tichonov functional is defined as:

J (T ) = 1 2 OST -u 0 2 + ε 2 T 2 . ( 4
)
Here u 0 = (u 1 0 , . . . , u N 0 ) are the known displacements at x 1 , . . . , x N respectively. The Tichonov functional J is the sum of two parts: the first term is the discrepancy between the measured displacement u 0 and calculated displacement for a given force T (i.e. ST ), evaluated at the beads location (i.e OST ); the second one is the force magnitude. The two additive contributions are weighted by the positive constant ε, the regularization parameter.

This traction field T is thus defined as the (unique) minimizer of the Tichonov functional J in the set of the admissible tractions F adm . It therefore satisfies:

P J ′ (T ) = εT + P (OS) T (OST -u 0 ) = 0.
(5)

Adjoint Equation

Although the stationarity condition ( 5) is in principle sufficient to define the optimal T , it is convenient to reformulate the problem in terms of a differential equation. It turns useful to define the adjoint state p as

S -T p = O T (Ou -u 0 ), (6)
where S is the solution operator and u is the displacement field, as defined in equation ( 1). Substituting ( 6) into (5), we find the relationship between the optimal force T and the adjoint state:

εT + P p = 0 (7)
The presence of the projection P requires, when solving the problem in practice, the introduction of Lagrange multipliers as done in [START_REF] Vitale | Force traction microscopy: an inverse problem with pointwise observation[END_REF][START_REF] Michel | Mathematical framework for traction force microscopy[END_REF]. Since it can be shown that S is self adjoint [3], equation ( 6) can be rewritten in a more familiar way as follows:

     -div(A∇p) = O T (Ou -u 0 ), in Ω, (A∇p)n = 0, on Γ N , u = 0, on Γ D . (8)
The latter equation is specific for the three dimensional problem. It is easy to deduce the adjoint equation relative to the two-dimensional averaged case [START_REF] Ambrosi | Traction patterns of tumor cells[END_REF], which can be found in [START_REF] Ambrosi | Cellular traction as an inverse problem[END_REF][START_REF] Michel | Mathematical framework for traction force microscopy[END_REF].

Content from this work may be used under the terms of the Creative Commons Attribution 3.0 licence. Any further distribution of this work must maintain attribution to the author(s) and the title of the work, book title and publisher. Published under licence by IOP Publishing Ltd.

Numerical Approximation

Summarizing, our differential problem in strong form is represented by the system of equations ( 1) and ( 8) supplemented by the definition of the adjoint state [START_REF] Legant | Measurement of mechanical traction exerted by cells in three dimensional matrices[END_REF] and the constraints (3). Well-posedness of the aforementioned system and further details are avaible in [START_REF] Vitale | Force traction microscopy: an inverse problem with pointwise observation[END_REF]. The numerical implementation using Finite Elements is thus straightforward. The choice of the regularization parameter ε is made by the use of the L-curve method [START_REF] Engl | Regularization of Inverse Problems[END_REF]. A preliminary validation test in the three-dimensional case can be seen in [START_REF] Vitale | A numerical method for the inverse problem of cell traction in 3d[END_REF].

Our approach differs from [START_REF] Legant | Measurement of mechanical traction exerted by cells in three dimensional matrices[END_REF] where the three dimensional TFM problem is tackled using numerically approximated Green functions. However, no details are provided on the mathematical well-posedness. Also, the latter method does not permit to incorporate the biomechanical constraints mentioned here. In addition, the Green function-based method is proved to be computationally more expensive as compared to ours [START_REF] Vitale | A numerical method for the inverse problem of cell traction in 3d[END_REF].

Applications

The well established tool for evaluating traction forces on a plane substrate [START_REF] Ambrosi | Cellular traction as an inverse problem[END_REF] has been succesfully applied when studying the migration of a living cancer cell [START_REF] Ambrosi | Traction patterns of tumor cells[END_REF]. In Fig. 1a) a living cancer cell is shown while migrating on a plane substrate and is tagged with fluorescence. Fig. 1b) shows the corresponding traction map. Relevant biophysical facts that may be deduced from two-dimensional TFM are detailed in [START_REF] Peschetola | Time-dependent traction force microscopy for cancer cells as a measure of invasiveness[END_REF].

The more difficult case corresponding to migration in a three-dimensional environment is currently under study. An image of the cell embedded in a collagen matrix is shown in Fig. 2a). In this case, the displacement field is measured using an auto-correlation technique and can be reconstructed pointwise. An interpolation has been required in order to visualize the latter as a colormap on a section plane in Fig. 2b), whereas it is not necessary for the numerical calculations. The force reconstruction is, therefore, possible since all the ingredients are available. Together with a discussion of the biological implications of such data, that will be the subject of forthcoming work by our lab.

Figure 1 .

 1 Figure 1. a) GFP actin fluorescent living cancer cell during migration on a plane substrate; b) Traction map.