Anatole Chessel
email: achessel@ifremer.fr

Frederic Cao
email: fcao@irisa.fr

Ronan Fablet
email: rfablet@ifremer.fr

Interpolating Orientation Fields: An Axiomatic Approach

We develop an axiomatic approach of vector field interpolation, which is useful as a feature extraction preprocessing step. Two operators will be singled out: the curvature operator, appearing in the total variation minimisation for image restoration and inpainting/disocclusion, and the Absolutely Minimizing Lipschitz Extension (AMLE), already known as a robust and coherent scalar image interpolation technique if we relax slightly the axioms. Numerical results, using a multiresolution scheme, show that they produce fields in accordance with the human perception of edges.

Introduction

Given a set of edgels (i.e. a set of points with an assigned direction), what are the most invariant and stable ways to reconstruct an orientation field in the whole plan? Because orientation live on the unit circle, an everywhere smooth interpolation is not always possible, due to global topological arguments. However, if we now use local arguments, a analysis similar to [START_REF] Alvarez | Axioms and fondamental equations of image processing[END_REF][START_REF] Caselles | An axiomatic approach to image interpolation[END_REF] is possible and leads also to similar necessary conditions, showing that only very few differential operators have good properties. Since the functions that will be considered in this paper are vector valued or have values in the unit circle, only little is known about existence, uniqueness or classification of the singularities of the solutions to equations we single out. This contribution is an insight of what could be those results and their interest from a low-level vision point of view.

Detecting what we intuitively call "edges" is a first step towards low level feature extraction and integration and has been the focus of a lot of work since the beginning of computer vision. But as noted by psycho-visual experiment and models [START_REF] Kaniza | La grammaire du voir[END_REF], that concept has appeared to be more difficult to define than simply "contrasted image part". Psychovision experiments by the Gestaltists [START_REF] Kaniza | La grammaire du voir[END_REF][START_REF] Wertheimer | Untersuchungen zur Lehre der Gestalt II[END_REF] has given us an acute and unified framework to analyse those effects, and many grouping laws are often involved in the recognition of what we call an edge. The so-called subjective contour effect in particular let us see edges which are strictly speaking not even actually present. It uses amodal completion (reconstruction of occluded edges due to the 2D projection of a 3D world) and modal completion (leading to illusory contour, where the object and the background have the same color). In both cases it rely on a curve interpolation process of unknown data according to the input.

The main origin for those subjective contours is the good continuation principle, which states that if two edgels (edge elements, i.e. a collection of points together with the orientation of the curve which should pass through it) are not too far apart and have compatible directions, we tend to see the curve to which they are both tangent as an edge. Many studies have aimed at computationally implementing this phenomenon. To this end, it is generally assumed that a filter has given us an image of edgels from which we want to extract the curves. Two classical approaches are the curve detector of Parent and Zucker [START_REF] Parent | Trace inference, curvature consistency, and curve detection[END_REF], which uses a discrete co-circularity measure to extract potentially interesting point, and Sha'ashua and Ullman saliency network [START_REF] Sha'ashua | Structual saliency: The detection of globally salient structures using a locally connected network[END_REF], where dynamic programming is used to exhaustively search for the "best" curves under curvature minimisation and length maximisation constraints. More recently, interesting approaches are Medioni's tensor voting [START_REF] Medioni | Grouping .,-,-¿ into regions, curves, and junctions[END_REF][START_REF] Medioni | A computational framework for segmentation and grouping[END_REF], where curves emerge from votes of sparse edgels, and Zweck and al. stochastic completion fields [START_REF] Zweck | Euclidian group invariant computation of stochastic completion fields using shiftable-twistable basis function[END_REF], an Euclidean group invariant implementation of the advection-diffusion model of Mumford [START_REF] Mumford | Elastica and computer vision[END_REF].

Related problems include image inpainting and restoration, and the operators described here are also applied in those cases. In particular, recent developpement extended them to the case of non-scalar image (vector or tensor valued images) [START_REF] Ballester | Disocclusion by joint interpolation of vector fields and gray levels[END_REF][START_REF] Kimmel | Orientation diffusion or how to comb a porcupine[END_REF][START_REF] Osher | Image decomposition and restoration using total variation minimization and the H 1 norm[END_REF]. The aim however is not exactly the same, as this work does not seek to recover the image itself, but an orientation field that would capture its geometrical features.

The good continuation principle states conditions on tangent vectors, and most of the approaches mentioned earlier rely, explicitly or not, on vector or orientation fields. The present work aims at finding out the most invariant interpolation methods based on partial differential equations (PDE). Experiments will be shown using artificial and natural images.

Section 2 states some generalities about interpolating angle, and in particular that singularities are often unavoidable. Section 3 is devoted to the actual axiomatic approach. The last two sections present in more details the two singled out operators, along with experiments. All the proofs are omitted and are given in [START_REF] Chessel | Orientation interpolation : an axiomatic approach[END_REF].

Interpolating Angles

Let Ω ⊂ R 2 and ∂Ω its boundary. Let S 1 be the unit circle of R 2 . We consider the extension problem: knowing I = (I 1 , I 2 ) : ∂Ω → S 1 ⊂ R 2 , how to extend I to the whole domain Ω? The circle represents angular data modulo 2π. In addition, we may also consider the directions of unoriented lines (i.e. angles modulo π). All the argumentation below will also apply to this case.

Topological Restriction

The first problem we encounter when extending vectors as opposed to scalars, is that singularities in the field may be unavoidable. Given a data to be interpolated when can we hope for a singularity free extension? A necessary and sufficient condition is the following. Proposition 1. Let f be a continuous vector field over ∂Ω. There exists a continuous extension of f to Ω if and only if f satisfies condition C.

∃α ∈ S 1 , α / ∈ f (∂Ω), (C)
that is to say only if f is not surjective.

These topological results mean that a singularity free extension is impossible for orientation field when the bounding data cover the whole unit circle. This classical result is equivalent to the Brouwer fixed point theorem [START_REF] Granas | Fixed point theory[END_REF].

A Fundamental Ambiguity

Another problem arising in orientation interpolation (compared to scalar interpolation) is that, due to the periodicity of the data (hence the absence of a total ordering) there are always two manners to interpolate between two fixed values. For instance, on Fig. 1, one goes through zero, while the other goes through For a singularity free (non-surjective) field, at least one value is excluded hence only one of the two possibilities is available. The operators we study are numerically solved by iterative methods, hence, in absence of uniqueness result on the solution of the stationary problem, the final result may depend on the initialisation. Moreover, the topological condition (C) might not be fulfilled in practice for the whole domain but only in sub-domains. This could lead to instabilities that do not exist in the scalar case. A possible workaround is a multiresolution scheme, which will be detailed in a forthcoming work [START_REF] Chessel | Orientation interpolation : an axiomatic approach[END_REF].

