Ichrak Toumi
email: toumi@lma.cnrs-mrs.fr

Olivier Derrien
email: derrien@lma.cnrs-mrs.fr

SPARSE DECOMPOSITION OF AUDIO SIGNALS USING A PERCEPTUAL MEASURE OF DISTORTION. APPLICATION TO LOSSY AUDIO CODING

Keywords: Audio coding, Sparse approximation, Iterative thresholding algorithm, Perceptual model

. In future work, a more sophisticated quantization scheme would probably allow our method to challenge standard codecs e.g. AAC.

INTRODUCTION

Actually, state-of-the-art lossy audio coders (e.g. MP3, AAC, OGG) globally share the same structure [START_REF] Spanias | Audio Signal Processing and Coding[END_REF]: The signal is first processed using a time-frequency (TF) transform, the transform coefficients are quantized according to a psycho-acoustic model and finally lossless binary coded. The TF transform is usually adaptive (i.e. it can switch between two pre-defined resolutions), perfectly invertible and introduces no redundancy. The most popular choice is the Modifed Discrete Cosine Transform (MDCT) [START_REF] Princen | Analysis/synthesis filter bank design based on time domain aliasing cancellation[END_REF]. Optimizing the quantization is performed in the frequency domain (i.e. independently in each time slot corresponding to one analysis/synthesis long window or a block of consecutive short windows). This process takes into account a perceptual weighting of frequency components computed by a psychoacoustic model.

However, this structure introduces two main limitations: 1) the invertible and non-redundant (i.e. orthogonal) TF transform takes advantage of an aliasing-cancelation property [START_REF] Princen | Analysis/synthesis filter bank design based on time domain aliasing cancellation[END_REF]. But when some spectral components are removed, because they are considered perceptually less relevant that others, some disturbing artifacts may become audible (called birdies). This occurs usually at low bitrate, where many components have to be removed. 2) Optimizing in the frequency domain does not allow to take into account the time-domain properties of audition (e.g. the temporal masking effect), which implies sub-optimal performance. Some propositions have been made to compensate for these drawbacks (e.g. [START_REF] Derrien | Reduction of artefacts in mpegaac with mdct spectrum regularisation[END_REF][START_REF] Najaf-Zadeh | Use of auditory temporal masking in the mpeg psychoacoustic model 2[END_REF]), but it only resulted in marginal improvements.

The most important breakthrough in the last few years consists in performing a sparse approximation of the signal on an overcomplete dictionary of waveforms instead of using an orthogonal transform [START_REF] Ravelli | Union of MDCT bases for audio coding[END_REF]. This does not guarantee perfect reconstruction, but this property is not mandatory in a lossy codec. This approach was proved to significantly improve the audio quality on some audio files, especially at low bitrate. The proposed method relies on a Matching-Pursuit (MP) derived algorithm, which was modified to reduce artifacts called pre-echo. In its basic version, the MP is optimal with respect to the minimum Mean Square Error (MSE) criterion. Some improvements have been proposed in order to introduce perceptual weights in MP, but this was not selected in [START_REF] Ravelli | Union of MDCT bases for audio coding[END_REF] because it significantly increases the complexity. Recently, a more sophisticated method for sparse audio signal decomposition was proposed [START_REF] Christensen | A perceptually reweighted mixed-norm method for spaser approximation of audio signals[END_REF]. It uses a variational algorithm and includes a perceptual weighting. However, the optimization is still performed only in the frequency domain and its implementation remains uncomplete (no quantization and binary coding were proposed).

In this paper, we describe a new method that is partially inspired by [START_REF] Ravelli | Union of MDCT bases for audio coding[END_REF] and [START_REF] Christensen | A perceptually reweighted mixed-norm method for spaser approximation of audio signals[END_REF]. Our algorithm performs a sparse approximation on a over-complete dictionary, more precisely a union of MDCT bases, as in [START_REF] Ravelli | Union of MDCT bases for audio coding[END_REF]. The optimization uses a variational algorithm, which appears to be more flexible than MP with respect to the distortion constraint. Our main contribution is that this algorithm performs the optimization in the TF plane, using a new TF audition model based on recent studies [START_REF] Necciari | Time-Frequency Masking: Psychoacoustical Measures and Application to the Analysis-Synthesis of Sound Signals[END_REF]. Then, we apply a quantization and binary coding scheme from the literature to evaluate the potential of this method in an audio-coding context. In [START_REF] Ravelli | Union of MDCT bases for audio coding[END_REF], two codecs are described: #1 uses a simple bit-plane algorithm as a quantization and binary coding stage and #2 uses a socalled Perceptual bit-plane algorithm. It was shown that codec #1 usually does not perform as well as AAC, but codec #2 is able to outperform AAC, especially at low bitrate. In this paper, we only implement the simple bit-plane algorithm, because using the perceptual version is not straightforward in our coding scheme. For the same bitrates, we compare the audio quality for our codec and codec #1. In further works, we plan to improve our method by implementing a perceptual quantizer. This paper is organized as follows: In section 2, we present the general method for sparse decomposition. In section 3, we motivate and describe our implementation. And in section 4, we compare our method to codec #1 described in [START_REF] Ravelli | Union of MDCT bases for audio coding[END_REF].

THE METHOD FOR SPARSE DECOMPOSITION

Let's consider a block of audio samples noted as a vector x of size 1 × N , N being the number of samples. x can indifferently represent the whole signal, or a time-segment. We define the coding dictionary as a matrix S of size M × N . Each row can be interpreted an elementary waveform, often called atom. The reconstructed signal is a linear combination of atoms, which can be written as:

x = a S (1)
x is the reconstructed signal and a is a vector of coefficients of size 1 × M , M being the number of atoms in the dictionary. This general scheme applies to state-of-the-art audio codecs: S is determined by the TF transform and a represents the transform coefficients to be quantized and coded. If we assume an invertible dictionary, M = N and the coefficients can be easily computed using:

a = x S -1
This is the case for instance with TF transforms derived from the Fourier transform (DFT, DCT). An over-complete dictionary implies that M > N , which means that S is not invertible. In the general case, there is no a such that x = x. Then, a must be computed using an iterative algorithm, with respect to a pre-defined distortion measure. For instance, the MP algorithm [START_REF] Mallat | Matching pursuits with timefrequency dictionaries[END_REF] finds a which minimizes the mean-square-error (MSE):

D(a) = a S -x 2
For audio signals, MSE is known to be poorly correlated to the auditory perception of distortion. In [START_REF] Christensen | A perceptually reweighted mixed-norm method for spaser approximation of audio signals[END_REF], a more efficient measure is proposed based on the concept of perceptual weights for each frequency component.

Here, we propose a more general formulation of the perceptual weighting technique. This comes from the observation that reasonably good perceptual distortion measures have already been proposed, e.g. the mean Noise-to-Mask Ratio (NMR) [START_REF] Brandenburg | NMR and Masking Flag : Evaluation of quality using perceptual criteria[END_REF]. But these measures are usually associated to an analysis filterbank that follows the characteristics of human perception. Thus, we introduce a second matrix P of size N × K that represents a perceptual transform. We assume that K > N , i.e. this transform is redundant. In the general case, P and S are not directly related, because a good perceptual filterbank (with respect to the accuracy of perceptual modeling) is generally not so good for coding purposes (with respect to the efficiency of coding). The coefficients of this transform corresponding to the signal x can be written as:

px = x P
We define the perceptual distortion measure as a weighted version of MSE computed in the perceptual-transform domain:

Dp(a) = (p x -px) ∆x 2 = (aS -x) P∆x 2
Where ∆x = diag (µ k (x)) is a diagonal matrix of size K × K containing perceptual weights µ k associated to the components of px. These weights depend on x and can be computed using a hearing model of masking. If we note T k (x) the masking threshold values corresponding to the signal x, we choose

µ k (x) = (T k (x)) -1 2 .
Then, the perceptual distortion Dp(a) can be interpreted as a mean NMR. Note that, if the perceptual transform performs a TF analysis with a sufficiently good time precision, both temporal and frequential masking can be modeled.

The general coding problem then consists in finding, for any input signal x, the optimal vector of coefficients a which minimizes the perceptual distortion measure Dp(a). The existence of a solution depends on the matrix K = SP called the mixture matrix. It is quite easy to prove that the maximum value for the rank of K is N . The minimization problem would have a single minimum if K would have a full-rank (i.e. min{M, K}), which is never the case here because N < min{M, K}. In other words, there are many equivalent solutions. We must add a sparsity constraint on a to select the best solution, i.e. we add a regularization term Ψ to the distortion measure and then we minimize the following objective function:

Φ(a) = Dp(a) + Ψ(a)
This class of problem can be solved using an iterative thresholding algorithm [START_REF] Daubechies | An iterative thresholding algorithm for linear inverse problems with a sparsity constraint[END_REF][START_REF] Beck | A fast iterative shrinkagethresholding algorithm for linear inverse problems[END_REF]. However, when K has not a maximum rank (i.e. N ), the convergence is not always satisfactory.

