Matthieu Perrault

Philippe Guéguen

Correlation between ground motion and building response using Californian earthquake records

Using data from the California Strong Motion Instrumentation Program, we studied the relationship between building response and parameters describing the noxiousness of ground motion. According to vulnerability methods that use structural drift as damage criteria, we estimated the building response on the basis of the normalized relative roof displacement (NRRD), considered as damage criteria. The relationships between the NRRD and the intensity measures of the ground motion are developed using simulated annealing method. Grouping buildings by typology (defined according to their main construction material and height) reduces the variability of the building response. Furthermore, by combining IMs, the NRRD can be predicted more accurately by a Building Damage Prediction Equation. A functional form is thus proposed to estimate the NRRD for several building typologies, calibrated on the building responses recorded in California. This functional form can be used to obtain a fast and overall damage forecast after an earthquake.

INTRODUCTION

In seismic vulnerability assessment methods, such as Hazus [START_REF] Fema | HAZUS Earthquake loss estimation methodology[END_REF], the vulnerability of a building category can be studied by assessing fragility curves. By relating a probability of damage to a ground motion intensity measure (IM), such as peak ground acceleration (PGA) or displacement response spectra (Sd), these curves account for variability due to the building model definition, the threshold damage parameter and the IMs, considered as being the engineering demand parameters providing the best representation of the noxiousness of ground motion. The relation between building response and IM can be linked to the damage prediction when associated with damage state thresholds. For example, [START_REF] Hancock | Numbers of scaled and matched accelerograms required for inelastic dynamic analyses[END_REF] proposed an equation relating earthquake parameters (magnitude, epicentral distance) and damage prediction. [START_REF] Luco | Probabilistic seismic demand analysis, SMRF connection fractures, and near-source effects[END_REF] defined two natures of IM: (1) an efficient IM defined as providing small variability of the conditional probability of damage (DM) given IM P(DM|IM); and (2) a sufficient IM defined as providing conditionally independent DM of the magnitude and distance. Coupled with real-time estimates of seismic ground motion (e.g., [START_REF] Wald | Relationships between peak ground acceleration, peak ground velocity, and modified Mercalli intensity in California[END_REF][START_REF] Worden | A revised ground motion and intensity interpolation scheme for shakemap[END_REF], [START_REF] Wald | ShakeCast: Automating and improving the use of ShakeMap for post-earthquake decision-making and response[END_REF] proposed a real-time assessment of expected damage for critical facilities. More recently, [START_REF] Krishnan | Rapid estimation of damage to tall buildings using near real-time earthquake and archived structural simulations[END_REF] estimated damage for tall buildings present in a wide area by coupling ground motion predictions with a numerical analysis of building models. Most of these studies use numerical methods to model building response, and a [building-specific\ prediction for a given scenario can provide accurate estimation of damage. However, numerical methods may require costly resources not suited for a global scale analysis. Moreover, they are often not fully representative of the response of existing buildings because of the lack of information and description of the structural design that finally require some assumptions (quality of the materials, boundary conditions, etc]).

Our study aims at assessing the correlation between efficient IMs according to the definition of [START_REF] Luco | Probabilistic seismic demand analysis, SMRF connection fractures, and near-source effects[END_REF] and building response using experimental data from a database of strong motions recorded in Californian buildings. From this database, the variability of building response is assessed through the normalized relative roof displacement (NRRD) and a combination of IMs is proposed as predictive equations of the expected response of building. These equations can be considered as a tentative of providing a building damage prediction equation (BDPE).

DATABASE OF STRONG MOTION RECORDED WITHIN BUILDINGS

The US Geological Survey (USGS) and the California Geological Survey (CGS) set up the Center for Engineering Strong Motion Data (CESMD) to provide earthquake recording data. An instrumentation program was launched in California in 1972 (CSMIP: California Strong Motion Instrumentation Program) to equip buildings, dams, bridges and empty fields with earthquake sensors. Only one set of earthquake data recorded in buildings was used in our study. These data include records of 115 ground motions from 234 buildings (Figure 1). These buildings are mainly steel frame (114) or reinforced concrete structures (76). There are also 13 masonry buildings and 8 wooden structures. Data from buildings that were built a) ISTerre, University of Grenoble, CNRS/IFSTTAR, France.

or retrofitted with earthquake isolation systems (19) were not used in this study. In order to compare with Hazuse building classification, three height categories were considered: lowrise (LR), mid-rise (MR) and high-rise (HR) buildings, for buildings comprising 1 to 3 floors, 4 to 7 floors, and more than 8 floors, respectively. The LR class counts 70 buildings, the MR class 96 buildings and the HR class 68 buildings.

The entire database comprises 772 {building/earthquake} pairs, a pair representing a building in which at least one accelerogram was recorded. The corresponding magnitudeepicentral distance diagram is represented in Figure 1. Magnitude (M W or M L as provided by the CESMD website) ranges from 3.5 to 7.2, including well-known earthquakes included such as the 02/09/1979 San-Fernando earthquake (M W = 6.6), the 10/17/1989 Loma-Prieta earthquake (M W = 6.9), the 06/28/1992 Landers earthquake (M W = 7.2) and the 01/17/1994

Northridge earthquake (M W = 6.5). The epicentral distances range from 2 km to 430 km. [START_REF] Borzi | Inelastic spectra for displacement-based seismic design[END_REF] and [START_REF] Ambrasey | Equations for the estimation of strong ground motions from shallow crustal earthquakes using data from Europe and the Middle East: horizontal peak ground acceleration and spectral acceleration[END_REF] identified the engineering significance of earthquakes for M W > 5 and a source-to-distance less than 150 km. Without additional parameters on magnitude and distance, the percentage of data respecting the magnitude (M L or M W )/distance (epicentral) criteria is 40 % (as noted on Figure 1).

The sensors are 1D accelerometers, with sampling frequency varying from 50 Hz to 200 Hz. Velocities and displacements were computed from the acceleration records, by removing the trend and mean value and filtering with an acausal Butterworth filter with 4 poles between 0.1 and 25 Hz and one or two signal integrations. The 25 Hz cut-off frequency was selected for providing a homogenous database, respecting the Nyquist condition. The peak ground accelerations recorded range from 0.017 to 8.161 m/s 2 , peak ground velocities from 0.001 to 1.119 m/s and peak ground displacements from 0.003 to 26.740 cm. In order to increase the number of data used for the prediction of damage, we considered separately the two horizontal directions of the buildings, although they are not independent.

Finally, we processed 772 pieces of data as described herein. First, considering the displacements at the top and ground floors, total relative displacement was computed and normalized with building height to obtain the Normalized Relative Roof Displacement (NRRD) considered herein as damage criterion. The correlation between damage and NRRD may be subject to caution but this parameter is easily computed with the available data in this study. Moreover, this value is a global parameter corresponding to the mean inter-storey drift observed during earthquakes and it can also be compared to the Hazus inter-storey drift limits [START_REF] Fema | HAZUS Earthquake loss estimation methodology[END_REF] given for several building typologies and four damage levels (slight, moderate, extensive and complete). The four damage states of the Hazus methodology are similar to those defined in EERI (1994), but for each Hazus building type based on the structural system. If the Hazus threshold is exceeded, a damage level is expected. In the following, we will deal with the variability between ground motion parameters and building response, considering NRRD as the damaged or undamaged criterion for characterizing building integrity.

INTENSITY MEASURES AND BUILDING RESPONSE