Axiomatic Approach

This section details the axiomatic approach exploited to defined operators for the interpolation of orientation field.

Let Γ be a continuous Jordan curve bounding a simply connected domain Ω. We look for an extension operator E, which associates with each directional data θ 0 : Γ → S 1 a unique extension E(Γ, θ 0 ). Throughout all the discussion to follow, it is assumed that θ 0 satisfies the global condition (C). The set of all those functions will be denoted by F(Γ ). It is necessary to parameterise S 1 , to be able to deal with numerical functions. A slight difficulty will arise, since it is not possible to describe the whole circle by a unique chart. Let φ be such a local parameterisation, that is to say, a bijective function U → V , where U is an open subset strictly included in S 1 and V an open subset of R. Let us now consider the extension operator E φ interpolating real valued boundary data u, defined by

E φ (Γ, u) = φ•E(Γ, φ -1 •u). Since E φ (Γ, u) is a numerical function, it
is easier to formulate conditions on the operator E φ . However, since the parameterization φ is arbitrary, the result should be independent on the parametrisation and condition on E φ should stand for any φ. Following [START_REF] Caselles | An axiomatic approach to image interpolation[END_REF], we ask that E φ satisfies the following axioms: Axiom (A1): Comparison principle Let θ 1 , θ 2 ∈ F(Γ ) such that they can be described by a common chart φ. Then φ(θ 1 ) ≥ φ(θ 2 ) implies

E φ (Γ, φ(θ 1 )) ≥ E φ (Γ, φ(θ 2 )).
( 1)

Axiom (A2): Stability principle Let Γ ∈ C, θ 0 ∈ F(Γ ), and Γ ∈ C such that D(Γ ) ⊆ D(Γ ). Then, E(Γ , E(Γ, θ 0 )| Γ ) = E(Γ, θ 0 )| D(Γ ) (2)
Axiom (A3): Regularity principle Let us denote by D(x, r) the disc with center x and radius r. Let Q : R 2 → S 1 such that there exists a parameterisation φ such that

φ(Q)(y) = 1 2 A φ (y -x, y -x) + (p φ , y -x) + c φ
where A φ ∈ SM (2) the set of two dimensional symmetric matrices, p φ ∈ R 2 \{0}, x ∈ R 2 and c φ ∈ R. Then there exists a continuous function

F : SM (2) × R 2 \{0} × R × R 2 , independent of φ such that lim r→0 + φ(E(∂D(x, r), Q| ∂D(x,r) ))(x) -φ(Q)(x) r 2 /2 → F (A φ , p φ , c φ , x). (3)
Axiom (A4): Translation invariance Let τ h θ 0 (x) = θ 0 (x-h), θ 0 : R 2 → S 1 , h ∈ R 2 . Then for all h, E(Γ -h, τ h θ 0 ) = τ h E(Γ, θ 0 ). ( 4
)
Axiom (A5): Domain rotation invariance For any plane rotation R,

E(RΓ, θ 0 • R -1 ) = E(Γ, θ 0 ) • R -1 . (5)
Axiom (A6): Zoom invariance Let H λ θ 0 (x) = θ 0 (λx), for λ > 0. Then,

E(λ -1 Γ, H λ θ 0 ) = H λ E(Γ, θ 0 ). ( 6
)
Once the parametrisation is taken care of, all result obtained in the scalar case are extended to orientation fields. This extension is nearly straightforward, complete proof can be found in [START_REF] Caselles | An axiomatic approach to image interpolation[END_REF].

Theorem 1. Assume that the interpolation operator E satisfies (A1)-(A3). Then F (A, p, x, c) does not depend on c. Moreover, if θ 0 ∈ F(Γ ), then φ(E(Γ, θ 0 )) is a viscosity solution of

F (D 2 u, Du, x) = 0 in D(Γ ) u = φ(θ 0 ) on Γ. ( 7
)
Remark 1. In the scalar case [START_REF] Caselles | An axiomatic approach to image interpolation[END_REF], grey scale shift invariance is assumed to prove this result. Of course, it does not make sense since angles do not add. However, since the result must be invariant with respect to the parameterisation, we get an equivalent property for free.

Theorem 2. Assume that E satisfies axioms ( 1)-( 6) and that F is differentiable at 0. Then, for all parameterisation φ, φ(E(Γ, θ 0 )) is solution of

D 2 u(Du ⊥ , Du ⊥ ) = 0 in D(Γ ). u = φ(θ 0 ) on Γ. (8)
Remark that this operator is the curvature of the level lines of u, up to a |Du| 3 factor. These level lines are independent of the parameterisation, which makes the result possible. Indeed, the independence with respect to parameterisation implies that, for all admissible φ and ψ,

E(Γ, θ 0 ) = φ -1 • E φ (Γ, φ • θ 0 ) = ψ -1 • E ψ (Γ, ψ • θ 0 ).
By noting u = φ • θ 0 and g = ψ • φ -1 , this equation becomes

g • E φ (Γ, u) = E ψ (Γ, g • u).
This condition is closely related to invariance with respect to contrast change for scalar data, and the arguments developed in [START_REF] Alvarez | Axioms and fondamental equations of image processing[END_REF] indeed apply.

As noted in disocclusion experiments [START_REF] Masnou | Level-line based disocclusion[END_REF], this operator interpolates the level lines of the data with straight lines. A well known problem is that the solution of this equation may be not unique, and we will see in the next section that if it manages to keep the discontinuities structuring the image, it fails to give a field smooth enough to recover subjective contour. Thus we may drop the full independence over parameterisation and slightly relax Axioms (1) and (3).

Proposition 2. Assume that Axioms (1) and (3) only holds for parameterisation that are Euclidean, up to a multiplicative factor. Then φ(E(Γ, θ 0 )) is solution of

aD 2 u(Du, Du) + bD 2 u(Du, Du ⊥ ) + cD 2 u(Du ⊥ , Du ⊥ ) = 0, (9)
where ac -b 2 ≥ 0.

The condition ac -b 2 ≥ 0 ensures that the equation is elliptic, and that the maximum principle can hold. As expected, a solution of ( 9) is invariant with respect to an affine reparameterisation of the circle, but not to any general parameterisation. Among all those operators, the case b = c = 0 is the Absolutely Minimizing Lipschitz Extension (AMLE)