In the literature, Ψ is often assimilated to an α norm on the coefficients, where α is usually equal to 1. Although 1 norm does not directly measure the amount of zero coefficients, one usually assume that the minimum 1 norm corresponds to a sparse solution [START_REF] Daubechies | An iterative thresholding algorithm for linear inverse problems with a sparsity constraint[END_REF]. Then,

Φ(a) = (aS -x) P∆x 2 + λ a α (2)
where the λ parameter allows to set the tradeoff between the distortion constraint and the sparsity constraint.

OUR PROPOSED IMPLEMENTATION

Figure 1 shows a block-diagram of our coder. It is mainly composed of 3 parts: The adaptive analysis over a coding dictionary which is a union of MDCT bases, a psycho-acoustic model and a quantization and coding section. The main point consists in finding a suitable coding dictionary S and a perceptual transform P such that K = SP has a maximum rank.

Union of MDCT bases

First, we assume a coding scheme that splits the audio file in timesegments of reasonably long duration, called macro blocks, in order to be compatible with audio streaming applications (e.g. over the internet). Thus, we assume that the signal x corresponds to a macro-block. This approach requires an overlap between macroblocks, and we wish that it does not introduce redundancy. Thus, a suitable choice for S is a union of MDCTs of different sizes, in the same way as in [START_REF] Ravelli | Union of MDCT bases for audio coding[END_REF]. In this paper, the authors proposed a union of 8 MDCTs: 64, 128, 256, 512, 1024, 2048, 4096 and 8192 bands.

In this study, we restrained our dictionary to a union of 2 MDCTs with medium numbers of bands: 256 and 2048. In this implementation, one macro-bloc is composed of 16 short-windows (for 256band) and 2 long-windows (for 2048-band). Choosing the number of MDCT sizes is tricky: More MDCTs, especially longer ones, would increase the efficiency of the decomposition, but would also degrade the efficiency of the coding stage: More MDCTs means more coefficients in a, and even with the same number of non-zero coefficients, the additional zero coefficients require more coding bits.

The vector of MDCT coefficients a is made of two parts: Figure 2: MDCT windows in the coding dictionary for one macroblock, and corresponding analysis segment.

Perceptual model

For the perceptual transform, we seek for a filterbank that follows the audition, but also with nice mathematical properties. The rank constraint on K requires that P is a full-rank matrix, i.e. its rank is N . A good choice is the ERB-MDCT, which is a near-orthogonal transform that follows the ERB frequency-scale [START_REF] Derrien | A quasi-orthogonal, invertible and perceptually relevant time-frequency transform for audio coding[END_REF]. In this transform, the time-resolution is adapted to the frequency resolution in each frequency band: Low frequencies have a high frequencyresolution and a low time-resolution, whereas high frequencies have a lower-frequency resolution and a higher time-resolution.

We choose an ERB-MDCT with one band per ERB, which corresponds to K = 15126. One can check that we actually have K > N .

The square-value of ERB-MDCT coefficients is interpreted as the temporal and spectral density of energy. To obtain a masking threshold, we convolve the TF energy image with the TF masking kernel described in [START_REF] Necciari | Time-Frequency Masking: Psychoacoustical Measures and Application to the Analysis-Synthesis of Sound Signals[END_REF], and plotted on figure 3. We assume the additivity of masking patterns in the energy scale, which is known to sometimes under-estimate the masking level. Finally, we keep the minimum value of the masking threshold and the absolute threshold of hearing in quiet as described in MPEG #1 psycho-acoustic model [13], and get an estimation of the µ k (x).

Adaptive decomposition

For the optimization in the adaptive analysis step, we used the FISTA algorithm proposed by Beck et al. in [START_REF] Beck | A fast iterative shrinkagethresholding algorithm for linear inverse problems[END_REF]. The algorithm includes a gradient step followed by a shrinkage/soft-thresholding step. The general update rule in FISTA is

ai = T λγ (G(yi))
where i is the iteration index, γ is the gradient stepsize, G(.) is the gradient of the distortion term Dp, yi is a specific linear combination of the previous two iterations {ai-1, ai-2} and T λγ : R M → R M is the shrinkage operator:

T λγ (a)m = sign(am) max (0, |am| -λγ)

A typical condition which guarantees the convergence to a single minimiser a is γ ∈ 0, 1/ K T K∆x .

In the following, we assume for each macro-block:

• α = 1, which amounts to perform a soft thresholding on the coefficients a.

• γ = 1/ K T K∆x .

• For i = 0, a0 = 01×M .

A crucial parameter for the optimization process is λ. We found out that a constant value for λ does not always produce the same sparsity ratio for a. Thus, the suitable value for λ has to be set for each macro-block in order to get a constant sparsity ratio. Our goal was to get approximately the same rate of non-zero coefficients in a as with quantized MDCT coefficients in a MPEG-AAC bitstream. We measured that, for a monophonic AAC operating at 48 kbps, approximately 30 % of the MDCT coefficients are equal to zero. Then, we target a sparsity ratio of 15 %, since our dictionary has a redundancy factor of 2.

Furthermore, we found out that setting the same value of λ for a256 and a2048 was not the optimal choice. Thus, we rewrite equation ( 2) as

Φ(a) = (aS -x) P∆x 2 + λ1 a256 α + λ2 a2048 α (3)
where λ1 and λ2 are respectively the regularization parameters corresponding to a256 and a2048. From informal listening tests, it appears that the best audio quality is obtained with λ1 > λ2, which means more non-zero coefficients in the long MDCT than in the short one.

In the literature, the choice of a suitable value for λ is known to be a difficult problem [START_REF] Hansen | Rank-Deficient and Discrete Ill-Posed Problems[END_REF][START_REF] Tibshirani | Regression Shrinkage and Selection Via the Lasso[END_REF], that would deserve a more specific study by itself. In this work, we use a simplified formula from [START_REF] Daubechies | An iterative thresholding algorithm for linear inverse problems with a sparsity constraint[END_REF] which is valid in the following conditions: (i) K is the identity operator and (ii) α = 1:

T λ (x) =    x + λ 2 if x ≤ -λ 2 , 0 if |x| < λ 2 , x -λ 2 if x ≥ λ 2 . ( 4
)
Then, assuming that ∀k ∈ {1

• • • K} , |px(k)∆x(k)| ≥ λ 2 , we get: λ = 2Ys K(1 -sp) + 1 ;
where sp is the amount of sparsity, Ys are the values of |px∆x| sorted in increasing order, and . denotes the rounding operator towards -∞. Although K is not the identity operator in our case, this formula still gives a good approximation of the optimal λ in each macro-block. Thus, controlling sparsity only requires to a adjust the parameter sp.

Coefficients quantization and coding

When the coefficients a are computed, we carry out the quantification and binary coding preceded by an interleaving step which aims to group together the coefficients that are close in the timefrequency plane, and to make the coding algorithm more efficient. We chose the same binary coder as in [START_REF] Ravelli | Union of MDCT bases for audio coding[END_REF]: a bit-plane coder, which is a special case of adaptive Rice-Golomb codes. It also performs an implicit quantization: the number of coding bits corresponds to the number of bit-planes used for coding. This method is particularly efficient when the sequence to be encoded exhibits long ranges of zeros, and thus is well suited for encoding sparse vectors. However, our interleaving scheme, as illustrated on figure 4, is different from the one described in [START_REF] Ravelli | Union of MDCT bases for audio coding[END_REF]. The idea is to group the low frequencies at the beginning and the high frequencies at the end. We found out that our method reduces the bitrate, because it favors long ranges of zeros at the end.

RESULTS AND DISCUSSION

To test our method, we used the same audio material as in [START_REF] Ravelli | Union of MDCT bases for audio coding[END_REF]. This collection of 4 musical pieces, sampled at 44.1 KHz, is described on table 1. Note that all the audio signals mentioned in this paper (unprocessed, re-synthesized and coded/decoded signals) can be downloaded from the companion webpage of this paper: http://potion.cnrs-mrs.fr/dafx2015.html. Table 1: Test signals used for evaluation.

Performance of the adaptive decomposition algorithm

In a first step, we analyze the results of the adaptive decomposition algorithm on the harp signal. For evaluation purpose, we re-synthesize the audio signal without quantization by applying equation [START_REF] Ravelli | Union of MDCT bases for audio coding[END_REF]. As explained in section 3.3, the sparsity ratio is about 10% on 256-band MDCT coefficients and 20% on 2048band MDCT coefficients. We can see on the spectrograms plotted on figure 5 that the original and re-synthesized signals are very similar. However, the reconstructed signal has less energy on the regions of the TF plane that are between the partials, which should correspond to masked regions. One can also notice that partials sometimes cross vertical structures corresponding to attacks, which can be associated to pre-echo. This is most probably due to the fact that our short MDCT has 256 bands, whereas in [START_REF] Ravelli | Union of MDCT bases for audio coding[END_REF], the shortest MDCT has only 64 bands.