The variability of building responses with respect to the ground motion parameter is analyzed by testing the relationship between NRRD and standard or specific IMs. Peak Ground Acceleration (PGA), Velocity (PGV) and Displacement (PGD) are computed from the acceleration, velocity and displacement time histories, corresponding to the maximum peaks. From the acceleration time history a(t), the Arias Intensity (AI; [START_REF] Arias | Seismic design for Nuclear Power Plants[END_REF] and the Cumulative Absolute Velocity (CAV; EPRI, 1988) are defined as:

AI = π 2g a 2 (t) 0 t E dt ( 1 ) CAV = a(t) 0 t E dt ( 2 )
where t E is the total duration of the recording. Two other definitions are derived from the original CAV definition:

CAV 5 = π 2g H a(t) -0.05 ( ) a(t) 0 t E dt ( 3 ) CAV std = H PGA i -0.025 ( ) a(t) dt t i t i+1 ( ) i=1 N ( 4
)
where N is the number of consecutive 1s intervals and H is the Heaviside function (i.e., zero for negative argument and one for positive arguments). CAV 5 only takes into account amplitudes greater than 0.05 m/s 2 , and CAV std does not take into account 1s intervals with amplitudes less than 0.025 m/s 2 (EPRI, 1991).

Additional IMs, such as spectral acceleration, velocity and displacement, are defined for a given frequency and damping ratio. Different solutions exist for extracting the frequency of existing buildings. Six frequency definitions are considered herein, using empirical relationships or signal processing of the recordings at the top of the structures (Figure 2):

(1) the frequency corresponding to the peak in the Fourier transform of the top recording, noted f FFT ;

(2) the frequency corresponding to pre-seismic behavior, noted f 0,ST , and evaluated from a time-frequency distribution at the beginning of the recording. When the top and the bottom recordings are synchronized, an S-transform [START_REF] Stockwell | Localization of the complex spectrum: the S transform[END_REF] is performed on the relative displacement between the top and the bottom floors, as proposed by [START_REF] Todorovska | Earthquake damage detection in the Imperial County Services Building I: The data and timesfrequency analysis[END_REF]. If synchronization was not confirmed, the frequency assessment was carried out on the top acceleration recording, as suggested by [START_REF] Michel | Time-frequency analysis of small frequency variations in civil engineering structures under weak and strong motions using a reassignment method[END_REF].

Frequency f 0,ST corresponds to the value of the building frequency before the earthquake occurs. This frequency can be extracted directly from ambient vibrations, in the framework of identification surveys of modal properties for existing buildings, as developed by [START_REF] Michel | Quantification of fundamental frequency drop for unreinforced masonry buildings from dynamic tests[END_REF][START_REF] Michel | Seismic vulnerability assessment to slight damage based on experimental modal parameters[END_REF].

(3) the frequency corresponding to the minimal value of the time-frequency distribution, noted f min,ST . Due to the opening and closing of pre-existing cracks, building frequency shifts during earthquakes, and the more the structure is damaged, the greater the frequency reduction [START_REF] Todorovska | Earthquake damage detection in the Imperial County Services Building I: The data and timesfrequency analysis[END_REF][START_REF] Michel | Time-frequency analysis of small frequency variations in civil engineering structures under weak and strong motions using a reassignment method[END_REF][START_REF] Michel | Quantification of fundamental frequency drop for unreinforced masonry buildings from dynamic tests[END_REF].

(4) the frequency corresponding to post-seismic behavior, noted f end,ST, and evaluated from the time-frequency diagram at the end of the recording. A shift between the f 0,ST and f end,ST frequencies provides information on the structural integrity of the building (e.g. [START_REF] Dunand | Utilisation du bruit de fond pour leanalyse des dommages des bâtiments de Boumerdes suite au séisme du 21 mai 2003[END_REF][START_REF] Clinton | The observed wander of the natural frequencies in a structure[END_REF].

(5) the frequency corresponding to the empirical frequency f emp =10/N, where N is the number of floors.

(6) the frequency corresponding to the European Eurocode 8 code (CEN, 2005) and the Canadian code (NBCC, 2005), and simplified, in term of periods, as: T = 0.085H 3/4 for steel structures, T = 0.075H 3/4 for RC structures, T = 0.025H 3/4 for masonry structures and T = N/10 for the other structures. This frequency is noted f emp,EC8 .

The damping ratio was obtained by computing the impulse response of the structure (Figure 2), corresponding to the inverse Fourier transform of the ratio between the top and the bottom recordings. The damping ratio is then assessed by applying the logarithmic decrement method to the impulse response function [START_REF] Clough | Dynamics of structures[END_REF] and is proportional to the width of the function at -3dB. Some ratios present values exceeding 15 %.

Higher values are expected with higher ground motions, but these values may come from the method, which may not be sufficiently robust. [START_REF] Perrault | Using experimental data to reduce the single-building sigma of fragility curves: case study of the BRD Tower (Bucharest -Romania)[END_REF] showed the significant effect of the variation of damping values on fragility functions. Without additional information on the relevancy of this estimate, it is difficult to adopt other solution for setting damping values. However, for all buildings, the median value of the damping ratios is approximately equal to 5 % (Figure 3), which corresponds to the standard value provided in the seismic codes. In spite of the large scatter of estimated damping, as observed in Figure 3, a constant value of 5 % is assumed in this paper.

We first analyzed the ground motion parameter distributions and the NRRD distributions.

We divided PGA in 5 classes of 150 pieces of the 772 building-earthquake pairs, removing extreme values, i.e. the 11 lowest values (from the minimum values 0.017 m/s 2 to 0.032 m/s 2 )

and the 11 highest values (above 4.0 m/s 2 ), and we plotted (Figure 4) the distribution of NRRD versus PGA. We observed that, for a uniform distribution of ground motion parameters, such as PGA or PGV, the NRRD follows a lognormal distribution, or a normal distribution in a logarithmic axis. Even if the fifth class corresponding to the strongest values lacks data to have a perfect uniform distribution of PGA, this figure confirms the log-normal distribution assumption, considering to represent the vulnerability of existing buildings through fragility curves, i.e. giving the probability to observe a certain level of damage (through NRRD for example) with respect to ground motion intensity (IMs) (e.g., [START_REF] Fema | HAZUS Earthquake loss estimation methodology[END_REF][START_REF] Milutinovic | Risk-UE: an advanced approach to earthquake risk scenarios with applications to different European towns[END_REF][START_REF] Mcguire | Seismic hazard and risk analysis[END_REF][START_REF] Calvi | Development of seismic vulnerability assessment methodologies over the past 30 years[END_REF][START_REF] Michel | Seismic vulnerability assessment to slight damage based on experimental modal parameters[END_REF][START_REF] Perrault | Using experimental data to reduce the single-building sigma of fragility curves: case study of the BRD Tower (Bucharest -Romania)[END_REF]. The same result is observed whatever the IM. Thus, in order to assess the variability of the building response due to ground motion variability, we introduce two functional forms between the IM and the structural response NRRD:

ln(NRRD) = a + b ln(IM) +ε ( 5 ) ln(NRRD) = a + b ln(IM) + c IM + ε ( 6 )
where a, b and c are coefficients that will be adjusted on the basis of ground motion records and ε corresponds to residual variability. The ln function corresponds to the natural logarithm.

VARIABILITY OF BUILDING RESPONSE WITH STANDARD INTENSITY

MEASURES

We first worked on the standard IMs, considering the two functional forms for each IM (Figure 5). On this figure, we observe that PGV provides the smaller standard deviation (σ = 0.88 for the two functional forms) while PGA is less efficient, suggesting that earthquake damage statistics give a much closer correlation with the PGV than with the peak ground acceleration, already reported by [START_REF] Lesueur | Linking ground motion measurements and macroseismic observations in France: a case study based on accelerometric and macroseismic databases[END_REF] and [START_REF] Wu | Near real-time seismic damage assessment of the rapid reporting system[END_REF][START_REF] Wu | Relationship between peak ground acceleration, peak velocity, and intensity in Taiwan[END_REF]. All results are presented in Table 1, where only the standard deviations are given.

As expected, the second functional form (Equation 6) provides smaller variabilities, even if the differences compared to the variabilities provided by Equation 5 are not significant (see Table 1) (see Electronic Supplements). This is also confirmed in Figure 5, comparing the two functional forms (Equations 5 and 6) considered PGA and PGV as IMs. However, the c IM term in Equation 6 causes a difference in the slope of the curves, as seen in the Figure 5 for the higher value of PGV. This is expected for the greater IM, which may be representative of the non-linear behavior of the structural response. If we consider all the buildings, we observe significant variabilities, which are all greater than 0.88 whatever the IM used and for the two functional forms. This results in a factor exp(2σ) ≈ 6 between the ratios of the median values plus or minus one standard deviation, respectively µ+σ and µ-σ. However, this variability is strongly reduced by specifying more the building category. It is generally assumed that overall variability is composed of random and epistemic components.

Introducing more information on the structural behavior of the buildings can reduce epistemic uncertainty used for building fragility curves, as suggested by [START_REF] Spence | Comparing loss estimation with observed damage: a study of the 1999 Kocaeli earthquake in Turkey[END_REF], and confirmed by [START_REF] Perrault | Using experimental data to reduce the single-building sigma of fragility curves: case study of the BRD Tower (Bucharest -Romania)[END_REF]. For example (Table 1), when building material is introduced, the variability in wooden building response is highly reduced (e.g., = 0.49 with AI). However the reduction is only slight for the other building categories. By adding information on height, we obtain variability values lower than 0.63 for HR-RC buildings for PGV, CAV or AI. For HR steel building category according to the Hazus definition (more than 8 floors), the variability of the response remain quite important (more than 1.0), as consequence of the wide range of number of floors. The year of construction, elevation or plan irregularities could be also introduced but this information is not available for all buildings tested herein. Even if CAV 5 and CAV std are defined with thresholds, they provide standard deviations very similar to those associated with CAV. PGA provides standard deviations greater than 0.75, except for the wooden buildings category. For the six categories defined by concrete and steel materials, PGA provides a mean standard deviation equal to 0.99. PGV or AI will thus be used in order to reduce the variability in building response (see Electronic Supplement).

INTENSITY MEASURES DEFINED BY SPECTRAL VALUES

In this section, IMs are defined by the spectral values of acceleration, velocity and displacement, computed at a damping ratio of 5% and for the different frequency definitions.

In order to take into account the frequency shift during seismic loading, [START_REF] Bommer | The influence of strong-motion duration on the seismic response of masonry structures[END_REF] proposed to compute the spectral values between two periods: T 0 and F.T 0 , where T 0 corresponds to the initial period and F is a fixed parameter to be determined in order to account for the frequency drop due to the non-linear response of civil engineering structures.

The use of a frequency/period interval has already been proposed by [START_REF] Housner | Intensity of ground shaking near the causative fault[END_REF] in the definition of the Response Spectrum Intensity, defined as the average spectral value in velocity between 0.1s and 2.5s. In our case, we computed the mean values of spectral acceleration, velocity and displacement between frequencies f 0,ST and f min,ST , as follows, for Mean(Sa): Adding dynamic parameters, such as frequency and damping of the structures, allows further reduction of the relational variabilities between NRRD and IMs (Table 2 and Figure 6). The two functional forms (Equations 5 and 6) were tested and we report in Table 2 only standard deviations associated with the second functional form (Equation 6) that provide the smallest variabilities. For spectral displacements and spectral velocities, lower variabilities are obtained with frequency f = f FFT (e.g., σ = 0.36 for HR RC buildings and σ = 0.52 for MR steel buildings) and frequency f = f 0,ST (e.g., σ = 0.35 for HR RC buildings and σ = 0.51 for MR steel buildings). However, these two frequencies provide higher variabilities when Sa(f, 5%) is considered, even higher than the standard deviations given by the two empirical frequencies (Table 2). In fact, for spectral accelerations, the empirical frequencies provide reduced variabilities. For spectral velocities Sv and spectral displacements Sd, f emp , which is defined only from the number of floors, is associated with high values of variability. The Eurocode 8 empirical frequency provides values that are as low as those obtained previously with PGV or AI.

Mean Sa f min,ST ; f 0,ST [ ] ,5% ( ) [ ] = Sa( f ,5%)df f min,ST f 0,ST f 0,ST -f min,ST df +1 ( 7 )
Taking into account the frequency reduction during seismic loading over the interval [f min,ST ; f 0,ST ] results in lower variability values (Figure 6), particularly for spectral accelerations and spectral velocities (e.g., σ = 0.59 and σ = 0.52 for low-and mid-rise steel buildings, respectively, and σ = 0.34 for high-rise RC buildings). Variabilities of 0.35 provide a ratio between µ+σ and µ-σ equal to exp(2σ) ≈ 2.

In conclusion, the introduction of dynamic parameters (frequencies and damping) of the structures for IMs and identification of building class improve the relationships between ground motion and building response. The Eurocode 8 empirical frequency provides higher 6) for each class of buildings. The thick and thin dashed red lines correspond to the median +/-σ and +/-2σ, respectively.

uncertainties, equal to those obtained with standard IMs. The NRRD, which is equivalent to the mean inter-storey drift, can thus be predicted from knowledge of the building material and the number of floors of the structures. Then, for a given ground motion scenario, the average (+/-σ) value of the inter-storey drift can be computed for each building class.

Additional parameters such as the year of construction, the quality of material, the plan and elevation irregularities could be introduced to reduce more the variability of the damage prediction. Herein, to be coherent with Hazus, we decided to consider only the material and the height of the building, the two parameters that were available for all buildings, ensuring thus a maximum usable data.

BUILDING DAMAGE PREDICTION EQUATION

According to the IMs used to represent ground motion, the prediction of building response can be computed empirically considering NRRD, i.e. by predicting a certain level of damage. [START_REF] Hancock | Numbers of scaled and matched accelerograms required for inelastic dynamic analyses[END_REF] proposed a functional form in which DM is related to moment magnitude, Joyner-Boore distance and parameters dependent on fault characteristics.

However, these parameters are not available immediately after an earthquake. We therefore decided to develop a functional form providing NRRD as a function of IMs, consistent with the functional model proposed by [START_REF] Luco | Probabilistic seismic demand analysis, SMRF connection fractures, and near-source effects[END_REF] and [START_REF] Hancock | Numbers of scaled and matched accelerograms required for inelastic dynamic analyses[END_REF] for damage prediction. By combining several IMs we aim to reduce variability even further, by predicting the building response according to the formula:

g(DM) = a + f 1 (IM 1 ) +K + f i (IM i ) +K + f n (IM n ) ( 8
)
where DM is a damage measurement and IM i intensity measure. The g, f 1 , D, f n functions are the functions that represent the damage parameter and the associated IMs.

The main objective of this paper is to define empirical models for damage prediction, based on IMs. They do not require lengthy processing and within the framework of contributing to earthquake warnings, we must account for efficient IMs according to [START_REF] Luco | Probabilistic seismic demand analysis, SMRF connection fractures, and near-source effects[END_REF] definition, such as PGA or Sd. In view of the above, we chose to represent ground motion by the following IMs: PGA, CAV, Sd and the mean value of Sv computed on a frequency interval. These IMs can be obtained directly from field data and used for buildingspecific damage prediction, or predicted using ground motion prediction equations (GMPEs)

that can be integrated into warning systems. The first three parameters are equivalent to acceleration, velocity and displacement respectively, and this combination covers a large frequency band of ground motion. According to [START_REF] Wald | Relationships between peak ground acceleration, peak ground velocity, and modified Mercalli intensity in California[END_REF], these IMs are complementary, although these authors deal with the three peaks (i.e. PGA, PGV and PGD) rather than PGA, CAV and Sd. The CAV is selected because it is computed directly on the accelerogram, without integration and then faster for warning systems. The last IM used in the functional form is Sv averaged over the interval [f min,bldg ; f 0,bldg ], in order to take into account the building frequency decay during seismic loading and providing the smallest variability of building responses. For each IM, the form f i = a + b ln(IM) + c IM + ε provides lower uncertainties in the prediction of ln(NRRD), except for CAV, for which the first form (i.e. f i = a + b ln(IM) + ε) provides the same uncertainties (see Electronic Supplement).

Finally, the following functional form is considered:

ln(NRRD) = a + b PGA + c ln(PGA) + d ln(CAV ) + e Sd( f bldg ,5%) + f ln Sd( f bldg ,5%) ( ) + g Mean Sv f min,bldg -f 0,bldg [ ] ,5% ( ) +h ln Mean Sv f min,bldg -f 0,bldg [ ] ,5% ( ) + ε (9)
where f bldg is the frequency of the building, f 0,bldg corresponds to the pre-seismic frequency and f min,bldg is the minimal co-seismic frequency. ε is the term corresponding to the residuals, i.e. equivalent to the variability in building response. Coefficients a, b, c, d, e, f, g and h are assessed from the database of recordings from Californian buildings, by performing a 50,000iteration simulated annealing process. Simulated annealing is a Monte Carlo method developed to find the global minimum of challenging functions having many local minima [START_REF] Kirkpatrick | Optimization by simulated annealing: Quantitative studies[END_REF][START_REF] Chen | Optimal placement of active/passive members in truss structures using simulated annealing[END_REF]. This method requires a global minimum of a scalar function called cost function in order to converge to a local minimum. Optimization consists in finding the optimal state that minimizes the cost function describing in our case the difference between the NRRD observed and calculated, i.e. aimed at minimizing the difference between the left and right parts of Equation 9. In this study, the cost function is the standard deviation of the residuals (Res) between observed and calculated NRRD, i.e.:

Cost = Ε Res 2 [ ] -Ε Res [ ] 2 ( 10
)
where E is the expectation operator and Res = ln(NRRD obs )-ln(NRRD mod ).

We chose to use a variant of simulated annealing, the Very fast simulated Annealing [START_REF] Ingber | Very fast simulated re-annealing[END_REF]. This method generates a new state from a reference state, associated with an initial cost. The first step is to validate the inversion process. The 8 coefficients are initially set to randomly determined values a 0 , b 0 , c 0 , d 0 , e 0 , f 0 , g 0 and h 0 . A standard deviation term σ 0 is also introduced whereby the residuals ε of Equation 9 are described by a normal distribution N(0, σ 0 2 ). Using the ground motion parameters (IMs) of Equation 9(PGA, CAV, Sd and Mean(Sv)) and considering a 0 , b 0 , c 0 , d 0 , e 0 , f 0 , g 0 , h 0 and σ 0 , synthetic values of ln(NRRD) are computed by applying Equation 9. Once the run finished, the synthetics NRRD and the four observed IMs are used for assessing the coefficients a to h and σ using the simulated annealing and compared to the initial values a 0 , b 0 , c 0 , d 0 , e 0 , f 0 , g 0 , h 0 and σ 0 (Figure 7). For a slight variability of the residuals ε (σ 0 ), the coefficients are well assessed, whereas for higher σ 0 values, coefficient determination is less accurate (Figure 8). However, as shown in Figure 8i, the inversion of the variability of residuals (σ 0 ) minimizes errors of the residuals (3.6 %) whatever the variability of the original residuals. These two tests show that even if the accuracy of the inversion of coefficients a-h depends on the original variability of the data, the observed variability does not introduce trade-off in the inversion of coefficients.

In Equation 9, different frequencies are mentioned. Several definitions of frequencies were assessed in the previous part, with or without the records. Two cases are considered here:

Case 1: the real frequencies of the buildings, such as:

f bldg =f FFT and [f min,bldg ; f 0,bldg ]=[f min,ST ; f 0,ST ];

Case 2: the empirical frequencies, such as: For each building class, a series of coefficients is provided (see Electronic Supplement), and associated with variability σ, as shown in Table 3. As expected, the functional form given in Equation 9 allows a reduction of variabilities for each class. For example, for midrise RC buildings with a single IM, the lowest variability was 0.84, which was obtained with Sd(f FFT ,5%), whereas the enhanced functional form provides a variability of 0.69 and 0.76 for cases 1 and 2. As expected, the introduction of real rather than empirical frequencies gives the best standard deviations, e.g. for the high-rise concrete building class (see Figure 9).

However, for the simulation of earthquakes losses at a global scale, it is difficult to obtain experimental frequencies. The use of empirical frequencies is a relevant solution that does not require building recordings. Moreover, the variabilities associated with empirical frequencies present low values, which can be lower than with real frequencies (e.g. for all buildings in the masonry building class and for high-rise steel buildings).

CONCLUSION

On the basis of the CSMIP database of earthquake recordings from Californian buildings, the correlation between building response, described by the normalized relative roof displacement (NRRD), and intensity measures IMs, such as PGA, PGV, CAV and Sd, was studied. It appears that classifying the buildings by typology may reduce the epistemic variability of building response. Moreover, IMs that use the dynamic parameters of buildings provide the lowest standard deviations (e.g., Sd(f FFT ) provides a mean standard deviation σ mean = 0.71 for the six building classes defined according to construction material and when used alone (Tables 1 and2) to describe ground motion and when combined according to Equation 9 (two last columns), using real or empirical frequencies. Variabilities are also given when IMs are added one by one in Equation 9. The sigma values of the fifth column (Sd) were obtained with f bldg= f FFT , and those from the sixth column (Mean(Sv)) with [f min,bldg ; f After testing several solutions, we find that the functional form using several IMs (Equation 9) gives the smallest variability of results, especially when the residual distribution is analyzed (Figure 9). This functional form and its coefficients assessed by simulated annealing are only compatible with the data used in this paper. These data are for California buildings, mainly in steel and reinforced concrete. The data correspond to earthquakes of magnitude between 3.5 and 7.0, and epicentral distances between 2 and 400 km. Rather than speaking of earthquake parameters, we believe that it is preferable to set the domain of validity of the empirical model for drift between 2.10 -6 and 2.10 -2 m/m. This model must be confirmed, especially by testing it on other types of construction. It can also be refined by reducing the uncertainties through the introduction of additional information for the building classification (e.g. date of construction, design, etc ...). Moreover, building-specific analysis (i.e., considering one building and several earthquakes) could improve information on the uncertainties related to the IM versus drift relationships, the number of these data being by now to small for carrying out an extensive analysis. When associated with empirical frequencies, it provides good estimates of NRRD (σ mean = 0.66 when empirical frequencies are used), which can be used in an earthquake early warning system, by comparing estimated NRRD with reference values for inter-storey drift defined for different damage levels (e.g. [START_REF] Fema | HAZUS Earthquake loss estimation methodology[END_REF]. Since empirical frequencies depend on the building height, this information can be easily filled from field survey, national census or remote sensing. Moreover, the damage prediction equation can also be used to build fragility curves for the associated building classes [START_REF] Perrault | Évaluation de la vulnérabilité sismique de bâtiments à partir de mesures in situ[END_REF].

AKNOWLEDGMENT

Figure 1 .

 1 Figure 1. Left: Map of California showing building location (squares) and earthquake epicenters (stars) from the CSMIP database. Right: magnitude-distance diagram corresponding to the entire database. The black rectangle contains 306 data (on a total of 772 data), which correspond to the criteria: M 5.0 and d 150 km.

Figure 2 .

 2 Figure 2. Left: Example of an S-transform[START_REF] Stockwell | Localization of the complex spectrum: the S transform[END_REF] performed on the 1971 San Fernando earthquake recording at the top floor of the Millikan Library building. Right: Impulse response used to assess the damping ratio.

Figure 3 .

 3 Figure 3. Damping ratios of the Californian buildings assessed by applying the logarithmic decrement method to each building impulse response.

Figure 4 .

 4 Figure 4. NRRD distributions according to PGA distributions for five categories of PGA (from left to right, and from top to bottom): [0.032 ; 0.106], [0.106 ; 0.513], [0.514 ; 0.803], [0.808 ; 1.266] and [1.266 ; 3.902] m/s 2 . Each category includes 150 pieces of data from the Californian buildings database.

Figure 5 .

 5 Figure 5. Comparison of the two functional forms (Equations 5 and 6) in order to represent NRRD as a function of intensity measures. The solid lines correspond to the median value, and the thick and thin dashed lines correspond to the median +/-σ and +/-2σ, respectively. Examples with PGA (first line) and PGV (second line). The standard deviations corresponding to each form are given in each title.

Figure 6 .

 6 Figure 6. Normalized Relative Roof Displacement (NRRD) as a function of Mean(Sd) evaluated on the frequency interval [f min,ST ; f 0,ST ], for the following typologies: all buildings, steel buildings, and the three steel classes defined with the number of floors. σ corresponds to the standard deviations provided by the second form (Equation6) for each class of buildings. The thick and thin dashed red lines correspond to the median +/-σ and +/-2σ, respectively.

Figure 7 .

 7 Figure 7. Example of the simulated annealing processing, to test the inversion of coefficients a and b for a given value of variability of residuals ε. The values a 0 and b 0 are initially set (red line) and 50,000 iterations of the simulated annealing process are used for converging. Red dots correspond to high values of the cost function and the more the dots are blue, the more the cost function is low.

Figure 8 .

 8 Figure 8. Accuracy of the simulated annealing process for the inversion of coefficients a to h of Equation 9 (label a to h in the subplots), following the initial variability σ 0 introduced to represent the distribution of the normally distributed residuals ε (i). The histogram corresponds to the distribution of the σ values (differences between initial and final value in %) (σ of Equation9).

Figure 9 .

 9 Figure 9. Left: comparison of observed NRRD as a function of synthetic NRRD, computed using Equation 9 and introducing real frequencies. Right: associated residuals, described by a normal distribution. Example for all buildings (first line) and high-rise concrete buildings (second line). Standard deviations are given in the titles. The thick and thin dashed red lines correspond to the median +/-σ and +/-2σ, respectively.

  This work was supported by the ANR (National Research Agency) as part of its RiskNat program (URBASIS project, no. ANR-09-RISK-009). Data were provided by the California Strong Motion Instrumentation Program (CSMIP), through the CSMIP web data access (http://strongmotioncenter.org/).

Table 1 .

 1 Correlation between building response, represented by NRRD, and intensity measures: PGA, PGV, PGD, CAV, CAV 5 , CAV std and AI. The values are the standard deviations associated with the first functional form (Equation5) and with the second functional form (Equation6). The numbers in brackets indicate the number of pieces of data available for each category of buildings. In each column, i.e. for each building class, the bold value indicates the smallest standard deviation. LR: lowrise buildings; MR: mid-rise buildings; HR: high-rise buildings.

		Eq.	All	Steel	Conc.	Mason.	Wood	LR Steel	MR Steel	HR Steel	LR Conc.	MR Conc.	HR Conc.
			(772)	(312)	(317)	(47)	(13)	(86)	(81)	(145)	(71)	(116)	(130)
	PGA	5 6	1.101 1.100	1.160 1.158	0.957 0.956	1.241 1.180	0.698 0.605	0.763 0.759	0.998 0.987	1.375 1.374	1.023 1.018	0.947 0.924	0.900 0.899
	PGV	5 6	0.889 0.888	1.835 0.826	0.823 0.823	1.052 1.046	0.544 0.521	0.690 0.690	0.760 0.701	0.776 0.767	0.970 0.970	0.878 0.872	0.616 0.615
	PGD	5 6	0.945 0.944	0.846 0.846	0.931 0.929	1.213 1.213	0.695 0.677	0.872 0.867	0.809 0.793	0.657 0.653	1.023 1.017	0.936 0.928	0.802 0.789
	CAV	5 6	0.960 0.960	0.897 0.897	0.850 0.850	1.199 1.186	0.526 0.518	0.782 0.776	0.826 0.813	0.804 0.795	0.934 0.932	0.895 0.889	0.627 0.627
	CAV 5	5 6	0.991 0.969	0.932 0.907	0.905 0.862	1.065 1.044	0.520 0.503	0.778 0.720	0.813 0.809	0.904 0.889	0.970 0.956	1.004 0.916	0.711 0.679
	CAV std	5 6	0.955 0.955	0.852 0.846	0.875 0.868	1.009 1.009	0.547 0.545	0.744 0.716	0.874 0.862	0.775 0.739	1.023 1.018	0.971 0.958	0.691 0.685
	AI	5 6	0.924 0.921	0.890 0.877	0.803 0.803	1.131 1.131	0.491 0.489	0.666 0.663	0.797 0.765	0.913 0.897	0.921 0.913	0.826 0.822	0.632 0.629

where df is the discretisation frequency of Sa, i.e. Mean(Sa) is the arithmetic mean of Sa over the interval [f min,ST ; f 0,ST ].

Table 2 .

 2 Correlation between the building response, represented by NRRD, and the intensity measures defined from spectral acceleration, velocities and displacements. The values are the standard deviations associated with the second functional form (Equation6). The numbers in brackets indicate the number of pieces of data available for each class of buildings. In each column, i.e. for each building class, the bold value indicates the smallest standard deviation. LR: low-rise buildings; MR: mid-rise buildings; HR: high-rise buildings.

			All	Steel	Conc.	Mason.	Wood	LR Steel	MR Steel	HR Steel	LR Conc.	MR Conc.	HR Conc.
			(772)	(312)	(317)	(47)	(13)	(86)	(81)	(145)	(71)	(116)	(130)
		Sa	1.334	1.243	1.193	1.243	0.458	0.775	0.981	1.199	1.076	1.233	0.994
	f FFT	Sv	0.957	0.777	0.847	1.128	0.671	0.597	0.648	0.807	1.059	0.906	0.513
		Sd	0.971	0.827	0.895	1.221	0.810	0.622	0.519	0.792	1.139	0.836	0.361
		Sa	1.285	1.178	1.154	1.188	0.473	0.756	0.930	1.124	1.052	1.186	0.942
	f 0,ST	Sv	0.951	0.779	0.845	1.153	0.714	0.604	0.618	0.835	1.061	0.912	0.470
		Sd	0.980	0.837	0.909	1.243	0.856	0.647	0.508	0.798	1.131	0.876	0.348
		Sa	1.351	0.815	1.233	1.217	0.469	0.803	1.043	1.178	1.093	1.256	1.056
	f min,ST	Sv	0.952	0.792	0.837	1.109	0.675	0.607	0.682	0.821	1.020	0.881	0.546
		Sd	0.950	1.247	0.867	1.221	0.812	0.622	0.538	0.773	1.085	0.800	0.385
	[f min,ST ;	Sa	1.316	1.210	1.188	1.192	0.460	0.777	0.986	1.141	1.068	1.227	0.989
	f 0,ST ]	Sv	0.949	0.773	0.842	1.119	0.683	0.593	0.646	0.808	1.045	0.900	0.492
	Mean:	Sd	0.964	0.821	0.887	1.224	0.824	0.623	0.517	0.778	1.113	0.836	0.341
	f emp = 10/N	Sa Sv Sd	1.188 1.085 1.263	1.101 0.865 1.058	1.043 1.205 1.407	1.446 1.437 1.613	0.778 0.859 0.744	0.759 0.769 0.713	0.734 0.812 0.800	0.811 0.726 0.783	1.126 1.549 1.637	1.012 0.954 0.922	0.666 0.561 0.652
		Sa	1.134	0.950	0.924	1.285	0.599	0.752	0.713	0.657	0.880	0.937	0.679
	f emp,EC8	Sv	0.879	0.761	0.845	1.202	0.611	0.625	0.711	0.702	0.993	0.879	0.609
		Sd	0.941	0.792	0.989	1.156	0.577	0.616	0.683	0.703	1.109	0.843	0.684

Table 3 .

 3 Variabilities associated with IMs PGA, CAV, Sd(f bldg ,5%) and Mean(Sv[f min,bldg ; f 0,bldg ],5%))

  0,bldg ]=[f min,ST ; f 0,ST ]. LR: low-rise buildings; MR: mid-rise buildings; HR: high-rise buildings. See Table3). The combination of four IMs, i.e. PGA, CAV, Sd and Mean(Sv) for a predictive equation of the expected response derived from experimental data provides lower standard deviations in the building response (σ mean = 0.61 when real frequencies are used). Since structural drift is considered, this equation can be considered as the first step towards the building damage prediction equation (BDPE), considering NRRD as damage criteria according to Hazus methodology.

			Previous part results	by using only Equation 9 results	Equation 9 results
		PGA CAV Sd	Mean Sv	PGA	PGA & CAV	PGA & CAV & Sd(f FFT )	Real frequencies	Empirical frequencies
	1	All	1.10 0.96 0.97 0.95 1.10	0.93	0.87	0.86	0.83
	2	Steel	1.16 0.90 0.83 0.77 1.16	0.88	0.77	0.75	0.70
	3	Concrete	0.96 0.85 0.90 0.84 0.96	0.80	0.76	0.76	0.78
	4	Masonry	1.18 1.19 1.22 1.11 1.18	1.08	1.07	1.03	0.86
	5	Wood	0.61 0.51 0.81 0.68 0.61	0.50	0.49	0.31	0.39
	6	LR Steel	0.76 0.78 0.62 0.59 0.76	0.68	0.57	0.56	0.54
	7	MR Steel	0.99 0.81 0.52 0.65 0.99	0.78	0.52	0.50	0.64
	8	HR Steel	1.37 0.80 0.79 0.81 1.37	0.79	0.72	0.70	0.61
	9	LR Concrete	1.02 0.93 1.14 1.05 1.02	0.89	0.88	0.88	0.87
	10 MR Concrete	0.92 0.89 0.84 0.90 0.92	0.80	0.75	0.69	0.76
	11 HR Concrete	0.90 0.63 0.36 0.49 0.90	0.61	0.35	0.35	0.54
		σ mean : mean on classes 6 to 11	0.99 0.81 0.71 0.75 0.99	0.76	0.63	0.61	0.66
		number of floors.


Mathilde Rieu
email: rieu@emse.fr

Malick Camara

Guy Tournier

Jean-Paul Viricelle

Christophe Pijolat

Nico F De Rooij

Danick Briand

Fully inkjet printed SnO 2 gas sensor on plastic substrate

Keywords: SnO 2 sensor, Sol-gel synthesis, Inkjet printing, Polyimide foil

A tin dioxide (SnO 2 ) sensor was fabricated by inkjet printing onto polyimide foil. Gold electrodes and heater were printed on each side of the substrate. A SnO 2 based ink was developed by sol-gel method and jetted onto the electrodes. A final annealing at 400°C compatible with the polymeric transducers allowed to synthetize the SnO 2 sensing film.

Electrical measurements were carried out to characterize the response of the fully printed sensors under oxidizing and reducing gases. The device was heated up at a temperature between 200 and 300°C using the integrated heater. The proper operation of the full printed metal-oxide gas sensors was validated under exposure to carbon monoxide and nitrogen dioxide, in dry and wet air.

Introduction

Development of sensing devices onto flexible foils is receiving special attention with the advances made in printed electronics [START_REF] Mattana | Recent advances in printed sensors on foil[END_REF][START_REF] Briand | Making environmental sensors on plastic foil[END_REF][START_REF] Oprea | Temperature, humidity and gas sensors integrated on plastic foil for low power applications[END_REF][START_REF] Khan | Technologies for Printing Sensors and Electronics Over Large Flexible Substrates: A Review[END_REF]. Indeed, the development of sensors on plastic substrates is of growing interest, with the aims of achieving low-cost, flexible, and wireless devices [START_REF] Vasquez Quintero | Smart RFID label with printed multisensor platform for environmental monitoring[END_REF].

Compared to silicon technology, fabrication on plastic substrates can simplify sensors manufacturing by avoiding bulk micromachining and high temperature process steps required, for instance, in the realization of silicon micro-hotplates. However, it demands the development of material synthesis and deposition procedures compatible with the temperature plastic supports.

Previous work reporting on the realization of metal-oxide gas sensors onto plastic substrate is based on clean room fabrication processes [START_REF] Courbat | Drop-coated metal-oxide gas sensor on polyimide foil with reduced power consumption for wireless applications[END_REF][START_REF] Vallejos | Microfabrication of flexible gas sensing devices based on nanostructured semiconducting metal oxides[END_REF]. Nowadays, additive manufacturing using printing is of interest to reduce the infrastructure and simplify the sensors fabrication procedure. Among the different printing techniques, ink-jet printing is a contactless technique able to deliver digitally precise quantities of materials to form well-defined patterns. This technique is used for the preparation of different functional films [START_REF] Mattana | Recent advances in printed sensors on foil[END_REF][START_REF] Molina-Lopez | All additive inkjet printed humidity sensors on plastic substrate[END_REF][START_REF] Claramunt | Flexible gas sensor array with an embedded heater based on metal decorated carbon nanofibres[END_REF][START_REF] Danesh | Development of a New Generation of Ammonia Sensors on Printed Polymeric Hotplates[END_REF][START_REF] Courbat | Polymeric foil optical waveguide with inkjet printed gas sensitive film for colorimetric sensing[END_REF].

Until now, inkjet printed hotplates on polymeric foil were limited to an operating temperature of about 100°C due to the use of silver ink [START_REF] Claramunt | Flexible gas sensor array with an embedded heater based on metal decorated carbon nanofibres[END_REF][START_REF] Danesh | Development of a New Generation of Ammonia Sensors on Printed Polymeric Hotplates[END_REF]. However, metal-oxide gas sensors normally require higher operating temperatures and more inert materials, such as gold or platinum, must be implemented.

Regarding inkjet printing of metal-oxide sensing layers, the works published deal with the deposition onto rigid substrates, made of silicon or alumina [START_REF] Shen | Properties of SnO 2 based gas-sensing thin films prepared by ink-jet printing[END_REF][START_REF] Vilà | Influence of In and Ga additives onto SnO2 inkjet-printed semiconductor[END_REF][START_REF] Andio | Comparison of gas sensor performance of SnO 2 nano-structures on microhotplate platforms[END_REF][START_REF] Shen | The preparation of ZnO based gas-sensing thin films by ink-jet printing method[END_REF]. The ink can be composed of either precursors or nanoparticles to prepare the oxide film, but the sintering step occurred at more than 500°C, which is not compatible with the use of plastic foil.

We report here on the preparation and the characterization of fully inkjet printed SnO 2 sensors with integrated heater onto polyimide (PI) film able to operate at temperatures up to at least 300°C. The hotplate transducers were made of inkjet printed gold electrodes and heater. SnO 2 was chosen as gas sensitive layer as it has been widely used as gas-sensing material since many years and plays a central role in the field of gas sensors [START_REF] Das | SnO 2 : A comprehensive review on structures and gas sensors[END_REF][START_REF] Kohl | Function and applications of gas sensors[END_REF][START_REF] Riviere | Development of tin oxide material by screenprinting technology for micro-machined gas sensors[END_REF][START_REF] Pijolat | Hydrogen detection on a cryogenic motor with a SnO 2 sensors network[END_REF]. A custom made ink of tin oxide was formulated and its synthesis using a sol-gel process was studied. The gold transducers on polyimide foil were coated via inkjet printing with the custom made ink solution and, finally, the electrical response of the sensor to different gases was characterized.

Experimental

Substrate preparation

A flexible foil was used as substrate. Polyimide (PI) Upilex-50S, 50 µm in thickness, was selected because of its heat resistance property, with the capability of preserving its physical and chemical properties at high temperature (up to 500°C). Prior to the printing step, the substrate was cleaned by successive immersion in acetone, isopropanol and deionized water.

After the cleaning process, it was dried into an oven at 120°C for 1h. The substrate could be optionally subjected to an oxygen plasma treatment (Diener Zepto Plasma reactor) for 2 min at 100 W (40 kHz) to modify its surface energy.

Design and printing of the sensors

The transducer consists of gold electrodes on one side of the PI foil and of a gold heater on the other side as shown in Figure 1. Electrodes were designed to have a simple two parallel lines structure (Fig. 1 b). The design of the heater (Fig. 1 c) was not optimized to minimize the power consumption of the device but to achieve a high printing yield.

Commercial gold ink (NPG-J from Harima) was ink jetted onto the PI substrate. The printer used was a Dimatix DMP-2800 with 10 pL cartridges. After printing, gold coatings were placed in an oven at 250°C for 3 hours to sinter the ink. SnO 2 ink was then deposited onto the electrodes. The shape of the sensing element was a rectangle of 4 x 2 mm 2 in order to fully cover the electrodes and the space in between.

SnO 2 ink preparation

The SnO 2 material was prepared using sol-gel process as described in the literature [START_REF] Sujatha Lekshmy | SnO 2 thin films doped indium prepared by the sol-gel method: structure, electrical and photoluminescence properties[END_REF][START_REF] Zhang | Facile synthesis of highly ethanol-sensitive SnO 2 nanoparticles[END_REF][START_REF] Bin | Investigation of the hydrogen response characteristics for sol-gel-derived Pd-doped, Fe-doped and PEG-added SnO 2 nano-thin films[END_REF].

SnCl 2 was dissolved in pure ethanol with a concentration of 0.25 mol.L -1 . The mixture was stirred and kept at 80°C for 5 h in a closed vessel to form a transparent tin ethoxide solution (Sn(OEt) 2 ). The obtained sol was clear and homogeneous. Furthermore, the sol was stored in a closed vessel and was stable over some months if keeping out of contact of humidity.

To inkjet the sol, viscosity and surface tension have to be adjusted to fulfil the requirements of the inkjet printer. For the printer used, the range of viscosity and surface tension are 10-12 mPa.s and 30-35 mN.m -1 , respectively. Inks viscosity and surface tension measurements were carried out with a viscosimeter DV-II+ by Brookfield and a portable tensiometer AquaPi by Kibron, respectively. Sol viscosity was 1.5 mPa.s and its surface tension was 22 mN.m -1 . In order to adjust these properties, ethylene glycol (EG) was added to the sol. Ethylene glycol was chosen because of its boiling temperature point near 200°C, which is higher than the 100°C recommended for inkjet printing in order to prevent quick evaporation and blocking of the nozzle. By adding EG to the sol, viscosity and surface tension were increased. Various amounts of sol and EG were mixed. The solution that fulfil the best the requirements of the inkjet printer was composed of 25 vol% of sol and 75 vol% EG with a viscosity of 10.2 mPa.s and surface tension of 35 mN.m -1 .

The SnO 2 ink was then put into a cartridge. Dimatix software was used to set the printing parameters. The jetting waveform used for the deposit was adapted from a model given into the Waveform Editor software of the printer, with a firing voltage set at 20V and a cartridge temperature fixed at 30°C.

Characterization techniques

Thermal analyses of the xerogel, thermo-gravimetric analysis (TGA) and differential scanning calorimetry (DSC), were performed with a coupled TGA/DSC Star e System from Mettler Toledo. Analyses were performed from room temperature to 800°C under 5 L h -1 air flow.

X-ray diffraction (XRD) analysis was performed onto powder with a Siemens D5000 diffractometer using Cu Kα radiation at 40 kV and 40 mA. Environmental scanning electron microscopy (ESEM Philips 30XL ESEM-FEG) was used to investigate the surface morphology of the coatings. Transmission Electron Microscopy (TEM Philips CM 2000) was used to characterize the SnO 2 powder and the cross section of the sensor. In order to prepare the sensor's cross section, the lamella was prepared using a Focus Ion Beam (FIB) equipment (Equipex MANUTECH-USD). Carbon and platinum were deposited on the top of the SnO 2 layer during the preparation of the lamella.

The gas sensing performances of the sensors were evaluated thanks to a gas mixing system equipped with mass flow controllers. Carbon monoxide and nitrogen dioxide, gases representative of atmospheric pollution field, were used as reducing and oxidizing gases, respectively. It was possible to generate concentrations from 1 to 50 ppm of CO and from 0.6 to 20 ppm of NO 2 , either in dry air or wet air (1 vol % absolute humidity, which corresponds to 50% relative humidity at 20°C) as carrier gas. The sensor was put in a cell exposed to a constant gas flow of 10 L h -1 . Thanks to a power supply, the gold heater at the backside of the sensor allowed to heat up the tin dioxide material between 200 and 300°C. The electrical responses of the sensors were measured using a multimeter (HP 34410A). Resistance data were recorded as a function of time. Chemo-resistive gas responses, defined as R air /R CO or R NO2 /R air , of the sensor under dry and wet air were calculated to compare the printed sensors characteristics to those of sensors reported in the literature. The chemo-resistive responses mean value of two tested devices and the standard deviations were calculated and used to plot calibration curves. The limit of detection (LOD) was also calculated and is defined as 3 times the standard deviation of the background noise. Unfortunately, sensors electrical responses to NO 2 injections of 5 and 20ppm were not measurable at 200°C (out of range R > 10 8 ohm).

Results and discussion

Characterization of SnO 2 precursor

Before establishing the inkjet printing process of the solution, the sol was hydrolyzed to obtain a gel in order to characterize it. By adding water (50 vol %) to the sol and letting the mixture 24h into a close bottle vessel at 80°C-4h in an oven, the transformation sol to gel occurred. The inorganic network of the gel retains a large amount of liquid phase composed of water and organic solvent. This gel was then dried at 80°C for 24h to obtain a xerogel.

Changing from a humid gel to dry xerogel powder, the major part of the liquid phase was removed. Then, the decomposition in air of the xerogel was investigated by thermogravimetric analysis (TGA) and differential scanning calorimetry (DSC). Weight loss and heat flow are reported as a function of the temperature in Fig. 2. The TGA curve shows, with increasing temperature, a continuing loss of mass linked to the decomposition of the organic compounds retained by the inorganic network of the xerogel. Around 400°C, the total loss of weight constitutes 50 wt % of the initial xerogel weight. The shape of the DSC curve indicates endothermic effects in the low temperature range (20-300°C), which are most frequently attributed to organics decomposition. At 375°C, a sharp exothermic process occurs simultaneously with the highest weight loss. As crystallization process is generally associated with an exothermic peak, this can be associated to SnO 2 crystallization [START_REF] Chepik | Crystallization of SnO 2 produced by sol-gel technique from salts of tin in different oxidation states[END_REF][START_REF] Ravaro | Effect of pH of colloidal suspension on crystallization and activation energy of deep levels in SnO 2 thin films obtained via sol-gel[END_REF].

To confirm the crystallization of the compound at around 400°C, X Ray Diffraction (XRD) analysis was done onto SnO 2 powder annealed at 400°C-1h. The powder was crystallized in the tetragonal structure of SnO 2 . The average crystallites size was calculated using Scherrer equation and was 11 ± 1 nm. The average crystallite size was confirmed by transmission electron microscopy (TEM) done on the powder annealed at 400°C. Indeed, as it can be observed on the micrograph reported in Fig. 3, the SnO 2 powder is composed of crystallites of around 10 nm of diameter. These crystallites are agglomerated into larger particles. Fig. 3. TEM micrograph of the SnO 2 powder after annealing at 400°C-1h.

Inkjet printing of the whole sensor onto polyimide substrate 3.2.1. Printing gold transducers

A commercial gold ink was printed to pattern the electrodes and heater onto cleaned and dried PI foil without further surface treatment. The inkjetted spots exhibited a diameter of around 150µm. Thus, a drop to drop spacing of 100 µm was chosen to obtain a uniform gold layer. Three layers were successively deposited in order to get an electrically conductive track. Gold electrodes were first printed on one side of the PI foil and annealed at 250°C for 3h. The heater was then printed on the other side of the substrate and sintered also at 250°C-3h.

Pictures of electrodes and heater are reported in Fig. 1 b) and c), respectively. The resistance of the heater at room temperature (22°C) was 39 Ω. Its temperature coefficient resistance (TCR) was measured between 20 and 200°C thanks to a climatic chamber. The TCR of the printed gold layer was of 2.0 10 -3 K -1 and was used to set the temperature of the sensor when injecting power into the heater.

SnO 2 ink deposition

The reservoir of the printhead was filled with the SnO 2 ink and inkjet printer parameters (cartridge temperature, jetting waveform) were adjusted in order to obtain a drop without satellite.

Then, the effect of the surface treatment of the substrate and of the drop to drop space (DDS) on the printing pattern was evaluated. As reported in Figures 4 a) andc), without a surface treatment of the PI film, the deposited drops were around 30-40 µm in diameter; in this case, a drop to drop space of 20 µm was necessary to get a good coverage of the surface (Fig. 4 c).