D 2 u(Du, Du) = 0 in D(Γ ), u| Γ = φ(θ 0 ) on Γ, (10)
for which existence and uniqueness of viscosity solution is known. It gives continuous oscillation free solution. It will be studied in more detail in Sect. 5.

Angle Interpolation with the Curvature Operator

As a result of the previous section, the only operator satisfying the given axioms is the curvature operator. It is well known in the computer vision community as a scalar restoration operator via total variation minimisation and has been used for scalar interpolation to solve the disocclusion problem [START_REF] Masnou | Level-line based disocclusion[END_REF][START_REF] Rudin | Nonlinear total variation bsaed noise removal algorithms[END_REF][START_REF] Chan | Local inpainting model and TV inpainting[END_REF]. The argumentation above gives the equation which is locally satisfied by the orientation of the vector field. An alternate formulation [START_REF] Tang | Diffusion of general data on non-flat manifolds[END_REF] is to consider the variational problem min |u| so that the constraint is automatically satisfied. Careful calculations lead to a system of the two coupled PDEs

div DI p-2 DI i + I i DI p = 0, 1 ≤ i ≤ 2. ( 11
)
It is worth noticing that DI p may be interpreted as the Lagrange multiplier of the constraint |I| = 1. The case p = 1, corresponding to the total variation, leads to

div DI i DI + I i DI = 0, 1 ≤ i ≤ 2. ( 12
)
As a sanity check, elementary calculations lead to the following result, which holds thanks to the particular choice of norm DI .

Proposition 3. Let I = (I 1 , I 2 ) ∈ C 2 (Ω, R 2 ) with |I| = 1 everywhere. Let θ such that I = (cos(θ), sin(θ)). Then

div DI i DI + I i DI = 0, 1 ≤ i ≤ 2. ⇐⇒ 1 |Dθ| D 2 θ Dθ ⊥ |Dθ| , Dθ ⊥ |Dθ| = 0. ( 13
)

Numerical Resolution

Experiments (not detailed here) show that the numerical solutions of the two outlined methods are indeed the same. Hence in the following we use the parametrized to S 1 equation. To compute a numerical solution, we will use the associated evolution problem

∂θ ∂t = D 2 θ Dθ ⊥ |Dθ| , Dθ ⊥ |Dθ| , (14)
and let t → ∞.

To solve this equation, we used a non-linear over relaxation scheme (NLOR) similar to the one found in [START_REF] Caselles | An axiomatic approach to image interpolation[END_REF] implemented in Megawave2 [START_REF] Froment | Megawave2[END_REF]. As noted in [START_REF] Perona | Orientation diffusion[END_REF] and [START_REF] Cecil | Numerical methods for minimization problems constrained to S 1 and S 2[END_REF], we took some special numerical attention to work with angles modulo 2π and circumvent the problem caused by the discontinuities at 2kπ, k ∈ N. Moreover, as noted above, a multiresolution scheme was used to take care of the initialisation, see [START_REF] Chessel | Orientation interpolation : an axiomatic approach[END_REF].

Experiments on the curvature operator

Geometrically, the curvature extension operator tries to extend the level lines of the boundary data by straight lines. Obviously, there are cases for which that approach does not apply [START_REF] Masnou | Level-line based disocclusion[END_REF]. In particular it fails to compute any solution for the simple artificial cases we will see in the next section (see figure 3).

However, experiments carried out for larger images with a larger set of boundary points yields interesting results. Figure 2, displays an example with the Lena image. We visualise the field via its field line, using Line Integral Convolution (LIC) [START_REF] Cabral | Imaging vector field using line integral convolution[END_REF]. The initial field is given by the orientation of the tangents to the level lines (the orthogonal to the gradient) decimated with a thresholded Canny-Deriche edge filter [START_REF] Deriche | Using canny's criteria to derive a recursively implemented optimal edge detector[END_REF]. As expected, the curvature operator keeps discontinuities, as at the top of Lena's hat. Interestingly enough, it also manages to keep singularities adequately. In particular, T-junctions are preserved, which is particularly relevant in a perceptual point of view. Not only singularities that are present on the boundary data are preserved, but they can also be created in the interpolated area in a suitable way (see for instance at the interface of the cheek and the hair).

To sum up, the curvature operator is (as in the scalar case) able to preserve singularities when necessary. It may be considered as a drawback when the smoothest solution is sought. Moreover, there is no existence and uniqueness result in the general case.

AMLE on angle

In this section, we will provide more insight to the AMLE extension. A more detailed presentation can be found in [START_REF] Caselles | An axiomatic approach to image interpolation[END_REF][START_REF] Aronsson | A tour of the theory of absolutely minimizing functions[END_REF][START_REF] Jensen | Uniqueness of lipschitz extentions: minimizing the sup norm of gradient[END_REF] We know that a non surjective data can be smoothly interpolated inside a single parameterisation, and that AMLE is independent of affine change of parametrisation.

AMLE was introduced in [START_REF] Aronsson | Extention of function satisfying lipschtiz conditions[END_REF]. It was proved (see [START_REF] Aronsson | A tour of the theory of absolutely minimizing functions[END_REF][START_REF] Jensen | Uniqueness of lipschitz extentions: minimizing the sup norm of gradient[END_REF] and references therein) that it can be equivalently defined, in the scalar case, as -the extension in Ω of a data defined on ∂Ω whose Lipschitz constant is minimal in any Ω ⊂ Ω. -the viscosity solution of the PDE D 2 u(Du, Du) = 0.

-the limit for p → ∞ of p-harmonic maps, defined as the minimization of the p-harmonic energy min

W 1,p (Ω) |Du| p .
Those results heavily rely on a maximum principle (eventually proved by Jensen [START_REF] Jensen | Uniqueness of lipschitz extentions: minimizing the sup norm of gradient[END_REF]), which guarantees in particular that the solution has no oscillation inside the domain. More importantly, it yields the existence and uniqueness of the solution.

Equivalence of AMLE on Angle and AMLE on Vector

Restricted to S 1

Again, we can link the intrinsic formulation on angle used until now and the R 2 restricted to S 1 one. Let us consider [START_REF] Ballester | Disocclusion by joint interpolation of vector fields and gray levels[END_REF] again and let p go to +∞. We formally obtain the two coupled equations

2 i=1 D 2 I i (DI i , DI j ) = 0 j = 1, 2. ( 15
)
The definition of a solution of this system is, to the best of our knowledge, an open problem. However, we point out the two following interesting facts.