We also plot the TF maps of the coefficients a256 and a2048 on figure 6. In other words, we plot the magnitude of MDCT coefficients in the TF plane for both MDCT bases separately. First, one can see that the long MDCT exhibits a good frequency resolution, but a poor time-resolution. The short MDCT has opposite properties. One can also see that partials are essentially represented by long MDCT coefficients, and attacks by short MDCT coefficients. This proves that our optimization algorithm dispatches the energy between MDCT bases in an accurate way. However, this process is not perfect: some partials are represented in the short MDCT.

Finally, we plot on figure 7 the spectrograms of the re-synthesized signal after quantization and coding at 48 kbps and 24 kbps, for the same original test signal. At 48 kbps, one can see that the spectrogram of the re-synthesized signal is highly similar to the spectrogram before quantization and coding. However, at low

Evaluation of audio quality

In order to evaluate the impact of our method on the audio quality and to compare our codec (denoted in section codec A) to the codec #1 proposed in [START_REF] Ravelli | Union of MDCT bases for audio coding[END_REF], we organized a MUSHRA listening test [START_REF]Method for the subjective assessment of intermediate quality levels of coding systems[END_REF] on the four audio signals described in table 1, with scores ranging from 0 (very bad quality) to 100 (excellent quality). Six versions of the test signals were evaluated by a total number of 16 listeners:

• A hidden reference,

• A 4 kHz low-pass anchor,

• Two coded versions with codec A at 24 and 48 kbps,

• Two coded versions with codec #1 at 24 and 48 kbps.

The listeners were post-screened according to the scoring of the reference: If the score attributed to the reference is lower than 80, the listener was judged unreliable and discarded. Finally, only 10 subjects were kept.

The results of MUSHRA listening tests are plotted on figure 8 and9, respectively at 24 and 48 kbps (mean values and 95% confidence intervals). For both bitrates, one can see that results highly depend on the test signal. Signals with naturally sparse spectrums (harp, bagp) are associated to high audio quality at both bitrates. This can be interpreted as follows: The amount of perceptually relevant information in these signals is relatively low, and a good reconstruction in achievable at low bitrate. For bagp, the results obtained with both codecs are similar at 24 and 48 kbps but for harp, codec #1 is better rated. According to the listeners' feedback, this was mainly due to the fact that attacks are not as sharp with our codec as with codec #1, which is probably related to the length of the short MDCT. Signals with more dense spectrums (orch, popm) are associated to lower audio quality, especially a 24 kbps, probably because the amount of perceptually relevant information in these signals is higher. At 24 kbps, both codecs perform similarly on these two signals, but at 48 kbps, our codec is slightly better on popm, and much better on orch. Note that, at 24 kbps, both codecs were rated lower than the anchor, because a reduced band-pass was judged less disturbing than a high level of coding artifacts.

Finally, on average over the 4 test signals, our codec was rated slightly worse at 24 kbps, only because it does not perform as well on harp, and rated significantly better at 48 kbps, because it performs better on orch and popm.

CONCLUSION

In this paper, we described a method that performs an adaptive decomposition of audio signal on redundant coding dictionary (a union of 2 MDCT bases) using a variational algorithm, i.e. the optimization of a cost function taking into account both a perceptual distortion measure derived form a hearing model and a sparsity constraint. We applied a simple quantization and coding scheme from the literature (simple bit-plane coder) and compared the final audio quality to one achieved by codec #1 in [START_REF] Ravelli | Union of MDCT bases for audio coding[END_REF], which uses the same quantization and coding scheme. We show that at low bitrate (24 kbps), our codec performs worse on signal with many sharp attacks, and performs similarly on other signals. This can be explained by the fact that our short MDCT has only 256 bands. At medium bitrate (48 kbps), our codec performs better on audio signals with a dense spectrum. This point is particularly interesting since in [START_REF] Ravelli | Union of MDCT bases for audio coding[END_REF], codec #2 was able to outperform AAC but not on signals with dense spectrum (orch and popm), despite the use of a complex perceptual quantization scheme. Thus, these results are promising: As it is, our codec can not outperform codec #2 or AAC, but it may seriously challenge both if we design a suitable perceptual quantization scheme. This would also require to extend the coding dictionary by adding a shorter MDCT (128 or 64 bands).

Figure 1 :

 1 Figure 1: Diagram of the proposed coder.

a

  = [a256 a2048]where a256 stands for the coefficients of the short MDCT and a2048 for the long MDCT. The time-support corresponding to one macro-block is determined by the long MDCT, i.e. N mb = 3 × 2048 = 6144 samples. But practically, the optimization must be performed on an analysis-segment longer than a macro-block, in order to avoid sharp variations of coding parameters between successive macro-blocks. Thus, N corresponds to the length of one analysis segment. Here, we chose N = 7 × 2048 = 14336 samples. This is illustrated on figure2. Both a2048 and a256 are of size 1 × M as , with M as = 6 × 2048 = 12288, and a is of size 1 × 2M as , which means M = 2M as = 24576. But practically, only 8192 coefficients (in the middle) are actually coded. The first and last 8192 coefficients are discarded. This part of signal will be analyzed and coded respectively in the previous and the next macro-blocks.

Figure 3 :

 3 Figure 3: Time-frequency masking kernel.

Figure 4 :

 4 Figure 4: Interleaving scheme for MDCT coefficients.

Figure 5 :

 5 Figure 5: Spectrograms of the harp signal

Figure 6 :

 6 Figure 6: TF maps of coefficients for the harp signal. Up: short MDCT. Down: long MDCT.

Figure 7 :

 7 Figure 7: Spectrograms of the harp signal after quantization and coding.

Figure 8 :

 8 Figure 8: Results of MUSHRA listening test at 24 kbps.

Proc. of the 18 th Int. Conference on Digital Audio Effects (DAFx-15), Trondheim, Norway, Nov 30 -Dec 3, 2015

* This work was supported by the joint French ANR and Austrian FWF project "POTION", refs. ANR-13-IS03-0004-01 and FWF-I-1362-N30.


Mathieu Delangle
email: mathieu.delangle@irccyn.ec-nantes.fr

Jean-François François Petiot
email: jean-françois.petiot@irccyn.ec-nantes.fremilie.poirson@irccyn.ec-nantes.fr

Emilie Poirson

Delangle Mathieu

Irccyn

Assessing the differences between numerical methods and real experiments for the evaluation of reach envelopes of the human body

Keywords: Reach assessment, User-centered design, Experiments, Numerical evaluation, Comparative study

Assessing the differences between numerical methods and real experiments for the evaluation of reach envelopes of the human body

I. INTRODUCTION

Modeling normal human reach is widely used by the designer to design and assess the accessibility of environments. Many tools and practices are based on anthropometrics data to perform these ergonomics evaluations [START_REF] Roebuck | Anthropometric methods: Designing to fit the human body[END_REF]. Appropriate anthropometric data regarding body size from the established data bases are used to analyze and design the intended product or environment. Existing database (e.g. ANSUR [START_REF] Gordon | anthropometric survey of U.S. army personnel: Methods and summary statistics[END_REF], NHANES [START_REF] Nchs | Analytic and reporting guidelines: The third national health and nutrition examination survey, NHANES III[END_REF]) are generally chosen as the reference population and are used directly to assess the environment, instead of recruiting sample-users to test the product. However, the reach behavior on an individual depends on many factors. Anthropometry, age, gender, joint mobility and muscle strength are a few such factors related to the individual being modeled. This human variability might cause difficulty to meet the requirements of a conception, especially for a certain percentage of the user population. Although numerical methods are faster and less expensive than the involvement of sample-user to test prototypes, the only use of these anthropometric data [START_REF] Moroney | Empirical reduction in potential user population as the result of imposed multivariate anthropometric limits[END_REF], often old, are not always representative of the target users population. In fact, they often consist of specifics surveys (military...), and typically provide only very limited information concerning children and people who are older and disabled [START_REF] Marshall | Using Hadrian for eliciting virtual user feedback in design for all[END_REF]. Moreover, these data are generally used in univariate case study (to determine the appropriate allocation of adjustability to achieve a desired accommodation level), where most problems are multidimensional. That is why, although design methods based on external body dimensions don't need experimental tests or the building of prototypes, methods based on this principle still pose questions about their ease of use and their reliability compared to reality. This paper describes to which extent an evaluation of the reach only based on the structural data of the human body may differ from those obtained with an experimental task. Thus, an experimental reach assessment and a numerical evaluation were performed and compared to highlight differences.

II. COMPARATIVE STUDY

The present study proposes to compare two ways of accessibility assessment of the human body. First during an experiment, a sample-user is asked to perform an accessibility task. Second using numerical data and a kinematic model of the body, the structural data of the participant are directly used to numerically assess the accessibility (Figure 1). The task proposed is a 2 dimensional reaching task of the hand of the participant in the frontal plane (Figure 2). Experimental test was divided into two sub-tests; a constrained test (feet fixed to the floor and body fixed relative to the vertical axis of the center of the plate) and an unconstrained test (feet fixed to the floor and rest of body free to move). The unconstrained 2/6