When performing an oxygen plasma treatment onto the PI surface (Fig. 4 b andd), the deposited drops were spread out with a diameter of 160 µm and a drop to drop space of 100 µm was sufficient to get the coverage (Fig. 4 b). However, the amount of deposited matter was very low in this last case, so making it difficult to see the deposited films as they were almost transparent. Finally, the parameters used for the preparation of the sensor were no surface treatment (only cleaning) and a DDS of 20 µm in order to have a continuous visible deposit (Fig. 4 c). To complete the sensor preparation, the SnO 2 ink was printed onto the substrate and the gold electrodes. Four layers were deposited in order to get a visible layer of SnO 2 material. A dry step at 50°C-30 min was implemented between each printed layer. Finally, a thermal annealing at 400°C-1h was applied to form the SnO 2 film.

Microstructural characterization of the sensing layer

After annealing, the coating was transparent and difficult to observe. SEM micrograph of the sensor is reported in Fig 5 . The PI substrate under the active layers is black on the micrograph. The two gold electrodes in grey in the micrograph are 500 µm-wide and spaced by 600 µm. In between the gold electrodes, the SnO 2 coating was transparent (black color) except where the material got accumulated, which resulted in the grey spots shown in Figure 5. Therefore the SnO 2 coating was not homogeneous with some non-transparent amounts of matter locally accumulated. The cross section of the sensor was then characterised by TEM analysis, in the area where the transparent coating was located, which means between the two gold electrodes and the grey spots. For this, the lamella was cut thanks to a FIB equipment. The PI substrate and the thin SnO 2 layer can be seen in the micrographs presented in Fig. 6. The SnO 2 layer thickness was determined to be between 10 and 15 nm, which means that this transparent coating is made of only a monolayer of SnO 2 crystallites. Finally, the electrical responses of the sensor under gas exposure were evaluated. The integrated heater at the back side was used to adjust the sensors temperature. A power supply of 0.59 Watt and 0.26 Watt was necessary to heat the sensor up to 300°C and 200°C, respectively. Temperature homogeneity was checked thanks to a K-type thermocouple and a difference of less than 20°C was found all over the sensing element zone which was acceptable. Two sensors prepared using the same procedure were tested and leaded to the same results. Furthermore, the sensors were tested for more than 30 days with the observation of a slight drift of their baseline but the amplitude of the sensor response remained the same.