This means that I is a vector AMLE on the circle if its argument is a scalar AMLE.

The second point is that the term corresponding to the constraint |I| = 1 has vanished from [START_REF] Ballester | Disocclusion by joint interpolation of vector fields and gray levels[END_REF] to [START_REF] Granas | Fixed point theory[END_REF]. Now, a method to solve the stationary problem [START_REF] Granas | Fixed point theory[END_REF] is to solve the corresponding evolution system

∂I j ∂t = 2 i=1 D 2 I i (DI i , DI j ) j = 1, 2. ( 17
)
If I is a continuous solution of [START_REF] Rudin | Nonlinear total variation bsaed noise removal algorithms[END_REF] such that |I| = 1 everywhere at time t = 0, does it remain true for t > 0? At this step, we cannot tell, but we have the following hint.

Lemma 1. Let I be a a C 2 vector field with |I| = 1 everywhere. Then the vector with coordinates 2 i=1 D 2 I i (DI i , DI j ) is everywhere normal to I.

Numerical Resolution

Again we can numerically check for the equivalence of the two formulation, and we used the simpler parametrized to S 1 one. It is proved, in the scalar case, that when t → ∞ the solution of the evolution problem

∂θ ∂t = D 2 θ(Dθ, Dθ), (18)
tends to the solution of the stationary problem, because this solution is unique.

In the vector case, we do not have such a result, but we display experiments showing that this is still reasonable. We used a multiresolution NLOR scheme similar to the one used for the curvature operator (see Sect. 4.1).

Experiments

Figure 3 shows numerical results on artificial data. The first one simply consists of two vectors. The interpolated vector field is as expected tangent to the curve with which we would like to connect the two vectors, something close to Euler elastica [START_REF] Mumford | Elastica and computer vision[END_REF]. The next two figures show the same mechanism with more complex curves: a circle and a tube. The interpolated field is perceptually sound. As asserted by Prop. 1, we do find singularities in the center of the circle and the extension is there somewhat chaotic, as we are looking for a Lipschitz function where it cannot even be continuous. The situation below the tube (Fig. 3, last experiment) is interesting as it is an example of the ambiguity of Sect. 2.2: a smooth extension do exist, but due to the lack of information the algorithm extended the orientation field the other way round and put a singularity.

We see Fig. 4 an experiment on Lena. The initial field is again the orientation of the tangents to the level lines (the orthogonal to the gradient) decimated with a thresholded Canny-Deriche edge filter [START_REF] Deriche | Using canny's criteria to derive a recursively implemented optimal edge detector[END_REF]. The interpolation field is again tangent to the edges as requested. On the other hand, there is no control on the position of the unavoidable singularities. Moreover, singularities are smoothed out, which can be expected, regarding the properties of the AMLE in the scalar case.

Nonetheless, the AMLE is a good candidate for an interpolation operator as we have a complete theory in the scalar case stating existence and uniqueness of solution. Moreover, it gives smooth solution from which extracting subjective contour as curves is possible. Compared to the curvature however, it tend to lack the ability to keep discontinuities in the fields it produce. Fig. 3. Three artificial geometric tests, initial orientation field on the left, AMLE extension visualised with LIC on the right. As expected, where we would put a curve it do find the tangent to that curve. It's behaviour is less predictable where there is no information.

Fig. 4. Test with the Lena image, initialised with the orthogonal to the gradient orientation field decimated using a Canny-Deriche filter. Notice that the recovered field is tangent to the edges, in particular at the top of the hat, on the strands of hair around the face and on Lena's jaw and chin.

Conclusion

An axiomatic approach of orientation field interpolation has been presented to define extension operators. There is a unique operator satisfying a small set of axioms including geometrical invariance and stability: the curvature operator. This operator is able to preserve singularities. On the other hand, one may require a smoother solution. Moreover, an existence and uniqueness of a solution are not well established. If the independence of the interpolation with respect to reparameterisation of the unit circle is relaxed, another operator becomes interesting: the AMLE. Existence and uniqueness holds in the scalar case. The AMLE is, to some extent, dual to the curvature operator (it minimizes the L ∞ norm of the gradient, while the curvature minimizes the norm L 1 ), and somehow smoothes out the singularities.

The proposed extension operators provides the required basis for the extraction of meaningful curves in images as curves tangent to the orientation field: for instance using Fast Marching approaches [START_REF] Sethian | Level Set Methods and Fast Marching Methods[END_REF]. Other obvious applications include LIC-based interpolation or restoration [START_REF] Cabral | Imaging vector field using line integral convolution[END_REF][START_REF] Tschumperlé | LIC-based regularization of multi-valued images[END_REF]. The operators above are the more natural popping out from the required axioms. However, if some of them are relaxed or more prior knowledge from the image is introduced, some variations of these operators may lead to new types of interpolation model.

Fig. 1 .

 1 Fig. 1. Ambiguity of interpolation between two directions (modulo π in the present figure). Left: there are two ways of going from u to v. Right: example when u -v = π/2

W 1 ,

 1 p (Ω) DI p , under the constraint |I| = 1. In this case I = (I 1 , I 2 ) and |I| is the Euclidean norm |I| = I 2 1 + I 2 2 and DI = |DI 1 | 2 + |DI 2 | 2 . Inspired by the scalar case, we can compute the Euler-Lagrange equations for the energy above by setting I = u

Fig. 2 .

 2 Fig. 2. Interpolation with the curvature operator on Lena. The initialization is orthogonal to the gradient orientation field decimated using a Canny-Deriche filter. A general observation is that T-junctions are preserved.

Proposition 4 . 2 i=1D 2

 422 Let I = (I 1 , I 2 ) ∈ C 2 (Ω, R 2 ) with |I| = 1 everywhere. Let θ such that I = (cos(θ), sin(θ)). Then I i (DI i , DI j ) = 0 j = 1, 2 ⇐⇒ D 2 θ(Dθ, Dθ) = 0.


Rafail Psiakis

Angeliki Kritikakou
email: angeliki.kritikakou@irisa.fr

Olivier Sentieys

Emmanuel Casseau
email: emmanuel.casseau@irisa.fr

Run-time Coarse-Grained Hardware Mitigation for Multiple Faults on VLIW Processors