La Plagne (73) -31 Mars au 02 avril 2015 situation represents the functional reach, that is to say the maximal distance one can reach with the hand while maintaining a fixed base of support in the standing position [START_REF] Duncan | Functional reach: A new clinical measure of balance[END_REF], [START_REF] Duncan | Functional reach: Predictive validity in a sample of elderly male veterans[END_REF].

III. EXPERIMENT

A. Sampling

The experiments were conducted with 40 adult volunteers, all French students or teachers of the Ecole Centrale de Nantes. Twenty five males and fifteen females were sampled in the study, covering a wide spectrum of physical characteristics, from 1482 mm for the smallest stature, to 1930 mm for the highest. The average of stature of the subjects was 1735 mm (with a standard deviation of 95.4). The summary statistics of the anthropometry of the subjects are provided in Table 1. No one of these volunteers reported motor disabilities or particular physiological limitations. The sample is considered as representative of the general population. It was a deliberate decision to not skew the data by "excluding" persons in the panel (e.g. old or disabled), in order to not bias the comparison. Indeed, the accessibility study of persons with specific physical limitations will depend of more parameters that would make the comparison more difficult.

B. The tests Principle:

Marking of the boundaries of the hand reach in standing posture. The boundary of reach volumes based on the reach capabilities of both arms were measured based on a technique used in [START_REF] Molenbroek | Reach envelopes of older adults[END_REF], [START_REF] Marshall | A taskbased design for all support tool[END_REF] and [START_REF] Gyi | Representing older and disabled people in virtual user trials: data collection methods[END_REF]. Briefly, a white board was fixed vertically, and individuals were asked to produce an 'arc' averaging 90° (from the shoulder height to the top of the head). In order to check if the subjects does not lift their feet to increase their vertical reach, an electronic sensor position was positioned under the heels indicating if the feet are off the ground or not. When the visual signal was triggered, the test was stopped and the subject repositioned. For each individual, the reached envelope was drawn on a white board. Arcs were drawn with the body, providing the volume described by the reach capability for both arms. Envelopes were marked to represent the task as "finger touch" function (one finger touches an object without holding it) in order to avoid much as possible grasp effects (reducing the reach envelope). Because most anthropometric data presented in databases represent nude body measurements and to permit reliable comparison, experiments were performed with light clothing (nude dimension and light clothing being regarded as synonymous for practical purposes). produced were captured on photo for analysis; pictures were treated with a MATLAB® program to accurately determine the area of the envelope (Figure 3). Each subject performed this test firstly constrained (foot and body fixed) and secondly unconstrained (only foot fixed). A similar protocol where participants are in interaction with physical points may be found in [START_REF] Delangle | Helping the designer to choose his software: Comparative analysis of human modeling tools[END_REF].

IV. STATIC NUMERICAL EVALUATION

A. Creation of the anthropometric database When using methods for ergonomic evaluations, anthropometric modeling can be either directly observed from anthropometric characteristics of the current users, or statistically derived from characteristics for the intended target population. In order to predict the accessibility by limiting the statistical biases in the comparison, accurate anthropometric data from the participants were needed. Thus, the presented evaluations were all based on the anthropometric characteristics of the subjects who performed the experimental tests. External anthropometry being the type most frequently available and collected, it was decided to collect some "direct"' measurements of external link-length anthropometry. Data in Table 1 were recorded to predict the upper body accessibility for each participant, and were collected in a laboratory environment from the 40 individuals. It was expected that this number would provide a manageable database for the development and validation of the comparative study.

B. The numerical model

The aim is to evaluate the reach characteristics from the recorded external anthropometrics that might be found in anthropometric database (Table 1). This methodology is based on the design limits approach, which is a common method used in design problems, where data about human physical characteristics are directly applied to solve design problem. The maximal reach is directly correlated to the greatest distance between the shoulder acromion and the fingertips, corresponding to an outstretched arm situation. This is kinematically modeled as a simple link (arm length) with a unique revolute joint (shoulder). So, the maximal reach envelope is defined by an arc circle, with a radius equal to the arm length of the operator and the shoulder as point of rotation (Figure 4). All points within this envelope (shaded area) were considered as reachable by the subject. The total reach area was defined using Equation 1. H represents the tip of the hand, S the shoulder location and O the center of the shoulder width, thus OS the half shoulder width and SH the arm length. C1 represents the arc circle area between H and S, C2 the arc circle area between S and O, and C3 the area of the SOP triangle.

Knowing the anthropometric characteristics of each participant, a program was implemented (using MATLAB® R2012b) allowing to automatically determinate the corresponding reach. The results are presented in section V.

Eq.( 1) :

( ) ∫ √ ∫ √ .OP

V. RESULTS

Results from experiments and evaluations were collected for each of the forty participants. The aim is to compare the shape of the reach envelopes obtained from the experimentation with those obtained from numerical assessments. Obviously, the anthropometric characteristics of the participants being different, the shapes of the envelopes are all different too. So as to make an objective comparison of the general behavior for the forty profiles, a normalization procedure was performed to compare the reaches between all the participants. First, experimental reach envelopes were recorded and drawn for each subject. Second, each envelope was rescaled in a common domain of comparison (1 unit representing the arm length of each participant). Thus, for all of them, the envelope defined by the numerical model (section IV.B) is modeled by a common arc circle, where 0 corresponds to the shoulder position and 1 to the theoretical extremity of the hand for a static position. This normalized theoretical reached envelope was used as a basis for comparisons to highlight the differences between the real reaches and the numerical model. Finally, the forty normalized reaches obtained from the experiments were aggregated on the same graph and the extreme boundaries (corresponding to the extreme position of fingertips) were extracted. The average extreme reach obtained for the sample-users is depicted Figure 5.

The sample of population being considered as representative of a normal distribution, statistical evaluations were made to represent, depending on the position of the arm, the variability of reach. This variability is the differences between the theoretical curve based on the angle of the arm, and the calculation of an overall confidence reach interval (95% of the 4/6 La Plagne (73) -31 Mars au 02 avril 2015 population). Given the sample size "n" (>30), it is estimated that the mean x of the sample follows a normal distribution. Thus, the confidence interval of 95% could be determinate using Equation 2. Where n is the size of the sample, σ the standard deviation and x the mean. The results of the differences between the experimental continuous reaches and those numerically found are presented Table 2 for the constrained test and Table 3 for the unconstrained. The associated graphic representation is shown Figure 5.

Eq. ( 2): The 95% interval indicates that there is 95% chance that the interval contains a new observation. It can be seen that, for the constrained and unconstrained experiments, the means curves, not perfectly follow the theoretical envelope numerically defined. Looking to the shape of the constrained experiment, it can be seen that the vertical reach is perfectly coincident with the theoretical results. The mean (with the corresponding interval confidence) is 1±0.02. So, for this situation, the experiment is coherent with the numerical model based on the structural body dimensions. For the lateral reach, the results show that the experimental reach is greater than those of the theoretical model, with an average of 1.086±0.024.

Results obtained for the unconstrained show that the reach is overall more important, compared to the numerical model and to the constrained experiment. The normalized mean is 1.076±0.018 for the vertical location, and 1.363±0.034 for the lateral reach.

In the constrained case, the extended lateral reach is certainly due to small displacements during the test. Indeed, although this experiment was implemented in order to limit as much as possible movements of the body, displacements of the upper body could appear during the tests. For the unconstrained task, this extended lateral reach represents the realistic functional reach of the subject. The body being not constrained (excepted the feet), the maximal reach strongly depends on the lateral reach capacity, represented by the medio-lateral balance stability and the pelvis rotation ability of the subject, as depicted Figure 6. Thus, the maximal lateral reach is greatly higher than expected. Considering the average arm length, this average maximal lateral reach corresponds to a deviation of 26cm, which is consistent with results that can be found in the literature [START_REF] Brauer | Lateral reach: a clinical measure of medio-lateral postural stability[END_REF]. The vertical envelope obtained from the constrained experiment is perfectly coherent with the numerical model, which is not the case for the unconstrained test. This involve that the movement executed with a free body allows the participants to increase their vertical reach. The main reason of this increase is probably due to the complex model of displacement of the shoulder, and particularly of its flexion, allowing the subjects to increase their reach upwardly, beyond their external structural dimensions. With the foot fixed on the ground, the strategy to enhance the vertical reach is to increase the flexion of the shoulder and to change the alignment of the shoulder and the pelvis The discrepancy between the constrained and the unconstrained reach is the highest for the lateral location, and decreases continuously up to the vertical. This represents the differences of strategy used throughout of the envelope to improve the maximum reach. The highest variability of the reach appears for the unconstrained test at 0°, a standard deviation of σ=0.104. Indeed, the balance capacity of each subject might be very different, involving an important variation of the reach from one participant to another. This variations would be greater for a sample of people with very different physical characteristics (disabled persons, children…), which is not the case in this study. So, the biomechanics characteristics of the body directly impact the maximum reach capacity.

VI. CONCLUSION