Several gas concentrations of CO and NO 2 were injected, either using dry or wet air as carrier gas. The electrical conductance of the sensor versus time is plotted in Fig. 7, with injections of CO and NO 2 in dry air, at 300°C (Fig. 7 a); 200°C (Fig. 7 b) and in humid air at 300°C (Fig. 7.c). Under dry air, one can notice that the pure SnO 2 sensor responds to the gases even at low concentrations, such as 1 ppm CO and 0.6 ppm NO 2 . Obviously, the reaction kinetics are slower at 200°C compared to 300°C, which implies longer response time and longer recovery time at 200°C. Comparing dry and wet atmosphere, at 300°C (Fig. 7. a) and c)), the main difference is in the baseline. Indeed, water is known to increase the conductance [START_REF] Bârsan | Understanding the fundamental principles of metal oxide based gas sensors; the example of CO sensing with SnO 2 sensors in the presence of humidity[END_REF].

That also implies that the responses to low concentration of CO are largely hidden in the case of wet air. Chemo-resistive gas responses (R air /R CO or R NO 2 /R air ) of the sensor under dry and wet air were extracted from these measurements and were plotted versus gas concentration in Fig. 8. These calibration curves evidence that the presence of the humidity drastically reduced the sensor response to CO gas (Fig. 8.a). Responses were smaller under wet air, nearly divided by 2. In the case of NO 2 (Fig. 8.b), the presence of humidity reduced also the responses but they were still visible. The influence of temperature was less pronounced than humidity effect. The chemo-resistive responses at 300°C were slightly lower than at 200°C but lowering the temperature affected mainly the kinetics as mentioned previously. Finally, the responses to CO and NO 2 in the investigated temperature range were comparable to result reported in literature for SnO 2 sensors [START_REF] Courbat | Drop-coated metal-oxide gas sensor on polyimide foil with reduced power consumption for wireless applications[END_REF][START_REF] Shen | Properties of SnO 2 based gas-sensing thin films prepared by ink-jet printing[END_REF][START_REF] Hahn | CO sensing with SnO 2 thick film sensors: role of oxygen and water vapour[END_REF]. Under humid atmosphere and for both tested temperatures, 20 ppm of CO was the lowest concentration to obtain exploitable responses of the sensor compared to 0.6 ppm for NO 2 . Limits of detection (LOD) were also calculated and are reported in Table 1. It confirms that the performances of pure SnO 2 are weaker under humid atmosphere than in dry air especially for CO injections. For example, at 300°C, under dry air, less than 1 ppm CO would be detectable while under wet air LOD was higher than 20 ppm. The LOD for NO 2 injections were quite low, a few ppb, which is agreement with literature data [START_REF] Stanoiu | Low level NO 2 detection under humid background and associated sensing mechanism for mesoporous SnO 2[END_REF].