As transistors scale down, processors are more vulnerable to radiation that can cause multiple transient faults in function units. Rather than excluding these units from execution, performance overhead of VLIW processors can be reduced when fault-free components of these affected units are still used. In the proposed approach, the function units are enhanced with coarse-grained fault detectors. A re-scheduling of the instructions is performed at run-time to use not only the healthy function units, but also the fault-free components of the faulty function units. The scheduling window of the proposed mechanism is two instruction bundles being able to explore mitigation solutions in the current and the next instruction execution. Experiments show that the proposed approach can mitigate a large number of faults with low performance and area overheads.

I. INTRODUCTION

The energy transferred by radiation is known to cause transient faults on system hardware, which may severely affect the application execution. As technology size decreases, the hardware becomes more and more susceptible to radiation. Although multiple simultaneous faults have been neglected for a long time, they can no longer be negligible for technologies of 130nm and beyond [START_REF] Heijmen | Radiation-induced soft errors in digital circuits -a literature survey[END_REF]. A single particle hitting the silicon creates secondary particles, which can be emitted in several directions, and, thus, affect different nodes of a circuit [START_REF] Ferlet-Cavrois | Statistical Analysis of the Charge Collected in SOI and Bulk Devices Under Heavy lon and Proton Irradiation -Implications for Digital SETs[END_REF].

Approaches to deal with multiple faults either apply hardware redundancy or software redundancy. Hardware redundancy inserts spare resources to the original processor to execute in parallel the same instruction and to compare the obtained results [START_REF] Klecka | Error self-checking and recovery using lock-step processor pair architecture[END_REF]. Small performance overhead is usually observed (due to the need of comparison), whereas the area overhead is significant. Software redundancy modifies the program by inserting redundant instructions to be executed on the original processor. Although the area is not increased, the impact on the execution time is significant due to the execution of redundant and comparison instructions [START_REF] Reis | SWIFT: Software Implemented Fault Tolerance[END_REF]. However, when the system processor has parallel Function Units (FUs), the idle FUs can be used to execute these instructions. Such type of processors are the Very Long Instruction Word (VLIW) processors, which usually consist of complex FUs able to execute all types of operations and simpler ones that cannot execute sophisticated operations, such as multiplications and divisions. Existing approaches on VLIW processors that detect the faulty units before the application execution, such as schedule multi-versioning [START_REF] Karri | Computer aided design of fault-tolerant application specific programmable processors[END_REF] and instruction modification [START_REF] Scholzel | Combining hardware-and software-based self-repair methods for statically scheduled data paths[END_REF], cannot be applied for soft errors. Few approaches detect the faulty FUs during execution. However, instruction duplication and re-execution approaches [START_REF] Sartor | Exploiting idle hardware to provide low overhead fault tolerance for vliw processors[END_REF] and run-time rescheduling mechanisms [START_REF] Psiakis | Neda: Nop exploitation with dependency awareness for re-liable vliw processors[END_REF] have not been explored for multiple errors. The work in [START_REF] Psiakis | Run-time instruction replication for permanent and soft error mitigation in vliw processors[END_REF] applies duplication and triplication of the instructions for multiple errors, whereas a mechanism is proposed to exclude the faulty FUs permanently for the rest of the execution. However, the negative impact on performance is significant.

To reduce this performance degradation while keeping the area overhead low, we propose a coarse-grained hardware mechanism that detects the faults during execution and excludes only the faulty components of the affected FUs. Our main contributions are: i) extension of the FUs with Built-In Current Sensors (BICS) [START_REF] Viera | Validation of Single BBICS Architecture in Detecting Multiple Faults[END_REF], which monitor the induced transient currents to detect a fault; ii) the run-time instruction scheduling that excludes only the faulty components of the affected FUs exploring mitigation solutions in the current and the next instruction execution; and iii) extensive fault injection simulations varying the number of total and concurrently occurring faults. The obtained results show that the proposed approach mitigates several multiple faults with low overheads.

The organization of this paper is as follows. Section II presents an overview of the approach and the architecture of the proposed coarse-grained mechanism is detailed. Section III presents the experimental results and Section IV concludes this work.

II. PROPOSED RUN-TIME MITIGATION MECHANISM

We use the 4-issue heterogeneous VLIW data-path of Fig. 1 to schematically illustrate our approach. In this example, the VLIW consists of a 3-stage pipeline with Fetch (F), Decode (DC) and Execute/Memory/Write-Back (EX/MEM/WB). A number of instructions is formatted as one big instruction, named as instruction bundle, which is executed in parallel by the FUs of the processor. The first issue is configured with 1 Arithmetic Logic Unit (ALU ) and 1 Branch unit (BR), the second issue with 1 ALU and 1 Memory unit (MEM ), and the third and fourth issues with 1 ALU and 1 Multiplication unit (MUL). The components in gray color correspond to the basic architecture and the components with green color are the extra hardware components required by our approach. The proposed approach focuses on multiple soft faults occurring in the arithmetic FUs, as they have the largest area footprint of the combinatorial components (see section III.A)). The faults in the storage components, e.g. register file, memory and pipeline registers, are assumed to be protected with other methods, such as Error Correction Codes (ECC).

DC[0] DC[0] DC[0] DC[0] 0 0 DC[0] 0 DC[0] 0 0 0 0 0 s[0] 1 NOP s[0] 1 NOP s[0] 1 NOP s[0] 1 NOP s[0] 1 NOP s[0] 1 1 NOP NOP DC[1] DC[1] DC[1] DC[1] DC[1] 0 0 DC[1] 0 0 0 0 0 1 s[1] NOP 1 s[1] NOP 1 s[1] NOP 1 s[1] NOP 1 s[1] NOP 1 1 s[1] NOP 1 1 DC[2] DC[2] DC[2] DC[2] DC[2] 0 0 DC[2] 0 0 0 0 0 1 NOP 1 s[2] NOP 1 s[2] NOP 1 s[2] NOP 1 s[2] NOP 1 s[2] NOP 1 1 s[2] DC[3] DC[3] DC[3] 0 0 DC[3] 0 0 DC[3] 0 DC[3] 0 0 0 0 NOP s[3] 1 NOP s[3] 1 NOP s[3] 1 NOP s[3] 1 NOP s[3] 1 NOP s[3]

A. Overview and motivation