This study compared results of a reach assessment obtained from experiments and from numerical evaluations. In spite of the simple nature of the presented task (two-dimensional), the results show important differences between the two evaluations. The reach capacities of the participants were not only correlated to their anthropometric characteristics. The biomechanics of the body and the physical abilities of the human implies that the maximum reach change during the reach. The reach capacity does not represent a perfect arc circle but a curve depending of the hand position in the space. Anthropometric design problems associated with human physical characteristics depend on many other factors, related to the task, the body position effects and the human behavior in interaction. Several task considerations should be taken into account in order to construct a reach envelope, as the nature and requirements of the task to be performed, the body position while reaching, the whole body movement capabilities and restraints. Moreover, numerous human variability factors might affect the reach results, as the age, gender, body build, fatigue, disease, clothing or environment. These parameters have to be taken into account in the reach evaluation, especially for multivariate design problem. Perform accessibility evaluations only considering the structural limits data of the problem can rapidly leads to misfit design solutions.

VII. PERSPECTIVES

In continuation of this study, a methodology using threedimensional human simulation is being investigated in order to incorporate real human reach capacity in virtual simulations, without the participation of sample-users. Based on CAD environment and virtual reality tools, reach assessment of multi-dimensional design problem will be performed in an intuitive and interactive manner. Thus, the external factors influencing the reach and the differences of physical capacities may be taken into account to carry out robust reach assessments, especially in the universal design field.

VIII. REFERNCES

Figure 1 .

 1 Figure 1.SYNOPSIS OF THE STUDY

Figure 2 .

 2 Figure 2. PARTICIPANT DRAWING HIS REACH ENVELOPE DURING THE CONSTRAINED TEST.

Figure 3 .

 3 Figure 3. IMAGE PROCESSING USING MATLAB®. BRIGHT SPACE REPRESENTS THE REACH ENVELOPE AREA DRAWN BY THE PARTICIPANT. THE DARK AREA (RECTANGLE) REPRESENTS THE PAPER SHEET USED AS REPOSITORY TO CALCULATE AREAS FROM PIXELS TO SQUARE METERS.

Figure 2

 2 Figure2shows a subject performing the reach envelope task using the traditional reach envelope board method. The arcs

Figure 4 .

 4 Figure 4. STATIC MODEL OF REACH ENVELOPE (SHADED AREA) OF THE RIGHT ARM CALCULATED FROM EXTERNAL ANTHROPOMETRICS DIMENSIONS OF THE SUBJECT.

  RESULTS OF THE CONSTRAINED TASK: POSITIONS OF THE HAND FOR SEVERAL ARM ANGLES, AFTER AGGREGATION OF THE 40 CONTINUOUS REACH ENVELOPES; WITH X THE MEAN, Σ THE STANDARD DEVIATION, I + AND I -RESPECTIVELY THE SUPERIOR AND INFERIOR LIMITS OF THE 95% CONFIDENCE INTERVAL. Table 3. RESULTS OF THE UNCONSTRAINED TASK: RESULTS OF POSITIONS OF THE HAND FOR SEVERAL ARM ANGLES, AFTER AGGREGATION OF THE 40 CONTINUOUS REACH ENVELOPES; WITH X THE MEAN, Σ THE STANDARD DEVIATION, I + AND I -RESPECTIVELY THE SUPERIOR AND INFERIOR LIMITS OF THE 95% CONFIDENCE INTERVAL.

FIGURE 5 .

 5 FIGURE 5. NORMALIZED MEAN REACH LINES OBTAINED FROM THE EXPERIMENT, WITH THE 95 PERCENT CONFIDENCE INTERVAL. THE SHADED AREA REPRESENTS THE NORMALIZED THEORETICAL REACHED AREA OBTAINED FROM THE STATIC NUMERICAL EVALUATION.

FIGURE 6 .

 6 FIGURE 6. PARTICIPANT EXTENDED HIS LATERAL REACH USING THE MEDIO-LATERAL BALANCE STABILITY AND THE PELVIS ROTATION, DURING THE UNCONSTRAINED TEST.

Table 1

 1

	.	ANTHROPOMETRIC CHARACTERISTICS OF THE SUBJECTS (IN
	MM) PARTICIPATING TO THE EXPERIMENT, WITH THE AVERAGE AND THE
	STANDARD DEVIATION (EXCERPT).


C Brun

X Buet

B Bresson

M Ciccotti

A Ghomari

P Lecomte

J P Roger

M N Petrovich

F Poletti

D J Richardson

D Vandembroucq

G Tessier

Picometer-scale surface roughness measurements inside hollow glass fibres

Keywords:

A differential profilometry technique is adapted to the problem of measuring the roughness of hollow glass fibres by use of immersion objectives and index-matching liquid. The technique can achieve picometer level sensitivity. Cross validation with AFM measurements is obtained through use of vitreous silica step calibration samples. Measurements on the inner surfaces of fiber-sized glass capillaries drawn from high purity suprasil F300 tubes show a sub-nanometer roughness, and the roughness power spectrum measured in the range [5 • 10 -3 µm -1 -10 -1 µm -1 ] is consistent with the description of the glass surface as a superposition of frozen capillary waves. The surface roughness spectrum of two capillary tubes of differing compositions can be quantitatively distinguished.

Introduction

Hollow core-photonic bandgap fibres (HC-PBGFs), in which a periodic cladding made of air holes allows light to be guided in a central air core, are of great interest as a potential transmission medium in next-generation optical communication systems. Their nonlinearity and latency are much lower compared to classical solid-core [START_REF] Poletti | Hollow core photonic bandgap fibres; technology and applications[END_REF][START_REF] Sleiffer | High capacity mode-division multiplexed optical transmission in a novel 37-cell hollow-core photonic bandgap fiber[END_REF], and they offer reduced Rayleigh scattering and absorption by impurities or phonons. However, losses in these fibres are still one order of magnitude higher than in typical single mode fibres. Scattering at the numerous air/glass interfaces present in these structures, particularly those at the hollow core boundary, has been identified as a major loss mechanism [START_REF] Roberts | Ultimate low loss of hollow-core photonic crystal fibres[END_REF]. In particular, short spatial frequencies (i.e. below 1 µm -1 ) seem to have the highest impact on losses [START_REF] Fokoua | Analysis of light scattering from surface roughness in hollow-core photonic bandgap fibers[END_REF]. Longitudinally, roughness should therefore be characterized over distances of up to several hundreds of micrometers.

The optimization of multi-mode propagation offers a way to control and manage the interaction of light with interfaces [START_REF] Wheeler | Wide-bandwidth, low-loss, 19-cell hollow core photonic band gap fiber and its potential for low latency data transmission[END_REF], and reduce losses [START_REF] Fokoua | Analysis of light scattering from surface roughness in hollow-core photonic bandgap fibers[END_REF][START_REF] Poletti | Towards high-capacity fibre-optic communications at the speed of light in vacuum[END_REF]. Enlarging the air core or moving to lower wavelengths are known ways to reduce roughness scattering losses. Nevertheless, in order for PBGFs to compete with classical solid-core fibers in long haul communications, the ability to measure and reduce glass surface roughness is essential.

The amorphous structure of a glass reflects the liquid state. In the same way, the glass surface inherits its roughness from the fluctuations of the liquid surface. During the fibre drawing process, the vitreous liquid goes through a glass transition in which the surface capillary waves (SCWs) that are thermally excited in the liquid state are suddenly frozen [START_REF] Roberts | Ultimate low loss of hollow-core photonic crystal fibres[END_REF][START_REF] Ja | Intrinsic roughness of glass surfaces[END_REF][START_REF] Sarlat | Frozen capillary waves on glass surfaces: an AFM study[END_REF]. According to this scenario, the roughness statistics of a glass surface are therefore the same as those of a snapshot of the liquid surface taken at the glass transition temperature T g . This fundamental thermodynamical process imposes the ultimate limit to roughness in glasses, and in HC-PBGFs in particular. The resulting roughness is extremly low. It is controlled by the competition between the thermal noise at the glass transition k B T G (where k B is the Boltzmann constant) and the surface tension γ. Typical values of the RMS roughness amount to σ ≈ k B T G /γ ≈ 0.4nm for amorphous silica. This global estimate actually results from the superposition of surface capillary waves modes with spatial wavelengths that may extend from the nanometric to the millimetric range [START_REF] Sarlat | Frozen capillary waves on glass surfaces: an AFM study[END_REF].