Conclusions

We have demonstrated the first fully inkjet printed metal-oxide sensor onto polymeric foil, including the gas sensitive layer and heating transducer. The reported technology allows the simple preparation of tin oxide sensors with results comparable to what has been published so far in the literature on this type of sensors. We have reported on a sol-gel tin oxide film compatible with inkjet printing and the use of a polyimide substrate thanks to its synthesis temperature of 400°C. The different fabrication steps of the sensor were characterized and the gas sensing performances were evaluated. Further optimization work on the sensor is necessary, with the aims of lowering its power consumption and improving its sensing performance.

Fig. 1 .

 1 Fig. 1. Cross section schematic view of the device (a) and photographs of the top side (electrodes and SnO 2 ) of the device (b) and of the back side (heater) (c).

Fig. 2 .

 2 Fig. 2. Thermo-gravimetric analysis (weight loss) and differential scanning calorimetry (heat flow) of SnO 2 xerogel as a function of the temperature.

Fig. 4 .

 4 Fig. 4. Optical micrographs of the SnO 2 sol drop deposited onto PI substrate: a) and c) without surface treatment of the substrate, drop to drop spaced of 100 µm and 20 µm, respectively; b) and d) with O 2 plasma treatment of the substrate, drop to drop spaced of 100 µm and 20 µm, respectively. Red circles underline deposited drop size.