We will use the following example to present the idea of our mechanism. Fig. 2(a) depicts the original schedule. The instruction bundle is made of 3 additions (ADD 1 , ADD 2 , ADD 3 ) and 1 multiplication (M U L 1 ). To be able to detect an error, existing software approaches require to duplicate the instructions, such as in Fig. 2(b) and compare the results, and, thus, insert additional time slots increasing the application execution time. Existing hardware approaches are applied for permanent errors, whereas error detection is performed before the application execution, such as [START_REF] Scholzel | Combining hardware-and software-based self-repair methods for statically scheduled data paths[END_REF]. In case an error is detected, for instance if a fault occurs in the MUL unit of the third issue, hardware approaches insert an additional time slot and re-execute the instruction scheduled in the third issue to another compatible FU of another issue, as depicted in Fig. 2(c).

To reduce the negative impact on performance of existing approaches, we propose a mechanism that removes the need of replicated instructions and explores at run-time the rescheduling of the faulty instructions to explore the healthy and idle FUs of the current bundle and the upcoming bundle, as depicted in Fig. 2(d). In this example, the instruction ADD 3 is moved to the third issue whose ALU is still healthy, whereas the instruction M U L 1 is moved to fourth issue. Therefore, no need exists for an additional time slot in this example.

B. Error detector

The error detector keeps the faulty status of the FUs and takes care of miscalculated results. Each FU is internally

ADD 1 ADD 2 MUL 1 ADD 3 a) ALU BR ALU MEM ALU MUL ALU MUL ADD 1 ADD 2 MUL 1 ADD 3 b) ADD 1 ADD 2 ADD 3 NOP NOP MUL 1 c) ADD 1 ADD 2 ADD 3 MUL 1 d) t i t i+1 ADD 1 ADD 2 MUL 1 ADD 3 t i t i+1 t i t i
Fig. 2: Example: Motivation and overview of the approach.

enhanced with BICS sensors [START_REF] Viera | Validation of Single BBICS Architecture in Detecting Multiple Faults[END_REF], [START_REF] Bastos | Comparison of bulk built-in current sensors in terms of transient-fault detection sensitivity[END_REF]. A BICS is attached to a group of transistors. During normal operation, the current in the bulk of these transistors is approximately zero. Only the leakage current flows through the biased junction, which is still very low compared to the current generated by energetic particles. When an energetic particle generates a current in the bulk, the bulk-BICS captures that a transient fault occurs. The bulk-BICS has a reset mechanism that allows the fault detection to be active only as long as it takes to dissipate the transitory energy pulse. When the fault has vanished, the affected transistors can be used once again. A signal from the BICS from each FU of one issue is combined into a status signal, s[i] in Fig. 1, with a size equalling to the number of FUs of the issue i. Hence, if a bit in the signal s[i] is set, it means the corresponding FUs is affected by a fault.

In case of one or more faults at cycle k -1, the results of the corresponding instruction(s) of the execution stage are miscalculated, and, thus, they must not be committed. For this purpose, each VLIW issue is enhanced with a multiplexer controlled by the signal sel_WB (with a size equalling to the number of VLIW issues) computed by the error detector. When a bit in sel_WB is set, the corresponding multiplexer passes a NOP result (instead of the miscalculated result) and the WB and MEM enable of the corresponding issue is disabled.

C. Mitigation switch

At cycle k, the DC unsched shadow register stores the decoded instructions that couldn't be scheduled and the EX sched shadow register keeps the instructions executed but in which an error occurred during their execution. Therefore, at cycle k the instructions to be scheduled can potentially come from three inputs: 1) the decoded instructions at cycle k -1 (DC) 2) the remaining instructions not scheduled at cycle k-1 (DC unsched register) and 3) the executed instructions with a fault at cycle k -1 (EX sched register). A mitigation switch is required to move at run-time the instructions in different issues than the ones defined by the compiler's original schedule at compile time. We reduce the complexity of the mitigation switch by performing hierarchical decisions. The first decision is whether a shadow register must be used or the currently decoded instructions and the second decision is which shadow register should be used (EX sched or DC unsched).

D. Error Mitigator

The error mitigator extracts from the fetch stage the required information and the dependencies between bundles and decides which of the three potential inputs that provide instructions for execution (DC, DC unsched, and EX sched) will be scheduled based on the status of the FUs and the type of the instructions.

The mitigator initially performs an early decoding of the instructions at the fetch stage to find the opcode, the destination registers and the source registers of each instruction and to identify instruction dependencies between two consecutive bundles.

To perform the mitigation process, each potential instruction input to the mitigation switch is modeled by a bit mask. A bit mask corresponds to the instruction scheduled at position i, as depicted in Fig. 3. a) bit 3: when this bit is set, the FU component required for the instruction execution is an ALU FU; b) bit 2: when this bit is set, the FU component required for the instruction execution is a MUL FU; c) bit 1: when this bit is set, it means the instruction is an instruction from a previous bundle that has not been scheduled yet; and d) bit 0: when this bit is set, it means the instruction has a dependency with at least one of the instructions of the next bundle. Then, each of the three potential inputs of the switch is represented by an array of bit masks (EX sched mask, DC unsched mask and DC mask) that has a size equalling to the number of VLIW issues.

The status of all FUs components is represented by the array status, where each element is the s signal of an issue enhanced with an additional bit that is set when the issue is occupied as shown in Fig. 3 (bit 0). (S3): If no new fault occurred at cycle k -1, the mechanism schedules first the remaining decoded instructions from cycle k -2 for execution at cycle k. To do so, the array DC unsched mask is used as input to the mitigation process. Then, the current decoded instructions DC mask are used as input to the mitigation process.

DC_mask[i] DC_unsched_mask[i] EX_sched_mask[i]

(S4): If there are still remaining instructions from the previous bundle and/or if there is any dependent instruction in current decoded instructions that cannot be scheduled, the stall signal is activated to stall the fetch and the decode stage for one cycle, so as to schedule these instructions.

The inputs of the mitigation process are: 1) the mask of the input signal to be scheduled, input; and 2) the array that represents the status of the FU components, status. The outputs are: 1) the control signal sel switch of the mitigation switch; 2) the control signal sel PB of the multiplexer between the shadow registers; 3) the updated mask of the input signal; and 4) the updated status. The procedure is given by Alg. 1. For all the instructions i of the input (line 3) and for each issue j described by status (line 4), if the instruction i has not been scheduled and the j issue is unoccupied (line 5), check if the required FU is available (line 6). If this is true, the occupied bit of the corresponding status is set (line 7), the remaining bit of the input is cleared, since the instruction is scheduled (line 8), and the signal sel switch instructs the switch to pass the instruction currently at issue i to issue j (line 9) and the next instruction is explored.

Algorithm 1 Mitigation process 1: Inputs: input, status 2: Outputs: input, status, sel_switch 3: for i ∈ {0, n} do for each instruction 4:

for j ∈ {0, n} do for each issue if the instruction is not scheduled and the issue is free

5: if (input[i][1]&&status[j][0]) then if the required FU is healthy 6: if (input[i][3]&&status[j][2] || input[i][2]&&status[j][1]
) then The issue is occupied end for 14: end for

III. EVALUATION RESULTS

The original processor and the processor enhanced with the proposed approach have been developed in C++ and synthesized using the Catapult High Level Synthesis (HLS) tool to obtain the RTL design. The gate-level netlist was generated by the Design Compiler tool from Synopsys using 28 nm ASIC library. To evaluate our approach, we use ten benchmarks from the MediaBench suite [START_REF] Lee | Mediabench: a tool for evaluating and synthesizing multimedia and communications systems[END_REF]. The benchmarks are compiled with VEX compiler [START_REF]VEX Toolchain[END_REF] for each configuration. For the evaluation results, we used the VEX VLIW processor architecture [START_REF] Fisher | Embedded computing: a VLIW approach to architecture, compilers and tools[END_REF] with two different configurations:

• 4-issue VLIW configured with 4 ALU FUs, 2 MUL FUs, 1 MEM FU and 1 BR FU, • 8-issue VLIW configured with 8 ALU FUs, 4 MUL FUs, 2 MEM FUs and 1 BR FU.

For every benchmark, Table I shows the average number of useful instructions per bundle, ILP, for both the 4-issue and 8issue VLIW configurations. Because ILP inside bundles is low compared to the number of issues, idle issues of the current or next bundles can be used to execute the faulty instructions in addition to switching instructions with healthy components of units affected by a fault.

A. Area analysis of the unprotected VLIW processor

First of all, we present the logic area of each unit of the original unprotected VLIW architecture (as depicted in gray color in Fig. 1) in Table II and the area of each pipeline stage for the 4-issue and 8-issue VLIW in Table III. The results motivate the focus of the proposed approach on the execute stage of the VLIW processor, since it covers more area, and, thus, it is more likely to be exposed to faults.

B. Performance

In this paper, we present first experiments to evaluate the benefit of our approach based on coarse-grained FU exploration. The fault model is simplified as much as possible : we randomly injected multiple faults during the benchmarks' execution and, in order to reduce the number of scenarios, we consider the faults as persistent after they occur, i.e. they last for the rest of the execution. Performance overhead is thus the upper bound compared to a transient fault based simulation. Each benchmark is tested for: i) 0 up to 4 multiple faults for the 4-issue configuration; and ii) 0 up to 10 multiple faults for the 8-issue configuration 1 .

Table IV and Table V show the number of execution cycles required to execute the benchmarks considering :

• the original application on the non-protected VLIW (Original), • the proposed run-time coarse-grained mitigation approach, i.e a time and space idle slot exploitation method in addition to current sensors to use healthy components of faulty function units • the version where only detection can occur through time and space idle slot exploitation with duplication of the instructions during the execution (2I), • the version where correction and masking occur through time and space idle slot exploitation with triplication of the instructions during the execution (3I), • the approach where correction and masking occur by triplicating the FUs in each issue with respect to the original unprotected processor while running the original application (3FUs). The performance results are obtained by taking the mean value of 20 simulations running the same benchmark and the faults are injected at random execution cycles for each simulation. When no faults occur, the proposed approach have the same performance as the unprotected version, i.e. the number of execution cycles is the same as the original execution. When detection and correction occurs by triplicating the FUs in the issues, the number of execution cycles is the same as the original version one (not presented in Table IV andTable V).

From the Table IV, we observe that the average performance overhead of the proposed approach for the 4-issue configuration is 0% for zero fault up to 83% for four faults. For the 2I (3I) approach, the average performance overhead is 72% (152%) for zero fault up to 262% (439%) for four faults. From the Table V, we observe that the average overhead of the proposed approach for the 8-issue configuration is 0% for zero fault up to 141% for ten faults. For the 2I (3I) approach the average overhead is 38% (104%) for zero fault up to 378% (626%) for ten faults.

Overall, we observe that: 1) the proposed approach inserts significantly lower performance overhead than the approaches with instruction duplication or triplication and 2) in several benchmarks, when the number of faults is low, our approach's performance is very close to the original one, i.e. without faults. The obtained gain of the proposed approach comes from the fact that whenever a persistent fault is detected, the proposed hardware mechanism exploits healthy FUs in the current and the next bundle execution in addition to idle FUs.

C. Area and delay overhead

Table VI presents the impact in area and in delay of each component of the proposed run-time coarse-grained mechanism for the 4-issue VLIW. It should be stressed that the delay introduced by the error mitigator does not affect the clock of the VLIW processor, since it is executed in parallel with the VLIW data-path. Table VII shows the area overhead of the proposed mechanism, the approach that provides error masking by applying triplication of instructions and the approach that triplicates the FUs at each issue. Compared to the unprotected version, the proposed approach implies an area overhead of 18% and 28% for the 4-issue and for the 8-issue configurations respectively. For the 3I approach the area overhead is 24% for the 4issue configuration while for the 8-issue configurations it is 30%. It is quite similar to the area overhead of the proposed approach. However the proposed approach achieves much better performance as seen in previous section. For the 3FUs approach the area overhead is 45% for the 4-issue while for the 8-issue configurations it is 56%.

IV. CONCLUSION