Pushed close to their instrumental limits, AFM microscopes have proven their ability to measure sub-nanometer roughness on flat glass surfaces, and to obtain roughness density spectra which are well described by theory [START_REF] Sarlat | Frozen capillary waves on glass surfaces: an AFM study[END_REF]. In HC-PBGF, however, such measurements are much where f is the focal length of the microscope, γ o the magnification of the objective, and G oc the angular magnification of the ocular lens. After reflection on the sample, these beams are recombined in the Wollaston prism, exactly along their previous path if the Wollaston prism is in the focal plane of the microscope (see Fig. 1). The phase between the two components of the polarisation is then shifted by ϕ = 2(ϕ w + ϕ s ), where ϕ w = 2π λ αT w is the phase shift induced by the Wollaston prism, T w being the distance between the axis of the beam and the center of the Wollaston, and ϕ s = 2π λ nD(x) is the phase difference introduced in the propagation medium (index n) if the beams hit the sample surface at heights differing by D(x). After reflection on the beamsplitter, 1/4 of the initial laser power is sent to the homemade photoelastic polarization modulator, in which a birefringence modulated at a frequency F = ω/2π = 50kHz is created by a sinusoidal strain applied to a silica glass rod, inducing a phase modulation Φ(t) = Φ M sin ωt. Both orthogonal polarisations are then projected on the axis of an analyzer forming an angle of 45 • in order to allow interference. Finally, a pinhole, conjugated with the impact spot of both beams on the sample, filters out stray light and parasitic reflections before measurement with a photodiode (Thorlabs PD100). The resulting signal is processed by a lock-in amplifier.

In our analysis, the surface of the glass samples are assumed to be optically homogeneous, i.e. exempt of local refractive index variations, and the reflection coefficients, r, are assumed identical for both beams. The intensity on the detector can thus be written as:

I(t) = rE 2 0 [1 + cos ϕ cos Φ(t) -sin ϕ sin Φ(t)] (1) ≈ rE 2 0 [1 + J 0 (Φ M ) cos ϕ -2J 1 (Φ M ) sin ϕ sin ωt + 2J 2 (Φ M ) cos ϕ cos 2ωt + ...]
where J n are Bessel functions of the first kind. The third term, at the angular frequency ω, is isolated by the lock-in to obtain a signal proportional to sin ϕ thus to ϕ s and D(x). The value of the proportionality factor, which is necessary in order to quantitatively measure D(x), is easily obtained by translating the Wollaston prism to change ϕ w and measuring the amplitude signal for positions corresponding to sin ϕ = sin 2(ϕ w + ϕ s ) equal to 1 or -1.

If we assume no overlap between the two spots, this measurement therefore delivers the height difference D(x) between the two spots; the case of overlapping Gaussian laser spots is slightly more complex and will be described later.

In order to be able to measure inside hollow glass fibres, we use an oil immersion objective (Olympus X20, N.A. = 0.85). A section of the fibre is fixed on a glass slide and immersed in a silica index matching oil. The outer interface of the fibre is then optically invisible, and approximately 4% of the incident light is reflected off the first glass/air interface, i.e. the upper inner surface of the glass fibre.

An important question is whether a portion of the transmitted beam can be collected after multiple reflections inside the hollow cylinder, as this might affect the measurement. Since we use a pinhole to filter out non-confocal light, light can only be collected if it emerges from the illuminated point, which has a radius r = 1.22λ /(2N.A.) = 500nm. This is only possible for beams lying in a plane perpendicular to the axis of the fibre, so the problem is actually that of a 2D reflection at a circular surface. If a beam enters the circle at an incidence angle i, its reflections inside will all form the same angle i. After k reflections, the beam impinges on the inner surface of the cylinder at an azimuthal position given by the angle k.(π -2.i), with relative intensities given by [R s (i)] k and [R p (i)] k , where R s (i) and R p (i) are the Fresnel coefficients for the s and p polarisations, respectively. These beams can only be collected by the confocal microscope if they emerge, after k reflections, at the illuminated point (radius r). We have verified numerically that, for a fibre with an internal diameter of 200µm, the collected intensity for both polarizations is lower than the precision of our calculations, and is therefore negligible. This holds true down to fibre radiuses of the order of the wavelength, i.e. a few hundreds of nanometers, at which this geometrical approach is no longer valid.

The sample is then scanned relative to the beams using a motorized stage (PI M-126-2S1), to obtain a line scan D(x). The differential interferometer is mostly insensitive to unwanted vertical movements, but pitch or yaw angular tilts are crucial, since a sample rotation of e.g. ρ = 10 -6 rad. results in an apparent differential displacement ρd = 3pm. However, if the scanning direction is kept constant, most good quality stages can provide this kind of precision.

An important characteristic of this type of interferometer is that, since the angle α is small (5 ), both beams follow very similar paths and are affected in similar ways by air-flow and vibrations. This makes the instrument relatively immune to these sources of noise.

The ultimate noise limit, which is driven by photon and detector noise, is of the order of 10 -14 m for a lock-in integration time of τ = 100ms, but it is in practice limited to higher values by vibrations and thermal drifts. In order to estimate this noise experimentally, we have carried out repeated measurements on a fixed region for different types of samples, measuring the fluctuations of the signal. Figure 2 depicts the results obtained. The noise level is plotted as a function of one over the square root of the time constant of the lock-in amplifier. As in most lock-in detection systems, drifts also appear in the amplitude signal. In all the experiments presented in this report, we adopted the conservative assumption that these drifts were linear, substracting a linear fit over the raw measurements.

On a perfectly reflecting surface such as a mirror, the noise level reaches the practical noise limit, which is around 4pm for τ = 100ms. Since glasses only reflect about 4% of the incoming light, the available number of photons is notably decreased and the signal to noise ratio is therefore decreased. However by increasing the power of the laser diode until we obtain a reflected signal which has the same intensity as on a mirror, we find the same noise levels. However, immersion measurements on capillary fibres lead to higher noise levels which might be caused by movements in the oil, and by the very low depth of field of immersion objectives.

Within the time constant range used to measure differential profiles (a few ms), the noise levels are around six to ten picometers. This allows precise measurements of the glass roughness. Moreover, we note that, even though the lateral resolution of the differential profilometer is as expected lower than that of an atomic force microscope, the vertical noise level of immersion measurements is at least three times lower than that of an AFM that typically amounts to 30 pm for the low noise commercial instrument used in the present study [14].

Response function

In the following sections we present results obtained on two different kinds of samples. Test silica samples that exhibit a deterministic step height profile of nanometric amplitude will first be examined. Hollow silica tubes characterized by a sub-nanometric random roughness will then be studied. In both cases, results will be quantitatively compared either with independent experimental results (AFM measurements on the silica steps) or with theoretical expectations (frozen capillary waves for the inner roughness of capillary tubes). While the quantitative analysis of a deterministic height profile naturally resorts to the real space, the characterization of a random roughness will be more practical in the Fourier space.

Before coming to the experimental results we give here more details on the experimental response function and its effects on the measurements. The two cases above discussed are presented and both the real space and the Fourier space representations of the response functions are discussed.

As presented in the above section, differential profilometry gives access to a local phase difference ϕ s (x) = (2πn/λ )D(x) on the sample surface due to the difference of optical path between the two reflected beams. If the distance d separating the two beams is significantly larger than the beam waist w, d >> w we expect that the measured signal directly gives the height difference between the two points of refllection:

D(x) ≈ h (x + d/2) -h (x + d/2) ≈ d.h (x) (2)
In practice however, measurements are performed with a separation between the two beams of the same order of magnitude as the beam waist, d/w ≈ 1 -3 and the two beams can partially overlap. For small values of d, the response of the profilometer is thus smeared out, essentially due to the diffraction limit which sets the diameter of the spots.

Step geometry-Real space response function

Let us first consider the case of a step profile of height M. The height difference D measured by the profilometer writes:

D(x) = λ 4π arcsin A(x) sin 4πM λ . (3)
In particular the maximum of the signal D max = D(0) is corrected by a factor A(0) which is equal to unity when the two beams are fully separated and takes lower and lower values when the beams overlap. A Gaussian hypothesis for the shape of the beams leads to:

A(x) = 1 2 erf x + d/2 w -erf x -d/2 w . ( 4
)
For steps of small height M λ , linearization gives:

D(x) ≈ MA(x) = M 2 erf x + d/2 w -erf x -d/2 w . (5)
Here, A(x) is thus nothing but the linear response function of the experimental set-up to a unit step. The main effect of the response function in this experimental configuration is an underestimate of the step height M. The maximum of the signal D max = D(0) is corrected by a factor A(0) which is equal to unity when the two beams are fully separated and takes lower and lower values when the beams overlap.

Random roughness-Fourier space response function

Following Eq.( 4) the linear response function of the experimental set-up to a height Dirac is:

R(x) = A (x) = G (x + d/2) -G (x -d/2) , (6)
where G is the normalized Gaussian function of width w:

G(x) = 1 w √ π exp - x w 2 . ( 7
)
The Fourier transform of the response function thus writes:

R( f ) = 2ie -π 2 f 2 w 2 sin(π f d) , (8)
where the Fourier transform is defined as: At low frequencies, f << 1/d, one retrieves the Fourier transform of the differential operator R( f ) ≈ 2iπ f d. Note here that the attenuation factor A(0) discussed in the previous subsection for the measurements of the steps results from the low pass filter and not from a mechanism that would attenuate all frequencies. The differential nature of the measurement set-up has another important feature. In Fig. 3 we have plotted the square modulus of the Fourier transform of the response function, | R( f )| 2 , using the beam waist value w = 1.3µm (determined from experimental measurements). Larger values of w (i.e. larger laser spots in case of e.g. imperfect focusing) will decrease the sensitivity to high spatial frequencies. As can be seen in Fig. 3, | R( f )| 2 has a sharp minimum corresponding to the spatial frequency f d = 1/d, the frequency to which our measurement system is entirely blind. This frequency therefore provides a way to measure accurately the beam separation. We get f d = 0.29µm -1 for immersion measurements, the corresponding value d = 3.4µm was used throughout our analysis.