Fig. 5 .

 5 Fig. 5. SEM micrograph of the sensor. The two grey parallel lines correspond to the gold electrodes and the grey dots correspond to SnO 2 accumulation.

Fig. 6 .

 6 Fig. 6. TEM micrographs of the cross section of the sensor, a) overview and b) SnO 2 layer at higher magnification. Carbon and platinum layers were deposited during the lamella preparation with FIB technique.

Fig. 7 .

 7 Fig. 7. Evolution of SnO 2 electrical conductance with CO and NO 2 injections at operating temperatures of a) 300°C and b) 200°C in dry air; c) at 300°C in wet air.

Fig. 8 .

 8 Fig. 8. Calibration curves for CO and NO 2 . The chemo-resistive responses (R air /R CO or R NO2 /R air ) are the mean value of the two tested sensors and error bars represent the standard deviation.

Table 1 .

 1 Limits of detection of CO and NO 2 of the SnO 2 inkjet printed sensors.

	Limit of detection	300°C	300°C	200°C	200°C
		dry air	wet air	dry air	wet air
	CO (ppm)	0.5	24	0.4	46
	NO 2 (ppb)	6	9	1	6

Acknowledgements

This work was supported by the LABEX MANUTECH-SISE (ANR-10-LABX-0075) of Université de Lyon, within the program "Investissements d'Avenir" (ANR-11-IDEX-0007) operated by the French National Research Agency (ANR).