A hardware mechanism is proposed for multiple soft faults that characterizes the FUs of the VLIW processor and reschedules at run-time the instructions, exploiting in time and space idle FUs and healthy components of faulty function units so as to mitigate faulty instructions. From the obtained results, it is shown that several multiple faults can be mitigated with significant reduction in the performance and area overheads compared with approaches that duplicate and triplicate the instructions during execution or triplicate the FUs in the issues. Work in progress is about experiments with faults that last various durations.

  [0] DC/EX[0] DC/EX[0] DC/EX[0] DC/EX[0] DC/EX[1] DC/EX[1] DC/EX[1

Fig. 1 :

 1 Fig. 1: Extension of the proposed mechanism on the original VLIW (4-issue VLIW example).

Fig. 3 : 1 C0Fig. 4 :

 314 Fig. 3: Bit masks to model the inputs and FU component status of the mitigation switch.

TABLE I :

 I Benchmark characteristics

	Benchmark ILP 4-issues ILP 8-issues
	adpcm dec	1.77	2.28
	adpcm enc	1.82	2.41
	bcnt	2.49	3.62
	fft32x32s	2.85	4.19
	huff	1.51	1.75
	motion	1.94	2.39
	dct	2.22	3.31
	crc	1.76	1.8
	fir	2.09	2.5
	matrix mul	2.61	4.46

TABLE II :

 II Area of main VLIW combinatorial units (µm 2 ).

	DC	BR	ALU FU MUL FU	MEM/WB
	250	2,290	1,533	2,310	358

TABLE III :

 III Area of VLIW stages (µm 2 ).

	Stage	4-issues 8-issues
	DC	3,290	6,290
	EX	10,752	21,5094
	MEM/WB	1,432	2,864

TABLE VI :

 VI Area and delay overheads introduced by the the components of the proposed mechanism for 4-issue VLIW.

	Component	area	delay
		(µm 2 )	(ns)
	Mitigation Switch	3,918	0.16
	Error mitigator	2,112	2.44
	Error detector	3,268	0.84

TABLE VII :

 VII VLIW processor area footprint.

	Approach	4-issues area area	8-issues area area
		(µm 2 )	(%)	(µm 2 )	(%)
	Original	50,844	-	79,661	-
	Instr. triplication (3I)	62,812	23.54	103,598	30.05
	FUs triplication (3FUs)	73,523	44.60	124,137	55.83
	Proposed mechanism	60,143	18.29	102,201	28.3

TABLE IV :

 IV Performance comparison for 4-issue VLIW

	Benchmarks Original	0	Proposed mechanism 1 2 3	4	0	Instr. duplication (2I) 1 2 3	4	0	Instr. triplication (3I) 1 2 3	4
	adpcm dec	386	386	412	437	538	565	676	731	809	837	1124	996	1150	1153	1605	1682
	adpcm enc	409	409	462	469	503	611	722	831	862	902	1210	1062	1094	1230	1387	1817
	bcnt	478	478	511	512	784	1149	870	904	904	1199	2308	1264	1297	1297	1989	3364
	fft32x32	569	569	758	775	932	1360	1097	1290	1568	1687	2710	1654	1943	2052	2600	4002
	huff	1101	1101	1136	1137	1177	1491	1668	1815	1895	1896	2830	2348	2349	2920	2920	4167
	motion	344	344	350	392	619	621	548	549	656	716	1211	813	868	868	1110	1817
	dct	1288	1288	1353	1626	1844	2607	2533	2662	3314	3428	5173	3779	3972	4962	5231	7752
	crc	12228	12228	12274	12275	14851	20972	15978	15979	18896	31420 41431	22600 22646 22647	45847 61887
	fir	6852	6852	8053	8594	8715	11595	11780	12802	14903	15564 23063	16890 20491 21031	24154 34520
	mat mul	11142	11142	15015	16039	16359	21597	20968	24841	25961	31084 43178	30795 34668 43884	47400 64733
	average overhead	0.0%	12.2%	18.4% 41.6%	82.8%	72.5%	88.0%	111%	139%	262%	152%	175%	202%	277%	439%

TABLE V :

 V Performance comparison for 8-issue VLIW.

	Benchmarks Original	0	1	2	3	Proposed mechanism 4 5 6	7	8	9	10
	adpcm dec	302	302	309	312	312	370	407	433	443	437	471	558
	adpcm enc	323	323	327	339	342	433	405	478	481	489	572	590
	bcnt	333	333	334	335	335	336	544	525	674	545	664	1130
	fft32x32	400	400	407	424	479	659	655	738	737	875	1385	1388
	huff	951	951	946	946	1039	959	960	1039	1040	1078	1113	1466
	motion	280	280	277	278	284	357	361	358	354	378	415	610
	dct	872	872	881	898	968	1078	1335	1369	1388	1543	2573	2576
	crc	11955	11955	11956 11956 11956	11957	12217	15105	12232	12516	15106	20972
	fir	5709	5709	5890	5889	5892	6072	6854	6853	7750	8111	8708	11814
	mat mul	6533	6533	6534	6534	10619	6537	10630	10950	11717	14323	14833	20356
	average overhead	0%	1.0%	2.2%	4.6%	17.1%	36.7%	42.5%	46.2%	62.9%	94.8%	141%
	Benchmarks Original	0	1	2	3	Instr. duplication (2I) 4 5 6	7	8	9	10
	adpcm dec	302	407	432	443	563	702	710	765	762	773	845	1116
	adpcm enc	323	439	470	478	492	496	756	769	848	897	866	1278
	bcnt	333	547	548	549	549	936	549	933	938	1022	1226	2064
	fft32x32	400	656	705	751	783	839	1015	1230	1391	1509	1778	2789
	huff	951	1038	1039	1040	1075	1076	1610	1605	1764	1744	1830	2800
	motion	280	361	383	383	551	384	558	558	561	655	715	1226
	dct	872	1339	1432	1448	2450	2490	1776	2547	2685	2949	3429	5131
	crc	11955	12230	12231 12231 12515	12517	16219	16216	16231	18908	41344	41382
	fir	5709	6853	7754	7934	7935	8115	10515	13697	13812	14533	15488	23505
	mat mul	6533	10950	11719 11719 11719	19904	19899	20520	22076	22851	28517	40736
	average overhead	38%	46%	48%	71%	96%	108%	143%	155%	174%	225%	378%
	Benchmarks Original	0	1	2	3	Instr. triplication (3I) 4 5 6	7	8	9	10
	adpcm dec	302	612	619	616	836	761	1037	1127	1154	1587	1295	1692
	adpcm enc	323	662	675	725	682	1049	1064	1204	1294	1280	1486	1899
	bcnt	333	764	782	782	783	1389	781	1417	1413	1422	1954	3449
	fft32x32	400	996	1029	1066	1085	1399	1830	1862	2077	4073	2745	4311
	huff	951	1519	1520	1521	1533	1593	1635	2268	2274	2832	2828	4085
	motion	280	496	497	502	597	603	807	790	828	1112	1117	1750
	dct	872	2136	2225	2218	2506	3566	2727	3900	4024	4286	5149	7626
	crc	11955	15962	15963 15978 15980	15980	16979	20019	22906	22914	60038	61808
	fir	5709	10877	11179 11179 11360	16219	17121	19814	19939	20584	24064	35147
	mat mul	6533	16408	16409 16410 16411	29433	30912	30870	33218	33502	43379	61015
	average overhead	104%	108%	110%	124%	200%	205%	259%	277%	362%	404%	626%

4 and 10 faults for the 4-issue and 8-issue VLIW configurations respectively, is the maximum number of faults that the approach can sustain at the same time