g( f ) = ∞ -∞ g(x)e 2iπ f x dx . (9)
If we denote h(x) and h( f ) the random height profile under measurement and its Fourier transform respectively, we get for the differential signal:

D( f ) = R( f ) • h( f ) = 2ie -π 2 f 2 w 2 sin(π f d) h( f ) (10)
While the integration of the differential measurement would require a simple division by the response function in the Fourier space, the blindness of the set-up at f d = 1/d makes this operation numerically hazardous. Since it involves a division by zero, integration in the Fourier space is thus only permitted in a frequency range strictly below f d = 1/d [START_REF] Gleyzes | Interferometric polarization picometric profile .1. single detector approach[END_REF]. While this operation is possible it obviously reduces the range of accessible frequencies. Another procedure consists in accounting for the response function and trying to directly analyze the bare differential signal. Let us consider a random height profile h(x). Such a random signal can be characterized by its statistical distribution and by its spatial correlation C(x). The power spectral density (PSD) of the differential signal will thus simply writes:

| D( f )| 2 = | R( f )| 2 C( f ) (11) = 4| sin(π f d)| 2 e -2π 2 f 2 w 2 | h( f )| 2 (12)
where we have

C( f ) = | h( f )| 2 .
For the sake of simplicity let us distinguish two different cases: short-range gaussian correlations and long-range correlation. In the first case the width of the distribution is given by the RMS roughness σ and the height-height correlation function C(x) = σ 2 exp(-x 2 /τ 2 ) is additionally characterized by the correlation length τ. We thus get:

C( f ) = √ πσ 2 τe -π 2 f 2 τ 2 (13) | D( f )| 2 = 4 √ πσ 2 τ| sin(π f d)| 2 e -π 2 f 2 (2w 2 +τ 2 ) (14)
Here the RMS roughness level simply gives a multiplicative factor while the gaussian correlation induces a low pass filter only measurable if the correlation length τ is significantly larger than the beam separation d. In particular at low frequency we expect in such a case to get a parabolic dependence of the PSD to the frequency:

| D( f )| 2 ≈ 4π 5/2 σ 2 τd 2 f 2 .
Let us now consider a second case of roughness that exhibits long range correlations. Glass surfaces obtained by a quench from liquid state show an ultra-low (sub-nanometric) roughness. It appears that this roughness can be quantitatively described as a superposition of frozen capillary waves [START_REF] Roberts | Ultimate low loss of hollow-core photonic crystal fibres[END_REF][START_REF] Ja | Intrinsic roughness of glass surfaces[END_REF][START_REF] Sarlat | Frozen capillary waves on glass surfaces: an AFM study[END_REF]. In the same way as the bulk structure of glass can be regarded as a snapshot of the liquid at the temperature of glass transition T G , glass surface is here regarded as a snapshot of the thermal fluctuations of the liquid interface at T G . If we call h SCW the height fluctuations induced by surface capillary waves, the power spectral density expected over 1D profiles can be expresssed as:

| h SCW ( f )| 2 = k b T G 2πγ f , (15)
where k B is the Boltzmann constant, T G the glass transition and γ the surface tension. Here k b T G /γ has the dimension of the square of a length and its square root gives a good estimate of the width of the interface. After accounting for the response function we get for the PSD of the differential signal:

| D SCW ( f )| 2 = 4e -2π 2 f 2 w 2 | sin(π f d)| 2 k b T G 2πγ f , (16)
≈ 2πd 2 k b T G γ f if f 1 d . (17)
Here we see that while a short range correlated roughness does not alter the low freqency behavior of the differential signal, the situation is very different in the case of the long-range correlated roughness induced by frozen capillary waves. The PSD shows now a linear frequency behavior.

In the following sections we make use of Eq. ( 5) and ( 16) to fit the results obtained on a deterministic silica step sample and on the inner roughness of hollow capillary silica tubes, respectively.

Tests on textured silica surfaces

In order to investigate the performance of our setup on representative test samples we performed measurements of textured silica wafers. We used 100 µm thick amorphous silica wafers (Corning 7980 Fused silica UV grade), with a standard surface finish (Optically polished both sides < 20/10 scratch/dig). Steps of known height were fabricated in these samples by use of Reactive Ion Etching (RIE) with processing time ranging from 8s to 120s. The profiles of the steps were measured by AFM (see Fig. 4). The height of the steps was found to range from 4.5 nm to 85 nm. As can be seen from the AFM image, the steps are extremely sharp with a typical transition width estimated to be about 70 nm, i.e. far below the resolution of optical profilometry.

After AFM measurements, the test samples were measured by differential profilometry i) in the conventional air-reflection configuration, and ii) in the new immersion configuration.

These test samples are planar and the height of the steps is obviously significantly larger than the sub-nanometric roughness of optical fibres. Still, the use of thin samples made of the same material (silica) as the fibres enables us to test our instrument in conditions of low sample reflectivity, in the case of oil-immersion measurements, and after transmission in the bulk material that are similar to those met in the case of fibres. and through the fit on the experimental measurements respectively, yielding d = 6.5µm and w = 2.1µm in air and d = 3.4µm and w = 1.3µm in oil immersion. The sharpness of the steps gives us the possibility to use the equation ( 5) derived in the previous section to fit the experimental measurements. An example of a differential profile is given in Fig. 4 for a 10nm high step. One obtains an excellent description of the profile with the predicted analytical expression. The parameters d and w were obtained by measuring f d AFM and air and immersion measurements were performed on a set of 5 samples textured with steps obtained by RIE exposition of respectively, 8s, 15s, 30s, 60s, and 120s. For each samples 8 AFM images, 8 differential profiles in air, and 8 profiles in immersion were performed. The results are summarized in Table 1.

For AFM measurements, 10µm × 5µm images were acquired in contact mode on an Icon (Brücker) AFM with a tip stiffness of 0.2N.m -1 and a scan rate of 1Hz. Planar fits (or flattening fits i.e. line by line along the fast scan direction) are performed on either the lower or the upper area at some distance from the step. After fitting the height histogram, two well defined peaks are identified. The distance between the two peaks is computed.

Table 1 summarizes the mean value and the standard deviation obtained over the different realizations for AFM, air and immersion differential profilometry. The different techniques give very consistent estimates. RIE etching times are roughly proportional to the actual etching depth as confirmed, with a much better accuracy, by AFM or optical profilometry. Note that even if measurements are performed on deterministic samples, i.e. steps, the results still have to be interpreted within a statistical framework. The main reason is that the RIE process is not expected to produce a perfect homogeneous step. Two main effects can be invoked here. First, the Table 1. Height of silica steps made by RIE, as measured by AFM and by differential profilometry in air and in immersion . The heights obtained by optical profilometry correspond to the parameters M providing the best fits to the experimental data using expression [START_REF] Wheeler | Wide-bandwidth, low-loss, 19-cell hollow core photonic band gap fiber and its potential for low latency data transmission[END_REF] intensity of RIE etching may spatially vary (over large length scales) and owing to the different scales of lateral resolutions of the two instruments, AFM and differential profilometry can not probe the very same area. The intrinsic roughness of the silica sample and that induced by the etching process (altogether a RMS roughnes of 2nm was estimated on the samples) induces a significant statistical uncertainty in the estimate of the step height. The results presented above with the different experimental techniques are therefore extremely consistent, and provide a cross validation of the AFM and differential profilometry.

Measurements on inner surfaces of hollow glass fibres

Direct measurements on the inner surfaces of hollow glass fibres have been performed using immersion differential profilometry. Two different samples have been studied. A hollow pure silica fibre of outer diameter (OD) 1350µm and wall thickness (WT) 88µm and another capillary tube of OD 950µm and WT 75µm made of alumina-doped silica. The first sample was produced by drawing a Suprasil F300 silica tube from Heraeus using a conventional fibre draw tower. The second sample was prepared using a conventional solution doping technique and then drawn in the same way as the first one, aside from the slightly different diameter. 11 independent profiles have been measured on each of the two samples.

During the drawing process the samples experience heating at a temperature of about 2000 • C i.e. far above the glass transition temperature T G of the materials (T G ≈ 1250 • C for silica). After cooling below T G , the roughness of the glass surface is therefore expected to result from the subsequent freezing of surface capillary waves [START_REF] Roberts | Ultimate low loss of hollow-core photonic crystal fibres[END_REF][START_REF] Ja | Intrinsic roughness of glass surfaces[END_REF][START_REF] Sarlat | Frozen capillary waves on glass surfaces: an AFM study[END_REF] and their amplitude to be controlled by the ratio T G /γ where γ is the interface tension between the viscous glass and the surrounding atmosphere at T G . The two samples have different compositions and different values of T G and γ, and are thus expected to have different thermodynamic roughnesses.

The purpose of this investigation is i) to test whether the roughness of the inner surfaces of the fibers can be quantitatively charcaterized in terms of surface capillary waves in the accessible spatial frequency range and ii) to establish whether our technique is able to resolve the difference between the two materials. In the present case the composition change (presence of alumina in one of the two samples) affects these parameters and is a potential cause of roughness contrast between the two fibres.

In Fig. 5 we show a differential profile obtained on the inner surface of a hollow silica fiber over 200µm. The height fluctuations appear to be extremely low, in the range of a few hundreds of picometers. Such an ultra-low level of surface roughness is compatible with expectations from a frozen capillary waves scenario. A more quantitative confirmation requires to characterize the Fourier spectrum of the height profile to test whether it is consistent with the spectrum of Fig. 6. Differential profilometry measurements on inner surfaces of two hollow glass fibres (silica and alumina doped silica). The power spectral densities of the differential height signal are represented as continuous lines. The dashed lines are fits obtained with the analytical expression [START_REF] Thiebaud | Nonequilibrium fluctuations of an interface under shear[END_REF]. The linear trend observed at low frequency is consistent with a scenario of frozen capillary waves at the glass surface. capillary waves. In particular, as discussed above the latter is expected to show a 1/ f behavior (where f is the spatial frequency).

Here we first show the comparison of the power spectral density of the bare differential signal | D( f )| 2 with the theoretical expectation | D SCW ( f )| 2 obtained in Eq. ( 16) for a spectrum of surface capillary waves accounting for the experimental response.

The average power spectral densities of the differential signal obtained for the silica fibre and the alumina-doped fibre are displayed in Fig. 6. The linear trend (see Eq. ( 17)) expected for frozen capillary waves is clearly visible in the range [5 • 10 -3 µm -1 -5 • 10 -2 µm -1 ]. But as Fig. 7. Roughness measurements on inner surfaces of two hollow glass fibres (silica and alumina doped silica). The power spectral densities of the height signal are obtained after division by the response function. The dashed and dotted lines are fits obtained with the analytical expression [START_REF] Derks | Suppression of thermally excited capillary waves by shear flow[END_REF]. The clear 1/ f trend observed at low frequency is consistent with a scenario of frozen capillary waves at the glass surface.

shown in Fig. 6, beyond this linear trend, the expression ( 16) accounting for the experimental response function appears to capture quantitatively the experimental power spectral densities over most of the accessible range of spatial frequencies. In practice an excellent agreement is found in the range [5 • 10 -3 µm -1 -10 -1 µm -1 ] as shown by the dashed lines in Fig. 6. At higher frequencies the signal is filtered by diffraction. To our knowledge, this is the first time that frozen capillary waves are identified in this frequency range on glass surfaces. A larger magnification may slightly enlarge the accessible frequency range toward high frequencies at the price of a reduced signal-to-noise ratio. At low frequencies (below 5 • 10 -3 µm -1 ) the signal is dominated by contributions of the nanometric scale vertical fluctuations of the motorized stage.

In order to make the expected 1/ f behavior of the roughnes spectrum even more visible we have also performed a division by the response function in the frequency range [2 • 10 -3 µm -1 -5 • 10 -1 µm -1 ]. As shown in Fig. 7, a very clear 1/ f behavior is obtained at low frequency. At higher frequency the proximity of the zero of the response function at f d = 1/d induces a high level of noise.

The power spectral density of the differential profiles obtained on the inner surfaces of the hollow glass fibers can thus be quantitatively accounted for by a superposition of frozen surface capillary waves. It is then possible to extract a single parameter T G /γ quantifying the amplitude of the surface capillary modes.

From the quantitative point of view, we obtain as a central value T G /γ = 450 ± 200 for the silica fibre and T G /γ = 1400 ± 500 for the alumina-silica fibre. The error bars have been evaluated after identifying lower and upper values bounding the fluctuations of the power spectral density in the range [5 • 10 -3 µm -1 -5 • 10 -2 µm -1 ].

The measurement clearly indicates a substantial difference in surface roughness for the two analyzed samples, and we speculate that the difference is a consequence of their different compositions, since both the temperature of glass transition and the interface tension are expected to depend on composition.

Still we remark that the two values are well below the equilibrium value estimate for amorphous silica T G /γ ≈ 5000 obtained from the standard values for the temperature fo glass tran-sition T G ≈ 1500K and surface tension γ ≈ 0.3J.m -1 . Note that low values of T G /γ had already been reported in [START_REF] Roberts | Ultimate low loss of hollow-core photonic crystal fibres[END_REF]. In this work, AFM measurements on a fragment of photonic band gap fiber led to estimate the ratio T G /γ in the range [500 -5000] with a central value T G /γ = 1500. A possible origin of this lower than expected surface roughness of glass fibers may lie in the drawing process. Recent experimental and numerical results suggest indeed that capillary waves may be attenuated under flowing conditions [START_REF] Derks | Suppression of thermally excited capillary waves by shear flow[END_REF][START_REF] Thiebaud | Nonequilibrium fluctuations of an interface under shear[END_REF].

Conclusion

We have presented an optical, non-contact method able to characterize the roughness of ultrasmooth glass surfaces with sensitivities down to 10pm, and have validated its ability to deliver quantitative height measurements on model samples. On flat glass surfaces, the capillary wave-driven sub-nanometer roughness has been obtained in a new spatial frequency range [5 • 10 -3 µm -1 -10 -1 µm -1 ] with a quantitative statistical charactarization of the spatial correlations. These measurements compare well to glass physics models as well as to AFM measurements [START_REF] Sarlat | Frozen capillary waves on glass surfaces: an AFM study[END_REF]. Beyond its excellent sensitivity, which provides a significant improvement over most current AFM performances regarding height measurements, one strength of this method resides in its ability to obtain measurements over distances of the order of hundreds of micrometers. Using an oil immersion technique, we show that these measurements can be performed at a distant glass/air interface, either on the back side of thin glass plates or, more interestingly, inside hollow glass tubes. This immersion-based method allows the study of unopened hollow fibres, thus avoiding the necessity to mechanically open the fibre and therefore preserving the pristine cleanness of the surface, which is essential considering the extremely low roughness levels measured here. The ability to access the inner surface of glass tubes, associated to the ability to measure a previously unexplored range of spatial frequencies attaining very high sensitivity presented here is expected to be of great interest for the study of roughness and loss mechanisms in photonic bandgap fibers, providing valuable input to their understanding and optimization.

Fig. 1 .

 1 Fig. 1. Polarization-modulated differential scanning interferometer. The light of the laser diode is separated into two beams carrying orthogonal polarizations by a Wollaston prism. The immersion objective is focused on the top inner interface of the capillary tube. After reflection, these beams are recombined and sent in a photoelastic polarization modulator and sent on a laser diode for lock-in detection.

Fig. 2 .

 2 Fig. 2. Current performances Noise levels measured on different kinds of samples. τ is the time constant of the lock-in amplifier.

Fig. 3 .

 3 Fig. 3. Response function | R( f )| 2 for a beam waist w = 1.3µm and a beam separation d = 3.4µm. At low frequency, the profilometer acts as a differential operator and and we get the parabolic behavior | R( f )| 2 ∝ f 2 . In addition to the diffraction effect at high frequency, the zero induced by the differential operator is clearly visible at frequency f d = 1/d = 0.29µm -1 , and can provide an experimental determination of d.

Fig. 4 .

 4 Fig. 4. Tests on silica steps obtained in air by RIE etching (exposure time 15s). Left: AFM image (estimated step height 11.7nm) Right: Differential profilomety measurement obtained at the glass/air interface on the same step (estimated step height 11.1nm after fitting with analytical expression (5), using d = 6.5µm and w = 2.0µm). The displacement step size used for the measurement was 0.18µm.

Fig. 5 .

 5 Fig. 5. Differential profile obtained on the inner surface of a silica capillary tube over a length of 200 µm. Height fluctuations are extremely low, in the range of a few hundreds of picometers.

  .

	RIE exposure time	AFM	Air	Immersion
	8s	4.7 ± 0.6 nm	4.6 ± 0.4 nm	4.3 ± 0.3 nm
	15s	12.2 ± 1.0 nm 11.2 ± 1.0 nm 10.3 ± 0.8 nm
	30s	20.9 ± 1.7 nm 22.1 ± 1.5 nm 19.1 ± 1.3 nm
	60s	45.5 ± 2.5 nm 48.3 ± 1.8 nm 46.8 ± 3.7 nm
	120s	84.3 ± 2.5 nm 79.6 ± 4.0 nm 86.8 ± 2.4 nm

Acknowledgements

This work was supported by the EU 7th Framework Programme under grant agreement 228033 (MODE-GAP). X.B. acnowledges the support of Programme DIM OxyMORE, Ile de France. The authors thank the technology center MIMENTO, Femto ST institute, for the realization of the silica steps.