Simon Pageaud

Nicolas Ponthus

Romain Gauchon

Catherine Pothier

Christophe Rigotti

Anne Eyraud-Loisel

Jean-Pierre Bertoglio

Alexis Bienvenüe

François Gueyffier

Philippe Vanhems

Nicolas Leboisne

Jean Iwaz

Stéphane Loisel

Pascal Roy
email: pascal.roy@chu-lyon.fr

Adapting French COVID-19 vaccination campaign duration to variant dissemination

Background The outbreak of SARS-CoV-2 virus has caused a major international health crisis with serious consequences in terms of public health and economy. In France, two lockdown periods were decided in 2020 to avoid the saturation of intensive care units (ICU) and an increase in mortality. The rapid dissemination of variant SARS-CoV-2 VOC 202012/01 has strongly influenced the course of the epidemic. Vaccines have been rapidly developed. Their efficacy against the severe forms of the disease has been established, and their efficacy against disease transmission is under evaluation. The aim of this paper is to compare the efficacy of several vaccination strategies in the presence of variants in controlling the COVID-19 epidemic through population immunity.

Methods An agent-based model was designed to simulate with different scenarios the evolution of COVID-19 pandemic in France over 2021 and 2022. The simulations were carried out ignoring the occurrence of variants then taking into account their diffusion over time. The expected effects of three Non-Pharmaceutical Interventions (Relaxed-NPI, Intensive-NPI, and Extended-NPI) to limit the epidemic extension were compared. The expected efficacy of vaccines were the values recently estimated in preventing severe forms of the disease (75% and 94%) for the current used vaccines in France (Pfizer-BioNTech and Moderna since January 1 11, 2021, and AstraZeneca since February 2, 2021). All vaccination campaigns reproduced an advanced agebased priority advised by the Haute Autorité de Santé. Putative reductions of virus transmission were fixed at 0, 50, 75 and 90%. The effects of four vaccination campaign durations (6-month, 12-month, 18-month and 24-month) were compared.

Results In the absence of vaccination, the presence of variants led to reject the Relaxed-NPI because of a high expected number of deaths (170 to 210 thousands) and the significant overload of ICUs from which 35 thousand patients would be deprived. In comparison with the situation without vaccination, the number of deaths was divided by 7 without ICU saturation with a 6-month vaccination campaign. A 12-month campaign would divide the number of death by 3 with Intensive-NPI and by 6 with Extended-NPI (the latter being necessary to avoid ICU saturation). With 18-month and 24-month vaccination campaigns without Extended-NPI, the number of deaths and ICU admissions would explode.

Conclusion Among the four compared strategies the 6-month vaccination campaign seems to be the best response to changes in the dynamics of the epidemic due to the variants. The race against the COVID-19 epidemic is a race of vaccination strategy. Any further vaccination delay would increase the need of strengthened measures such as Extended-NPI to limit the number of deaths and avoid ICU saturation.

INTRODUCTION

The outbreak of SARS-CoV-2 virus is causing major national, European, and international health crises, with serious consequences in terms of public health and economy. The COVID-19 epidemic started in China in late autumn 2019 and, on December 31 [START_REF] Zhu | A novel coronavirus from patients with pneumonia in China[END_REF], the WHO China Office [START_REF]SARS-CoV-2 variants[END_REF] received reports on severe pneumonia cases in the city of Wuhan (Hubei Province). On February 14, 2021, more than 100 million cases and 2.4 million officially recognized deaths were recorded worldwide of which 34 million cases and 760 thousand deaths occurred in Europe and 3.4 million cases and 81.6 thousand deaths occurred in France. Europe has been one of the most severely hit continents.

In France, during 2020, two lockdown periods were decided to avoid the saturation of intensive care units (ICU) that would have led to an increase in mortality. The first period (March 17 to May 11, 2020) was decided to limit the magnitude of the first wave. The corresponding maximum 7-day moving average of hospital deaths was 532 on April 8, 2020. Less than 30 daily hospital deaths were observed between June 6 and September 16, defining a plateau. Daily values of 222 and 276 were observed at the beginning (October 30, 2020) and the end (December 15, 2020) of the second lockdown, with a maximum 7-day moving average of hospital deaths of 419 on November 19, 2020. Nevertheless, the resulting epidemic situation was less favorable after the second lockdown; there was a residual plateau with 240 to 340 daily deaths.

On December 14, 2020, the variant SARS-CoV-2 VOC 202012/01 (B.1.1.7) was identified in United Kingdom [START_REF] Volz | Transmission of SARS-CoV-2 Lineage B.1.1.7 in England: Insights from linking epidemiological and genetic data[END_REF]. The first case with this variant appeared in France on December 2020. It was deemed responsible of 3.3% of all new cases on January 8, 2021 [START_REF]COVID-19 : point épidémiologique du 28 janvier 2021[END_REF], 14% on January 27 [START_REF]COVID-19 : point épidémiologique du 11 février 2021[END_REF], and 35% on February 16. By the end of February 2021, this percentage level had reached 50%. Preliminary studies suggest that SARS-CoV-2 VOC 202012/01 increased transmissibility (infectiousness) but not disease severity.

Vaccines have been rapidly developed. Their efficacy against the severe forms of the disease was estimated between 75% [START_REF] Voysey | Safety and efficacy of the ChAdOx1 nCoV-19 vaccine (AZD1222) against SARS-CoV-2: an interim analysis of four randomised controlled trials in Brazil, South Africa, and the UK[END_REF] and 94% [START_REF] Lindsey R Baden | Efficacy and safety of the mRNA-1273 SARS-CoV-2 vaccine[END_REF][START_REF] Fernando P Polack | Safety and efficacy of the BNT162b2 mRNA Covid-19 vaccine[END_REF]. The efficacy of these vaccines on the disease transmission is still under evaluation. The objective of the current massive vaccination strategy in France is to achieve the extinction of the epidemic through herd immunity resulting from an extended vaccination coverage.

An important public health issue is to optimize the population vaccination schedule and ensure that necessary protective measures are adapted to the epidemic situation and its evolution. The dynamics of the epidemic is continuously assessed through changes in hospital admissions, intensive care units (ICU) occupancy rate, and cumulative number of deaths associated with COVID-19. In France, protective measures have been constantly adapted to the dynamics of the epidemic. These measures are monitored using the value of the effective reproductive number and the weekly incidence rate.

The effectiveness of a vaccination campaign depends strongly on the availability of vaccine vials and on their distribution strategy. A campaign should also adapt to the evolution of the epidemic, especially in the presence of variants. This study aims to analyze the expected dynamics of the COVID-19 epidemic after applying protective measures and taking into account the increasing proportion of more infectious variants (e.g. SARS-CoV-2 VOC 202012/01) and several vaccination strategies.

MATERIALS AND METHODS

Compartmental model

The model is a compartmental model, stratified by age groups and adapted from the model described in [START_REF] Di | Impact of lockdown on COVID-19 epidemic in Île-de-France and possible exit strategies[END_REF]. A graphical view of the model can be found in [START_REF] Gauchon | Lessons learnt from the use of compartmental models over the COVID-19 induced lockdown in France[END_REF]. The model compartments are the following:

• S (Susceptible): individuals who are not infected by the virus;

• The incubation period included two compartments:

when infected, an individual goes from S to compartment E (Exposed) that contains infected individuals who did not develop symptoms yet and are not contagious. The mean length of stay in E is (t i -t p ) where t i is the incubation period and t p the duration of the prodromal state (see next item);

when an individual from E starts to be contagious, it is transferred to compartment I p (Prodromal phase which is the short phase that follows contamination without pathological symptoms but with possible non-specific prodromes) that contains individuals still in the incubation period already contagious but without apparent symptoms. After a stay of average duration t p , the individuals enter one of the four compartments A, I ps , I ms , and I ss defined hereafter.

• A (Asymptomatic): individuals who completed the incubation period, became infectious, but have no disease symptoms (with probability p A ). The mean length of stay in A is t s ;

• The symptomatic infectious period includes three compartments for individuals developing symptoms (with probability 1 -p A ):

-I ps (Paucisymptomatic disease): individuals with weak disease symptoms;

-I ms (Medium symptoms): individuals with disease symptoms (e.g. fever or cough) who will not require hospitalization. The average durations in states I ps and I ms are the same as the duration in state A;

-I ss (Severe symptoms): individuals severely infected, who will require hospitalization. They stay in I ss before being hospitalized. The mean length of stay in I ss is t bh .

Once symptomatic, the probabilities to enter I ps , I ms , or I ss are respectively p Ips , p Ims , and p H (probabilities summing to one).

• The hospitalization period includes two compartments:

-ICU : hospitalized individuals who during part of their stay will require an ICU bed;

-H (Hospitalized): hospitalized individuals who will not require ICU care during their stay.

Individuals coming from compartment I ss are admitted to ICU with probability p ICU and to H with probability 1 -p ICU . The mean length of stay in ICU (resp. H) is t icu (resp. t h ).

• In this model, there are two absorbing compartments:

-D (Deceased at hospital): probability to enter D after ICU (resp. H) is p D|ICU (resp. p D|H );

-R (Removed): individuals who recovered (from any compartment).

The dynamics can be expressed according to the following system of ordinary differential equations where α denotes the incidence hazard rate:

                                                                                                 dS dt = -αS, dE dt = αS - 1 t i -t p E, dI p dt = 1 t i -t p E - 1 t p I p , dA dt = p A t p I p - 1 t s A, dI ps dt = (1 -p A )p Ips t p I p - 1 t s I ps , dI ms dt = (1 -p A )p Ims t p I p - 1 t s I ms , dI ss dt = (1 -p A )p H t p I p - 1 t bh I ss , dH dt = (1 -p ICU ) t bh I ss - 1 t h H, dICU dt = p ICU t bh I ss - 1 t icu ICU, dR dt = 1 t s (I ps + I ms + A)+ 1 -p D|H t h H + 1 -p D|ICU t icu ICU, dD dt = p D|H t h H + p D|ICU t icu ICU. (1)
Here, each compartment was split into 17 age groups, as in [START_REF] Gauchon | Lessons learnt from the use of compartmental models over the COVID-19 induced lockdown in France[END_REF], whereas only three age groups were considered in [START_REF] Di | Impact of lockdown on COVID-19 epidemic in Île-de-France and possible exit strategies[END_REF]. The incidence hazard rate α could vary over time and was different between each age groups (as defined in the next section). Parameter values are given in tables 1 and 2.

Modeling the disease contagiousness

For a given age group i, at a given time t, the incidence hazard rate α i t was calculated as follows:

α i t = β 1 β 2,t 17 j=1   C i,j Z∈{E,Ip,A,Ips,Ims,Iss} i Z Z j t N j tot   . (2)
where:

• C is the contact matrix where an element C i,j is the average number of individuals from age group j encountered per day by a single individual from age group i;

• Z j t is the number of individuals of age group j in compartment Z at time t; • N j tot is the total number of individuals from age group j (who had contracted the disease or not); • i Z is the relative infectiousness of individuals in compartment Z. For the historical strain, it is set to 1 in case of medium or severe symptoms; it ranges from 0 to 1 for the other compartments;

• β 1 is the daily incidence hazard rate. In can be interpreted as a hazard rate in case of one daily contact with a medium-or severe-symptom case infected by the historical strain;

• The relative hazard β 2,t estimates the effect of Non-Pharmaceutical Interventions (NPI), also named barrier measures, such as mask wearing, social distancing, lockdown, and auto-isolation after infectious contacts. The underlying hypothesis is that NPI effects are similar in all compartments and age groups, with a multiplicative effect on the hazard (without multiplicative interaction). β 2,t can be interpreted as a coefficient reducing the mean number of daily contacts and/or the mean duration of contact, and/or their infectiousness. The introduction of barrier measures leads to 0 ≤ β 2,t ≤ 1, by ignoring here an improbable exacerbation of risk behaviors (β 2,t > 1) once the epidemic is declared.

Reproductive numbers

At an epidemic start, almost all individuals are susceptible and no protective measure are already taken. The effective reproductive number R e (t) is then at its maximum and corresponds to the basic reproductive number R 0 which characterizes the virus infectiousness and the initial environmental conditions of viral spread (i.e. before the initiation of the NPI). The reproductive numbers R 0 and R e (t) give an indication of the epidemic dynamics: the highest they are, the faster the epidemic spreads. The effective number R e (t) decreases with time, following the decreased proportion of susceptible individuals combined with the introduction of NPI and is defined by

R e (t) = R 0 β 2,t S(t) N tot , (3)
where N tot is the total size of the population.

To compute R 0 as defined by Diekmann et al. [START_REF] Odo Diekmann | On the definition and the computation of the basic reproduction ratio R 0 in models for infectious diseases in heterogeneous populations[END_REF], it is necessary to introduce the concept of next generation matrix. In the following, the approach and notational conventions of Van den Driessche and Watmough [START_REF] Van Den Driessche | Reproduction numbers and sub-threshold endemic equilibria for compartmental models of disease transmission[END_REF] are used. It is assumed that there are m different compartments of infected individuals and n different age groups. Thus, there is a total of nm groups of infected (but not necessarily infectious) individuals. Let x(t) = (x 1 (t), . . . , x nm (t)), with quantity x i (t), 1 ≤ i ≤ nm, denote the number of individuals in the ith group of infected individuals at time t. Let F i (x(t)) be the rate of new infected individuals in the ith group of infected individuals at time t.

Let V i (x(t)) = V - i (x(t)) -V + i (x(t)) where V - i (x(t)
) is the rate of departures from group i, and

V + i (x(t)
) is the rate of arrivals in group i by all other means than new infections. Let

F = ∂+Fi(x 0 ) ∂xj 1≤i,j≤nm and V = ∂+Vi(x 0 ) ∂xj 1≤i,j≤nm
where x 0 is a disease free equilibrium and ∂ + denotes right partial derivatives. The next generation matrix is then defined as

G = F V -1 . (4)
As shown in [START_REF] Van Den Driessche | Reproduction numbers and sub-threshold endemic equilibria for compartmental models of disease transmission[END_REF], V -1 i,j can be seen as the average time an individual introduced in the jth group of infected individuals will spend in the ith group of infected individuals during the remaining lifetime. Thus, intuitively, an entry (i, j) of the next generation matrix F V -1 is the expected number of new cases in group i due to individuals introduced in group j. R 0 is then defined as the spectral radius of G. It is computed before initiation of the NPI (i.e., with β 2,t = 1), then, as shown in [START_REF] Gauchon | Lessons learnt from the use of compartmental models over the COVID-19 induced lockdown in France[END_REF], the expression of R 0 is of the form R 0 = kβ 1 with k depending on the other parameters. In the following, using Equation ( 3), this provides the relationship to calculate β 2,t from R e after β 1 has been estimated.

Agent-based simulation of the pandemic

The model was implemented as a stochastic Agent-Based Model (ABM). A discrete time step of one day was used for the simulations assuming that the derivatives in the equation system (1) were constant functions over this interval.

Parameters

The agent states corresponded to the compartments, and the probability to exit from a state was obtained from the mean length of stay in the associated compartment. When more than one destination state were possible, the transition is selected at random using the corresponding probabilities (i.e., p A , p Ips , p Ims , p H , p ICU , p D|H and p D|ICU ).

The age distribution of the population came from INSEE (French Institut National de la Statistique et des Etudes Economiques) studies; it considered 16 five-year age groups from 0 to 79 years plus an extra agegroup with individuals aged 80 or older. The contact matrix was available for France on the basis of these age ranges [START_REF] Prem | Projecting social contact matrices in 152 countries using contact surveys and demographic data[END_REF]. The number of daily contacts per individual was set at the beginning of the simulation and was drawn from a Poisson distribution with a parameter set to the corresponding element of the contact matrix. At each time step, the daily contacts between infectious and susceptible agents were drawn uniformly at random. For each such contact, the probability of transmission was determined by β 1 , β 2,t and the value of the relative infectiousness coefficient associated with the compartment of the infectious agent.

The proportion of asymptomatic forms was set at 20% and 35% (Table 2). The proportions of various symptomatic forms, paucysymptomatic, medium symptom and severe symptom were specified as depending on age (Table 1). The relative infectiousness of prodromal, asymptomatic, and paucysymptomatic forms of the historical strain was set to i Ip = i A = i Ips = 0.55; medium and severe symptom forms were taken as references with i Ims = i Iss = 1 (Table 2).

Parameter β 1 was estimated using the same method as [START_REF] Gauchon | Lessons learnt from the use of compartmental models over the COVID-19 induced lockdown in France[END_REF], by fitting hospitalization data before and during the first lockdown in France. The maximum likelihood estimation led to β 1 = 0.05052 when the asymptomatic rate was 20% and to β 1 = 0.05331 when the asymptomatic rate was 35%. On January 1, 2021, the beginning of the simulations, the estimated proportion of removed individuals was 13% [START_REF]Proportion de la population ayant été infectée par SARS-CoV-2[END_REF] and the proportion of variant considered negligible [START_REF]COVID-19 : point épidémiologique du 7 janvier 2021[END_REF].

Scenarios

Three different scenarios were considered: Relaxed-NPI, Intensive-NPI, and Extended-NPI. These scenarios took into account a mix of several protective measures taken to reach a specific effective reproduction number R e . The three typical groups of protective measures were designed from the French restrictions observed during the past year: i) Intensive protective measures consisted in a combination of strict protective measures, such as overnight curfew, social distancing, mask-wearing, increased distancing at work, closing of bars, restaurants and leisure facilities, etc.; ii) Extended protective measures consisted in very strict sanitary measures enhanced compared to the intensive ones, such as weekend or full lockdown, as well as closing of schools, universities, and shopping facilities; iii) Moderate protective measures consisted in relaxing the intensive measures (such as restaurant and leisure facility reopening) but maintaining some protective measures such as mask wearing and social distancing. The Intensive-NPI included constant intensive protective measures. The Extended-NPI included constant extended protective measures. The Relaxed-NPI associated systematic alternation of 45-days periods with Intensive protective measures and Moderate protective measures.

These three scenarios that corresponded to different protective procedures were specified through the value of the relative hazard β 2,t , under the hypothesis of a similar multiplicative effect of protective measures whatever the viral strain, without interaction. On January 1, 2021, the estimated effective reproductive number was R e = 1 [START_REF]COVID-19 : point épidémiologique du 7 janvier 2021[END_REF]. For R e = 1, the calculated relative hazard at this date was β 2,t = 0.2501. This value of parameter β 2,t was the value associated with the so-called Intensive protective measures, in order to match the effective situation of the sanitary protective measures adopted in France on January 1, 2021. The values of parameter β 2,t for the other protective measures were calculated with hypothetical measures leading to different values of R e : R e = 0.8 (Extended protective measures, β 2,t = 0.2001) and R e = 1.5 (Moderate protective measures, β 2,t = 0.3752) on January 1, 2021.

The expected effects of the three protective scenarios were compared. Due to the stochasticity of the agentbased simulation, 50 simulations were performed with each set of parameters. Each run involved 322,488 agents (i.e., 0.5% of the population of metropolitan regions of France, of about 64.5 million inhabitants). To obtain estimations for France, simulation results were multiplied by 200. The mean values and standard errors of the mean across the runs were plotted.

New variant influence

The so-called 'variant strain' herein corresponds to several variants considered more contagious than the historical strain. The relative infectiousness of the variant strains vs. the historical strain were supposed independent of the compartment. This corresponded to multiplying each i z by a parameter β 3 ; the latter parameter was the same whatever the compartment and was characteristic of the relative infectiousness of the variant strain. The value of β 3 fitting with a prediction of a proportion of variants close to the available data on January 8, 2021 (3.3%), January 27, 2021 (14%), and February 18, 2021 (36%) [START_REF]COVID-19 : point épidémiologique du 28 janvier 2021[END_REF][START_REF]COVID-19 : point épidémiologique du 11 février 2021[END_REF][START_REF]COVID-19 : point épidémiologique du 18 février 2021[END_REF] was found to be β 3 = 1.7, with at the beginning of the simulation (January 1, 2021), an intercept of the new variant rounded up to 2% of the population infected with the historical strain.

Vaccine effects

The expected effects of vaccination were set at their recently estimated levels of efficacy in preventing severe forms of the disease; that is, 94% [START_REF] Lindsey R Baden | Efficacy and safety of the mRNA-1273 SARS-CoV-2 vaccine[END_REF][START_REF] Fernando P Polack | Safety and efficacy of the BNT162b2 mRNA Covid-19 vaccine[END_REF] and 75% [START_REF] Voysey | Safety and efficacy of the ChAdOx1 nCoV-19 vaccine (AZD1222) against SARS-CoV-2: an interim analysis of four randomised controlled trials in Brazil, South Africa, and the UK[END_REF], depending on the vaccine. When the vaccine was effective, a vaccinated individual could develop only an asymptomatic form of the disease. Putative reductions of virus transmission were fixed at 0, 50, 75 and 90%, the transmission probabilities being considered independent of the vaccine. An age-based priority for being vaccinated was used. Three vaccination priorities were considered and all individuals of a given age group had to receive their first injection before giving the first injections to the next priority age group. The age priorities considered were: 75 and older, 65 and older, then the rest of the population. In the rest of the population, two vaccination options were compared: vaccinate at all ages vs. vaccinate individuals aged 20 years and over (vaccination is under investigation in the youngest population). The effects of four vaccination campaigns over 6, 12, 18, 24 months were compared. At the end of the vaccination campaign, all individuals eligible to vaccination had to have received their second injection. In all simulations, vaccines with 94% efficacy were attributed to individuals aged 75 or older, whereas vaccines with 75% efficacy were attributed to the rest of the population. To be considered as vaccinated, an individual had to receive two injections of any vaccine with a 21-day interval. In all simulations, vaccine supplies were considered sufficient to follow the study's vaccination schedules. Furthermore, all vaccines were considered as efficient against the historical strain as against new variants similar to the B.1.1.7 lineage [START_REF] Wang | Antibody Resistance of SARS-CoV-2 Variants B.1.351 and B.1[END_REF].

Outcome criteria

The main outcome was the cumulative number of individuals removed, the cumulative number of deaths at hospital, the daily prevalence in ICU beds and its saturation indicator. The saturation of ICU beds was calculated as the cumulative number of new cases requiring ICU when all beds were already occupied. Three thresholds were considered: 5,000 (the current number of ICU beds), 8,000 (adding beds in intermediate care units), and 12,000 (adding also more beds in conventional care units). The model specified that, on average, people in the intensive care compartment stayed the same time in an intensive care bed or a conventional hospital bed 1 . The cumulative number of individuals removed due to disease was calculated as of December 31, 2022.

RESULTS

Dynamics of the epidemic with the historical strain without vaccination

Table 3 presents the results obtained by restricting the simulations to the historical strain and for 20% and 35% prevalence of asymptomatic cases, without vaccination. After two years of Relaxed-NPI, the simulated numbers of removed individuals were 20.4 and 20.5 millions and the numbers of deaths were 62.6 and 51.2 thousand respectively. The simulated numbers of ICU admissions that exceeded the 5,000-bed threshold were 2,860 and 1,890, respectively. In this scenario, the cumulative distribution of deaths plateaued after 500 days (Figure 1a). The waves observed for the Relaxed-NPI scenarios show the slightly delayed effects of the systematic alternation of periods of 45 days with intensive protection measures and moderate protection measures (Figure 1b). After an initial decrease, three successive waves were observed at the end of the three first 45-days periods with Moderate protective measures. The second wave exceeded ICU capacities (Figure 1b), whereas the third wave showed a decrease before an extinction of the epidemic.

After two years of the Intensive-NPI, the simulated numbers of removed individuals was much lower (around 9.6 million) with 12.7 and 10.4 thousand deaths, respectively. The cumulative distribution of deaths plateaued after 250 days (Figure 1a). The ICU capacities were not exceeded. A rapid decrease was observed whatever the prevalence of asymptomatic cases (Figure 1b) leading to the extinction of the epidemic.

After two years of the Extended-NPI, the simulated numbers of removed individuals were around 8.9 million with 9.7 and 8.1 thousand deaths, respectively. The cumulative distribution of deaths plateaued after 100 days (Figure 1a). The ICU capacities were not exceeded. A rapid decrease was observed whatever the prevalence of asymptomatic cases (Figure 1b) leading to the extinction of the epidemic.

Dynamics of the epidemic with the introduction of new variant strains without vaccination

Table 4 presents the results obtained taking into account a competition between the historical strain and the more infectious new variants with 20% and 35% prevalence of asymptomatic cases, without vaccination. Applying the Relaxed-NPI, the simulated numbers of removed individuals were around 40.2 million, with 207.7 and 168.5 thousand deaths. The cumulative distribution of deaths plateaued after 250 days (Figure 2a). The simulated numbers of ICU admissions that exceeded the 5,000-bed threshold were 35.1 and 27.2 thousand respectively. The capacities of the ICU were also saturated beyond 8,000 or 12,000 beds. A huge single wave was observed, resulting from the fusion of the effects of the two first 45-days periods with moderate protective measures, leading to the extinction of the epidemic after this single wave (Figure 2b). Given these simulated results, the Relaxed-NPI was not further retained.

When the Intensive-NPI was applied, the simulated numbers of removed individuals were still important (around 32.5 million) with 137.6 and 111.3 thousand deaths, respectively. The cumulative distribution of deaths plateaued after 300 days (Figure 2a). The simulated saturation of ICU capacities over 5,000 beds concerned 24.3, and 18.9 thousand patients. The capacities of the ICU were saturated beyond 8,000 or 12,000 beds. A huge wave was observed, despite the keep up of the Intensive-NPI, leading to the extinction of the epidemic after a massive contamination of the population (Figure 2b).

When the Extended-NPI were applied, the simulated numbers of removed individuals were 21.8 million with 70.2 and 56.4 thousand deaths, respectively. The cumulative distribution of deaths plateaued after 550 days (Figure 2a). The simulated numbers of ICU admissions that exceeded the 5,000-bed threshold were 1.9 and 0.7 thousand. The capacities of the ICU were never saturated beyond 8,000 beds. Compared to that simulated under Intensive-NPI, the wave observed under Extended-NPI was delayed with a maximum at 300 days and an extinction of the epidemic after 550 days (Figure 2b).

Effects of a 6-month vaccination campaign

Table 5 shows the results obtained with a 6-month vaccination campaign.

Applying the Intensive-NPI with a 20% prevalence of asymptomatic forms, a 90% reduction in virus transmission halved the number of removed and reduced the number of deaths by 20-25% vs. no vaccine effect on virus transmission.

In comparison with vaccination of individuals aged 20 or older, the vaccination of all individuals reduced the number of deaths by 10%.

The ICU resources were never saturated. With a 6-month vaccination campaign of all individuals without transmission reduction (Table 5) vs. no vaccination and Intensive-NPI (Table 4), the number of deaths was divided by 9 with the Intensive-NPI and by 13 with the Extended-NPI.

Effects of a 12-month vaccination campaign

Table 6 shows the results obtained with a 12-month vaccination campaign.

Applying the Intensive-NPI, the number of deaths was three times higher than with a 6-month vaccination campaign under similar conditions.

In all simulations with Intensive-NPI, the ICU resources were saturated. The application of the extended-NPI made possible to avoid the ICU resources saturation.

With a 12-month vaccination campaign of all individuals without transmission reduction (Table 6) vs. no vaccination and Intensive-NPI (Table 4), the number of deaths was divided by 3 with the Intensive-NPI and by 10 with the Extended-NPI.

Effects of 18-month and 24-month vaccination campaigns

Tables 7 and8 show the results obtained with an 18-month and a 24-month vaccination campaign respectively.

With the Intensive-NPI, whatever the simulated conditions, the ICU resources were saturated. The level of intensive care saturation increased markedly between the 12 and the 18 month vaccination campaigns, and continued to increase but to a lesser extent between the 18 and the 24 month vaccination campaigns.

In all vaccination campaigns exceeding 6 months, the Extended-NPI was necessary to avoid ICU resource saturation.

Versus no vaccination and Intensive-NPI (Table 4), an 18-month vaccination campaign of all individuals without transmission reduction divided the number of deaths by 2 with the Intensive-NPI and by 6 with the Extended-NPI (Table 7), whereas a 24-month vaccination campaign of all individuals without transmission reduction divided the number of deaths by 2 with the Intensive-NPI and by 4 with the Extended-NPI(Table 8).

Determinants of the vaccination campaign efficacy

Figure 3 shows the results of the simulations of the epidemic from January 1, 2021, depending on the duration of the vaccination campaign, the reduction in the risk of viral transmission, the prevalence of asymptomatic forms and the Intensive/Extensive type of NPI. The alternative between the Intensive-NPI / Extended-NPI strongly impacted the mortality and the ICU saturation. For a 2-year vaccination campaign restricted to individuals aged 20 or older, without reduction in viral transmission, with a prevalence of asymptomatic forms of 20%, the Intensive-NPI led to 93.3 thousand deaths while 21.1 thousand patients exceeded the ICU resources, whereas the Extended-NPI led to 37.1 thousand deaths without ICU saturation. Limiting the duration of the vaccination campaign to 6 months would lead, under the same other conditions, to 24.7 thousand deaths without ICU saturation under the Intensive-NPI, and to 11.2 thousand deaths without ICU saturation under the Extended-NPI. In case the 6-month strategy was not realistic, it was mandatory to limit the duration of the vaccination campaign and apply the Extended-NPI to limit the number of deaths and to avoid the ICU saturation.

DISCUSSION

Historical strain without vaccination

Table 3 presents the cumulative number of individuals removed due to disease, the cumulative number of diseaserelated deaths and ICU overload after two years with the historical strain and without vaccination, under Relaxed-NPI, Intensive-NPI and Extended-NPI strategies. It shows that adopting a Relaxed-NPI would lead to an important increase in the number of deaths, whatever the probability of asymptomatic form (either 20% or 35%). However, without a variant and even without a vaccine, it would not have been mandatory to switch to Extended-NPI; the Intensive-NPI would have been acceptable while waiting for a sufficient herd immunity, although the number of deaths generated would have been naturally higher than the one expected with Extended-NPI. Figure 1 shows also that Relaxed-NPI would have been sufficient to induce a virus extinction within approximately 500 days, and limit the economic cost of the sanitary situation, but at a certain cost in terms of deaths and ICU beds, and provided these NPI were continuously maintained.

Introduction of new variant strains without vaccination

Figure 2 shows that the Relaxed-NPI are still sufficient to lead to an extinction of both historical and variant strain virus within approximately 250 days, but at very high costs in terms of deaths (170 to 210 thousand) and of ICU beds (35 thousand patients would be deprived). With the introduction of variants, Relaxed-NPI would not be a reasonable option.

Strikingly, as soon as the variant strain is included in the computations, the Intensive-NPI works no more in limiting sufficiently the number of deaths and maintaining an acceptable pressure on ICU beds. The emergence of more infectious variant strains of the virus makes ineffective both Relaxed-NPI and Intensive-NPI. Table 4 shows that it is even necessary to strengthen Intensive-NPI and to impose Extended-NPI to avoid an explosion of the number of deaths, and a total overload of the ICU resources. The level of Herd immunity then reached does not make possible to release the NPI because it does not protect against a new epidemic wave once these protective measures have been withdrawn.

Without vaccines, 550 days would be necessary to halt the epidemic, and Extended-NPI would be inevitable. Without vaccines and with new variants, using simply current preventive measures would be very costly in terms of deaths, ICU beds, and paralysis of the health care system. It is then important to accelerate the vaccination campaigns. Increasing herd immunity through vaccination is necessary to avoid the resumption of the epidemic after NPI removal.

Adapting the vaccine strategy to new variant strains

In the absence of vaccination, the emergence of a new variant strain changes definitely the dynamics of the epidemic. The vaccination strategy has to be adapted to a new epidemic dynamics. The present study aims to comparing the impacts of 4 vaccination strategies to face a new epidemic dynamics (6, 12, 18 and 24-month campaigns).

The main conclusion of modelling the vaccination strategy effect is that the best ways to deal with a historical and a new strain, without overloading ICUs is to vaccinate the entire population or all persons aged over 20 within 6 months (The minimum age for vaccination will depend on the results of ongoing pediatric clinical trials). This is shown in Table 5 where a maximal of 18 thousand deaths and no ICU overload are observed with applying the Intensive-NPI. Currently, vaccinating the whole French population within 6 months does not seem to be an achievable goal; the vaccination campaign has begun on January 11, 2021, and on March 13, 2021, 2,220,608 people have been already fully vaccinated (3, 3% of the population). Thus, despite an ongoing significant acceleration, the goal would still not be achievable. Tables 6, 7 and 8 show the potential effects of 12, 18 and 24 month vaccination campaigns. Despite intensive NPI, vaccination campaigns lasting 18 months or 24 months would generate up to 81 and 93 thousand deaths respectively, and significant pressure on ICU beds in both scenarios.

Taking into account the currently unlikely vaccination strategy within 6 months, it appears that vaccinating the whole population within 1 year would bring real benefits in terms of deaths and ICU saturation. As shown in Table 6, the 12-month vaccination strategy provides reasonable results with Extended-NPI. By contrast, with Intensive-NPI, the number of deaths may reach 61 thousand in the worst scenario and exceed 20 to 30 thousand in many others whereas the pressure on ICU beds would exceed the maximum capacities by several thousand. With the Intensive-NPI scenario, there are real risks of paralyzing the health-care system. The only solution to overcome the possible failure of a 6-month vaccination strategy seems to be the adoption of Extended-NPI and the vaccination campaign within a maximum of 12 months.

The need for Extended-NPI seems to be inevitable, at least during a transitory period. Relaxing Intensive-NPI will be definitely a crucial issue for economic and social reasons. It will highly depend on the number of vaccinated persons : i.e. as long as the number of vaccinated persons is insufficient to protect the population and to obtain herd immunity, it will be difficult to release the NPI without worsening additional deaths and ICU pressures. One notice that the situation would be far better than a vaccination strategy over more than 12 months. An important additional argument in favour of limiting the duration of the vaccination campaign is the difficulty to maintain the acceptation of the appropriate NPI at population level over a large period of time.

It is important to mention that in all extended-NPI scenarios, the duration of the vaccination campaign impacts the number of deaths: 9-10 thousand with the 6-month strategy, 11-16 thousand with the 12 months strategy, 14-27 thousand with the 18 month strategy and 20-37 thousand with the 24-month one. Although ICU resources are not overwhelmed applying the extended-NPI, the duration of the vaccination campaign strongly influences their use (Figure 3). The hospitalization of COVID-19 patients in intensive care significantly disrupts medical activity, forcing the deprogramming of much of usual medical and surgical activity. There is a real need to accelerate the vaccination of "the whole population in one year" to maintain care capacities, treat patients in intensive care and avoid excessive deaths.

The present study makes the strong hypothesis of a 100% acceptation level of vaccination. Virus extinction will happen if the population is widely vaccinated. This means that more efforts should be made in terms of vaccination strategy, and wise communication regarding the vaccines. It would be interesting to tackle the issue of vaccine acceptation and measure the impacts of partial acceptation levels on the extinction of the pandemic.The adverse effects of various vaccines do influence vaccination strategies and the vaccine acceptation level. The present modeling deliberately ignored the side effects of COVID-19 vaccines, these being considered negligible compared to the expected benefits of establishing herd vaccine immunity. Another choice the modelling was ignoring the duration of vaccine protection. This might be the object of future investigations

Conclusion

The main conclusion of this study is that the key point of this race against the COVID-19 historical strain and its variant strains is an issue of vaccination strategy. It is mandatory to vaccinate most of the population within a year, and preferably within 6 months. Should a 6 month vaccination campaign not be feasible, then reinforced NPI should be considered.

These conclusions are in agreement with the recommendations of Bosetti et al. [START_REF] Bosetti | A race between SARS-CoV-2 variants and vaccination: The case of the B.1.1.7 variant in France[END_REF] that reinforced NPI applied to the most contagious variant are necessary to achieve epidemic control similar to that obtained by applying less stringent measures to the historical virus strain. These authors also insist that keeping transmission rates of the virus low helps prevent a worsening of the epidemic requiring enhanced NPI.

One advantage of the modeling presented in this paper is that it takes into account the particularity of the French vaccine strategy: vaccinating the most advanced age groups first. If it appears that it will not be possible to vaccinate all the population within one year, alternative solutions will have to be examined. Regarding accelerating vaccination strategy would it be more impactful to increase the delay between the two doses or increase the pace of vaccination through increasing the numbers of vaccination centers, vaccination human resources or available vaccine supplies? Answering this question is not the aim of the present report. Nevertheless, it seems important to mention that solutions such as increasing the delay between the two doses of vaccines could have a direct impact on the success of the vaccination, on its efficiency.

Another advantage of rapid vaccination is to help reduce the worldwide multiplication of the virus, and thus to limit the possible emergence of vaccine-resistant variants, and/or more virulent ones. The present paper is only calibrated on the historical strain and on the variants with contagiousness similar to the SARS-CoV-

  Prevalence of COVID-19 patients in needs of Intense Care Units. Horizontal black lines represent an ICU capacity of 5,000 and 8,000 beds.

Figure 1 :

 1 Figure 1: Projection of the COVID-19 pandemic in France for the historical strain, without vaccine, for Extended-NPI, Intensive-NPI and Relaxed-NPI, and for 20% and 35% asymptomatic rates.

  Cumulative number of deaths.

  (b) Prevalence of COVID-19 patients in needs of Intense Care Units. Horizontal black lines represent an ICU capacity of 5,000, 8,000 and 12,000 beds.

Figure 2 : 3 :

 23 Figure 2: Projection of the COVID-19 pandemic in France with the historical strain and new variants, without vaccine, for Extended-NPI, Intensive-NPI and Relaxed-NPI, and for 20% and 35% asymptomatic rates.

Table 1 :

 1 Simulation parameters, sources see[START_REF] Gauchon | Lessons learnt from the use of compartmental models over the COVID-19 induced lockdown in France[END_REF].

	age	0-4	5-9	10-14 15-19 20-24 25-29 30-34 35-39 40-44 45-49
	p Ips	0.25	0.25 0.249 0.249 0.247 0.247 0.242 0.242 0.238 0.238
	p Ims	0.749 0.749 0.748 0.748 0.741 0.741 0.726 0.726 0.713 0.713
	p Iss	0.001 0.001 0.003 0.003 0.012 0.012 0.032 0.032 0.049 0.049
	t h	8.1	8.1	8.1	8.1	8.1	8.1	7.4	7.4	8	8
	t icu	19.8	19.8	19.8	19.8	19.8	19.8	20.7	20.7	20.5	20.5
	p ICU	0.285 0.285 0.285 0.285 0.285 0.285 0.323 0.323 0.343 0.343
	p D|H	0	0	0	0	0	0	0.023 0.023 0.058 0.058
	p D|ICU	0	0	0	0	0	0	0.2	0.2	0.31	0.31
		age	50-54 55-59 60-64 65-69 70-74 75-79 80+
		p Ips	0.225 0.225 0.209 0.209 0.189 0.189 0.182
		p Ims	0.674 0.674 0.626 0.626 0.568 0.568 0.545
		p Iss	0.102 0.102 0.166 0.166 0.243 0.243 0.273
		t h	9.4	9.4	11.4	11.4	15.9	15.9	26.6
		t icu	24.5	24.5	17.6	17.6	23.4	23.4	19.2
		p ICU	0.408 0.408	0.36	0.36	0.282 0.282 0.05
		p D|H	0.081 0.081 0.147 0.147 0.227 0.227 0.434
		p D|ICU	0.25	0.25	0.47	0.47	0.54	0.54	0.47

Table 2 :

 2 Simulation parameters.

	Parameter Sources historical strain	Historical strain	New variants
	t i	see [10]	5.1	5.1
	t p	see [10]	1.5	1.5
	p A	[9, 20]	0.20 or 0.35	0.20 or 0.35
	i A	[21]	0.55	1.7 × 0.55 = 0.935
	i Ip	[21]	0.55	1.7 × 0.55 = 0.935
	i Ips	[21]	0.55	1.7 × 0.55 = 0.935
	i Ims	[21]	1	1.7
	i Iss	[21]	1	1.7
	β 1	calibration	0.05052 or 0.05331 0.05052 or 0.05331
	t s	see [10]	7	7
	t bh	see [10]	3	3

Table 3 :

 3 historical strain without vaccination: two-years cumulative numbers of individuals removed and hospital deaths, ICU overload exceeding beds capacity under different restrictions.

	NPI	Removed Deceased	ICU overload if max available beds is
				5,000	8,000	12,000
		Probability to develop asymptomatic form: 20%
	Extended-NPI	8,936,356	9,748	0	0	0
	Intensive-NPI	9,666,180	12,676	0	0	0
	Relaxed-NPI	20,427,252	62,612	2,860	0	0
		Probability to develop asymptomatic form: 35%
	Extended-NPI	8,944,300	8,128	0	0	0
	Intensive-NPI	9,634,884	10,356	0	0	0
	Relaxed-NPI	20,546,988	51,224	1,890	0	0

Table 4 :

 4 historical strain and variant development without vaccination: two-years cumulative numbers of individuals removed and hospital deaths, ICU overload exceeding beds capacity under different restrictions.

	NPI	Removed Deceased	ICU overload if max available beds is
				5,000	8,000	12,000
		Probability to develop asymptomatic form: 20%
	Extended-NPI 21,788,476	70,176	1,914	0	0
	Intensive-NPI 32,460,452	137,556	24,326	21,326	17,326
	Relaxed-NPI	40,223,456	207,732	35,060	32,060	28,060
		Probability to develop asymptomatic form: 35%
	Extended-NPI 21,818,740	56,416	672	0	0
	Intensive-NPI 32,482,356	111,332	18,940	15,940	11,940
	Relaxed-NPI	40,250,780	168,536	27,192	24,192	20,192

(a) Cumulative number of deaths.

Personal communication from the EHESP.

2VOC 202012/01 (B. 1.1.7), a virulence similar to the historic strain, and without specific vaccine-resistance. The available data don't allow to examine the effect of the emergence of more strains. Nevertheless, as soon as enough data are available to estimate a realistic β 3 for a new variant, the present model is a useful tool to compare the effects of various vaccine strategies on new variant strains. Studying other variants requires measuring the efficiency of the vaccines on these variants. It is possible to integrate a vaccine-resistance parameter specific to the variant in the model (variants such as B.1.351 [START_REF] Wang | Antibody Resistance of SARS-CoV-2 Variants B.1.351 and B.1[END_REF]), as soon as studies on the subject could provide serious estimations of such parameter. An increase of the virulence of the virus with a variant strain could be treated by modifying the transition probabilities in the model (see [START_REF] Challen | Increased hazard of mortality in cases compatible with SARS-CoV-2 variant of concern 202012/1-a matched cohort study[END_REF]).

Funding

This work was supported by Région Auvergne-Rhône-Alpes, Initiative d'excellence (IDEXLYON) -Université de Lyon, Université Claude Bernard Lyon 1 and École Centrale de Lyon and by the PHRCi-2020 COVID-3S supported by GIRCI Auvergne Rhône-Alpes.

Competing Interests

The authors have declared no competing interest.

Data Availability

All data used during this study are available on data.gouv.fr or in the published articles listed in the references, except APHP related data obtained from a personal communication of the EHESP. All data generated during this study are available from the corresponding author on reasonable request.


Erwan Bocher
email: erwan.bocher@ec-nantes.fr

Thomas Leduc
email: thomas.leduc@cerma.archi.fr

Guillaume Moreau
email: guillaume.moreau@ec-nantes.fr

Fernando González Cortés

GDMS: An abstraction layer to enhance Spatial Data Infrastructures usability

The practical exploitation of SDI (Spatial Data Infrastructures) raises number of issues as far as it grows. Among them is the heterogeneity of data sources and thus the difficulty for GIS users not to depend on the data source format and of course to learn different systems. This a major flaw with respect to reuse and data sharing. The purpose of our work is to propose a new semantic layer derived from the SQL language that is independent of the underlying data source. This layer, called GDMS (Generic Data source Management System) can first be seen as an abstraction layer between data sources and the SDI tools. We will also show how this layer extends both SQL and spatial semantics and improves the exploitation of the SDI, by providing feedback both in terms of work and data reuse. A simple example mixing heterogeneous data sources will be presented.

INTRODUCTION

Recently, attention has been focused on the development of Spatial Data Infrastructures at different levels (local, national) and for various policy areas and priority of actions such as social, economic, and environmental issues [START_REF] Executive | Order 12906: Coordinating geographic data acquisition and access: The National Spatial Data Infrastructure[END_REF], INSPIRE 2007). All of them share the same goal: maximize user's access to spatial data, minimize the redundancy of efforts and investments [START_REF] Nebert | Developping Spatial Data Infrastructures: the SDI Cookbook[END_REF]. To implement it in European states, the INSPIRE directive defines a set of rules and specifications and some basic components (figure 1).

Theoretically, data is stored on distributed repositories (Database Management Systems that are extended with the capability of manipulating georeferenced data in most cases). A set of tools is coming in between the repositories and the user applications such as GIS applications or thin web clients. These tools provide the user with an access to a set of services like catalog, geoprocessing, dispatcher that permit to view, query or plan advanced analysis. The interoperability of the system is guaranteed by several standards for data exchange and processing (GML, WFS, WMS, WPS) [START_REF] Smits | INSPIRE Architecture and Standards Position Paper[END_REF].

In practice, the development and exploitation of an SDI is complex and faces many difficulties. Unlike the web client side of an SDI, the server side and desktop clients still deal with the complexity due to the multiple data sources types (WFS, spatial databases, local files...) of typical SDI scenarios. Additionally, this complexity leads to an increasing difficulty to reuse work between the users of the SDI. Finally, there are human factors (new practices) [START_REF] Rajabifard | Developing Spatial Data Infrastructures to Facilitate Industrial and Mining Decision-Makings[END_REF][START_REF] De Oliveira | User Interface Issues in Geographic Information Systems[END_REF] that reduce SDI acceptance by the users. (INSPIRE, 2007) This paper is focused on improving the exploitation of SDI through a semantic layer between data sources and desktop applications. This layer will allow the users to work independently of data source types in the SDI and will let them reuse the semantic functionalities of some other users through a SQL derived language. The goal of this language is to be simple and powerful so that the users can manage the SDI independently of the repositories complexity.

After presenting the aims of the work we will describe the architecture of the solution and its internal model. Finally we will show the semantic layer in action in a theoretical use case of classifying watersheds by its runoff risk.

RELATED WORK

The questions of spatial query language, process layer and spatial semantic repositories have been thoroughly studied over the last years. As reviewed by (Lopes de Oliveira, 1996), the need for spatial query formalism has been clearly identified by and answered in (Egenhofer, 1988a,b), (Güting, 1988), [START_REF] Goh | A Graphic Query Language For Cartographic and Land Information Systems[END_REF], with GEOQL [START_REF] Ooi | Efficient Query Processing in Geographic Information Systems[END_REF], with DML [START_REF] Calcinelli | The Management of the Ambiguities in a Graphical Query Language for GIS[END_REF], with SQL [START_REF] Svensson | Geo-SAL: A Query Language for Spatial Data Analysis[END_REF], GPL (Egenhofer, 1991), Spatial SQL [START_REF] Egenhofer | Spatial SQL: A query and presentation language[END_REF], geoPOM [START_REF] Nittel | geoPOM: a heterogeneous geoscientific persistent object system[END_REF][START_REF] Sha | a Spatial Query Language of Object-Oriented GIS[END_REF]. And there are already several efficient and industrial implementations that spatially enable Relational DataBases Management Systems (such as the PostGIS, a geographic extension for PostgreSQL).

Spatial and alphanumerical data may be located on volatile memory or persistent storage devices. As for databases, common geographic software applications require, in order to be considered usable, implementing a minimum set of features: ability to Create new entries, as well as to Read, Update and Delete existing entries. The Structured Query Language (SQL) is a world famous declarative programming language dedicated to CRUD data from Relational DataBase Management Systems (RDBMS). Beyond this, it is a data-oriented programming and query language that consists of three more command subsets:

•

A standard Data Definition Language (DDL) with CREATE, DROP, ALTER common command,

• A standard Data Control Language (DCL) with GRANT, REVOKE common commands,

• A Transaction Control Language (TCL).

The SQL language has already been spatially enabled and extended by the Open Geospatial Consortium in the OpenGIS "Simple Features Specification for SQL" (Herring 2006a,b). Indeed, this specification defines a spatial object model and fundamental geometric functions, including:

• Spatial predicates (based on the DE-9IM model),

• Overlay functions (intersection, difference, union, symmetric difference),

• Buffer,

• Convex hull,

• Area and distance functions.

It proposes a standard SQL schema that supports storage, retrieval, query and update of simple geospatial feature collections but adapted to RDBMS (Herring 2006a,b). For example, the following query returns a new geometry defined by buffering a distance d around geometry object, where d is in the distance units for the Spatial Reference of geometry.

SELECT Buffer(geometry, d) FROM roads;

In this case roads is a table stored in a Relational DBMS and geometry is a column that contains the value encoded in Well Known Binary format (Herring 2006a).

GDMS [START_REF] Anguix | gvSIG: Open Source Solutions in spatial technologies[END_REF] is related with a similar layer used in the gvSIG project to manage the alphanumeric access. This layer is called gdbms and can only used for alphanumeric purposes. The work presented here is a general refactoring that mainly adds spatial functionalities and a more powerful SQL processor to GDBMS.

Motivations

Frequently, the heterogeneity of source types makes difficult the reuse of algorithms that are tightly coupled with the specified file format, database vendor, etc. With an intermediate layer between the user and the information source, the work developed by the former will not be coupled with the specificities of each format but with the intermediate layer itself, letting the work to be reused in a much more wide set of scenarios and of course simplify the learning curve for new developers.

One of the problems that such a layer arises is that it has to be able to access any potential source type . As the current number of source types is huge (ESRI Shapefile, spatial postgis, oracle spatial, WFS...) and can grow indefinitely, the user has to be able to extend the access capabilities of the layer. By introducing the concept of driver, GDMS intends to potentially fit any situation, even if the types of sources are not well known or currently the layer does not have the driver to access it.

The development of an abstract layer to manipulate spatial data in SQL is not a recent innovation. Several vendors like Esri with ArcSDE or MapInfo with SpatialWare proposed such a middleware. However, these systems have some limitations, due to:

•

Interoperability across data repositories and spatial extensions of database technology

•

Interoperability about SQL semantics. MapInfo used the User DefinedTypes (UDT) framework to build SpatialWare while ArcSDE used Abstract DataTypes (ADT).

Practically, it has an impact on SQL syntax. For example, with ADTs to compute the area of a geometry the user writes "select geometry.area() from mylayer" while with UDTs he write "select area(geometry) from mylayer".

• Complexity of the their installation (need server like SQL Server 2000).

With GDMS, the idea is to postpone these problems by developing a high flexible and portable tool to build SQL queries. Currently, GDMS is dedicated for GIS clients but in the future it will be also a software component for SDI to process spatial data according to Web Processing Service.

GDMS provides an API that lets the user operate independently of source type. However, this API is not friendly enough for an end user and reduces the acceptance of the solution by final users. With the purpose of simplifying both access and manipulation of data sources, GDMS provides a SQL processor that lets the execution of the common Data Manipulation Language statements against any source mapped by a driver. To avoid introducing a new grammar, GDMS fully preserves the SQL-92 grammar and adds to this standard geometric concepts and spatial functions as in OGC simple features SQL specification. As an analogy to spatial SQL for RDBMS, GDMS provides an extended SQL query language on heterogeneous data types.

It is the main purpose of GDMS to improve data creation and sharing. As in an SDI the consumption of data is as important as the production and sharing of data, the SQL processor in GDMS allows data feedback as if they were new data sources (materialized views). This means that the result of SQL queries can easily be integrated into the SDI as a new data source. Those data will be ready to be used by further SQL statements as any other existing data source (figure 3). Finally, to improve further code reuse, GDMS introduces the concept of semantic repository. GDMS allows an extension to the semantics of the SQL language in terms of functions and custom queries. Functions and custom queries are artifacts that contain the implementation of some operations on the data and can be reused just by referencing its name into an SQL statement. This way, some user can implement a buffer operation and other can reuse it just by calling buffer in an SQL query: select buffer(the_geom, 20) from mydata. The collection of all this artifacts is the semantic repository. The purpose of GDMS is to maintain such a repository and encourage the feedback of new artifacts from the user to make it a growing knowledge base (figure 4). The semantic repository is more than a library of functions that can be extended by new user-defined add-ons; we aim to go a step further with this concept. Indeed, the goal is to provide a high level of spatial semantics:

• checking: a « just in time » validation process of input data type, range... made by a dedicated spatial SQL pre-processor

• sharing: Future releases will embed network-groupware capacity because user needs also to share spatial process and not only geographic data

• interoperability:
The GDMS core is fully OGC compliant, so it is possible to copy and paste spatial SQL scripts from/into PostgreSQL/PostGIS.

At last, our semantic repository provides not only spatial functionalities but also description and (sometimes) useful use-cases for each of them. In spite of its extensibility, the driver layer is absolutely heterogeneous and building a SQL processor on top of such a layer would be very complicated. To solve this problem we place the Datasource layer on top of the drivers one. The main responsibility of this layer is the creation of a common API to manage all source types.

The SQL processor is built on top of this Datasource layer so it uses the common interface it provides to access any of the available data source types. Consequently it is possible to mix different types of data sources in the same SQL statement. For example it is possible to join a shapefile with a CSV file and a PostgreSQL table in one single SQL query. To reference the sources in a SQL statement it is necessary to map a name to each involved data source before. Here we have associated a Shapefile with 'shape', a PostgreSQL table with 'postgres' and a coma separated values file with 'csv'. Note that the associated name is the last argument of each "select register..." statement. Sources are referenced with those names in the 'from' clause of the last statement.

SELECT REGISTER ('postgresql','localhost', '5432', 'mygisdb','postgres','postgres','watershed', 'postgres'); SELECT REGISTER ('/tmp/waternetwork.shp','shape'); SELECT REGISTER ('/tmp/codes.csv','csv'); SELECT s.the_geom FROM shape s, postgres p, csv c WHERE s.id = p.id AND p.id = c.id AND c.code = 3 AND p.code = 4;

On top of the stack there is the semantic repository that provides a functional package ready to be reused in new queries by the SQL processor.

Data Model

We have seen that the Datasource layer deals with the heterogeneity of the drivers and provides a common interface to the upper layers. We will now focus on the data model that is used to implement all features. It consists of one main entity called DataSource. A DataSource is an abstraction that allows the management of one unique data source. This means that it is necessary to have as many DataSource instances as the number of sources to be accessed. This abstraction has a typical tabular structure that uses rows to store each of the elements of the data source. This structure is close to the JDBC standard, where each row consists of a set of fields with a common structure for all rows in the DataSource. It is possible to see a relationship between this model and the GeoAPI model (currently being standardized in ISO19123) where FeatureCollection and Feature respectively match DataSource and its rows. Finally the field values of each row match the attribute values of a Feature.

All the information about field types is stored into a Metadata object that contains all field names, field types and restrictions for the types, such as length for string values, and coordinate reference system for spatial types. The model in this point allows a wide range of types for the fields that are compliant with the JDBC standard for alphanumeric types or with the OpenGIS Simple Features Specification for spatial types. Also note that there is no restriction in field types for a DataSource, it can have zero or more fields of any type so it is as possible to have several spatial fields as to have just alphanumerical.

The major components of the GDMS data model are shown in figure 6.

EXAMPLE CASE STUDY

The result of our development is the explained semantic layer. It's an open source project that can be easily embedded in third party applications. Since the usage of the GDMS layer would be difficult without a user-friendly interface, GDMS is part of a broader GIS and SDI project named OrbisGIS [START_REF] Leduc | OrbisGIS : A GIS for scientific simulation[END_REF]. As the rich client application of our SDI, OrbisGIS allows executing one or more queries and permits the presentation of the geometry of the data source and the visualization of issues that are related to spatial objects. For non-spatial features, the user can display the result into tables. In this section we will demonstrate how various SQL expressions can be applied including heterogeneous data source and different kinds of processing: spatial and non-spatial. We suppose a theoretical example where the aim is to classify watersheds comparing to runoff risk. The watersheds are classified on the basis of one arbitrary criterion: The percentage area of cereals being less than 20 meters away from a river network. The example is divided in the 6 steps illustrated in figure 7.

Data and processing schema

1. Accessing data, 2. Filtering data, 3. Identify the parcels that are located less than 20 meters away from a river feature, 4. Extract the parcels that are intersected with watershed, 5. Aggregate all parcels by watershed to compute the total area, 6. Compute density.

Application with GDMS

First, to start processing, the user must associate a name with all the input data sources.

SELECT REGISTER ('postgresql','localhost', '5432', 'mygisDbName','postgres','postgres','watershed', 'watershed'); SELECT REGISTER ('/temp/waternetwork.shp','waternetwork'); SELECT REGISTER ('/temp/landcover2000.shp','landcover2000');

The data are loaded in the OrbisGIS graphical environment and represented with specific colors and symbols (see Appendix for graphical results).

The first step of the methodology (step 2 in figure) involves restricting the data area. The user identifies and displays only the parcels where landcover type is equal to 'cereals' and only the waternetwork where the type is 'rivers'.

SELECT REGISTER ('/temp/rivers.shp','rivers'); SELECT REGISTER ('/temp/cereals2000.shp','cereals2000'); CREATE Before computing the area we will group all the parcels that are in the same watershed keeping their spatial information with an aggregated function called GeomUnion that will compute the union of all geometries.

SELECT REGISTER ('/temp/cerealsUnion.shp','cerealsUnion'); CREATE TABLE cerealsUnion as SELECT GeomUnion(the_geom) as the_geom, gid FROM cerealsLess20Intersects WHERE the_geom IS NOT NULL GROUP BY gid;

Finally, with the aggregated data the user can compute the percentage area of cereals which is less than 20 meters of a river network. The final result is a list of three watersheds with their corresponding cereals density that the user can display with a thematic map analysis.

Extending SQL semantic is one of GDMS main purpose and ability. Indeed, as a matter of fact, one of the first use cases was to use it to produce morphological indicators. Among around ten indicators, compacity [START_REF] Miguet | Urban bioclimatic indicators for urban planners with the software tool SOLENE[END_REF] as the ratio between the exchange areas of the envelope and the floor area -has been implemented as a GDMS function. Compacity can be used to analyze landscape structuring with respect to current agricultural landcover changes.

CONCLUSION AND FUTURE WORKS

With the implementation of GDMS the exploitation of the SDI is improved in terms of usability due to the fact that a simple universal SQL interface is provided to the user to potentially manage all the sources in the SDI. This interface is as close as possible to the SQL-92 standard and to the "Simple Features Specification for SQL" (Herring 2006a,b) to guarantee a smooth learning curve.

In spite of its simplicity from the user point of view, the SQL interface provided by GDMS allows the extension of the semantics of the SQL language. These extensions can be developed and easily reused by the users of the SDI. The set of all these extensions is called the semantic repository and it will improve the cooperation and the reuse of work between the users of the SDI. Additionally, we have seen that the SQL processor is able to execute queries to create data sources of different types. This will let the user of the SDI to enrich the set of sources it has providing a feedback in terms of data through the common SQL interface.

As a further step in terms of simplicity for the user, we have planned to embed the GeoProcess tool. It is a graphic process builder that will let the user to create his own complex data processing sequence. The processes managed by this tool are composed of data sources and SQL scripts. GDMS has the ability to process both of them (execute script against data sources), so it is only necessary to link the inputs and outputs of the scripts between them. The input of the SQL scripts are the names of all data sources involved in the queries while the output results are the names of the produced data sources. As shown in figure 12, labeled rectangles correspond to data sources, while labeled circles correspond to the scripts. The directed edges that connect the nodes symbolize source dataflows. As an example, we provide the process diagram for the previous use case specified in figure 7, and represented in figure 8. This user interface enhancement will not only simplify the definition and the implementation of complex process, but will also let the possibility to export this whole data processing sequence to the semantic repository. Thus, all this complex data processing sequence will be stored ("serialize") in order to be re-used as a "black box" in next processes ("unserialize"). This sequence of serialization and unserialization mechanism in a local or remote semantic repository, will be a new step towards reuse, workflow and even geoinformatic groupware. This will be step towards both data and methods sharing which is the ultimate goal of our scientific SDI.

Figure 1 :

 1 Figure 1: Architecture reference model for the INSPIRE directive (INSPIRE, 2007)

Figure 2 :

 2 Figure 2: GDMS ability to access Data and Semantic repositories

Figure 3 :

 3 Figure 3: Data sources feedback after processing: implementation of materialized views

Figure 4 :

 4 Figure 4: Feedback into the semantic repositories

Figure 5 :

 5 Figure 5: GDMS layered architectureIn spite of its extensibility, the driver layer is absolutely heterogeneous and building a SQL processor on top of such a layer would be very complicated. To solve this problem we place the Datasource layer on top of the drivers one. The main responsibility of this layer is the creation of a common API to manage all source types.

Figure 6 :

 6 Figure 6: GDMS data model

Figure 7 :

 7 Figure 7: Data processing sequence

  SELECT REGISTER ('/temp/density.shp','density'); CREATE TABLE density as SELECT ((area(a.the_geom) / area(b.the_geom))*100) As density, b.the_geom as the_geom, gid FROM cerealsUnion AS a, watershed AS b WHERE a.gid = b.gid;

Figure 8 :

 8 Figure 8: Use case translation in the graphic query builder context

Table 1 :

 1 Three vector spatial data sources are considered: watershed, waternetwork and landcover2000. They are described in table 1. Input data sources

	Name	Geometry	Description	Data repository
	watershed	polygon	A set of 3 watersheds used to compute density.	PostGIS table
	waternetwork linestring	Ditches and rivers in Sterenn watershed. Esri Shapefile
	landcover2000 polygon	Land cover classification in 2000.	Esri Shapefile

TABLE rivers

 rivers The third step calculates the spatial intersection between cerealsLess20Intersects and the watersheds.

	AS SELECT * FROM waternetwork where
	type_axe='rivers';
	CREATE TABLE cereals AS SELECT * FROM landcover2000 where
	type='cereals';
	The second step uses two spatial functions (Buffer and Intersects) to select the parcels that are less
	than 20 meters of river features.
	SELECT REGISTER ('/temp/cerealsLess20.shp','cerealsLess20');
	CREATE TABLE cerealsLess20 AS SELECT a.the_geom FROM cereals AS a,
	rivers AS b WHERE Intersects(Buffer(b.the_geom, 20), a.the_geom);
	SELECT REGISTER
	('/temp/cerealsLess20Intersects.shp','cerealsLess20Intersects');
	CREATE TABLE cerealsLess20Intersects AS SELECT Intersection(a.the_geom,
	b.the_geom) AS the_geom, b.gid AS gid FROM cerealsLess20 AS a,
	watershed AS b WHERE Intersects(a.the_geom, b.the_geom) ;

AVAILABILITY

GDMS alpha version is available for public download at http://sourcesup.cru.fr/projects/orbisgis/

APPENDIX: EXAMPLE RESULTS

The following snapshots are the results of the SQL statements of the example:

• up: initial data • left low: rivers and cereals2000 layer representation • center low: cereals2000 parcels close to the rivers • center right: cereals2000 parcels close to the rivers that intersect watershed


Thomas Leduc
email: thomas.leduc@cerma.archi.fr

Erwan Bocher
email: erwan.bocher@ec-nantes.fr

Fernando González Cortés

Efficient Constrained Delaunay Triangulation implementation in Java for spatial hydrological analysis

To quantify the impact of increased urbanization on the peri-urban hydro-systems, it is necessary to design adapted models. Currently, most of the hydrological models, based on classical raster approach, suffer from a lack of efficiency because they essentially do not take into account the urban artifacts and do not match complex geometries such as urban ones (heterogeneous objects, strong space's partition, roads, pits, buildings...).

That is the reason why, we have decided to choose a Triangulated Irregular Network (TIN) model.

This approach requires a Constrained Delaunay Triangulation (CDT) preprocessing phase. Also CDT algorithms are numerous, we must admit that there is a real lack of free, efficient and robust implementations adapted to our technological choices. Indeed, all the hydrological simulations we have developed, had been made in the context of the OrbisGIS platform (an open source GIS developed in Java in the French IRSTV Institute). The main objective, with OrbisGIS, is to remain coherent and give to the community a consistent set of simulations tools and methods for the spatial analysis of the urban environment in the context of a single and cross-operating systems GIS.

In this paper, we will first present the OrbisGIS platform. The remainder of this paper is dedicated to the CDT algorithm we choose, its implementation and some comparisons with other implementations. The last part of this paper is dedicated to the hydrological process and the first results we have obtained.

Introduction

Due to the increase of urbanization, 60% of the worldwide population is expected to live in towns in 2030. The most important land-use modifications occur in peri-urban zones. A higher concentration of built infrastructures has an impact on hydrology and ecosystems functions. It is therefore necessary to design models to quantify the impact of increased urbanization on the peri-urban hydro-systems.

The main objective of the AVuPUR project (Assessing the Vulnerability of Peri-Urban Rivers) is to handle properly the heterogeneity of peri-urban surfaces. Existing hydrological models will be adapted to include the discontinuities induced by natural features such as the hydrographic network and anthropogenic elements such as roads… using a combination of an object-based description developed for urbanised area and a TIN-based approach, used traditionally for rural areas.

In the hydrology community, there has been a substantial amount of research and development using the classical raster approach. The two basics concepts in terrain analysis for hydrology processing are flow routing (model globally the water flow through the terrain) and flow accumulation (quantify how much water flows through each point). Spatial models are often realized with 2D regular orthogonal lattices, using an uniform grid with an elevation given for each cell. Based on this raster DEM, a flow direction grid can be created with several different methods. The simplest and most effective method is the D8 one [START_REF] O'callaghan | The Extraction of Drainage Networks From Digital Elevation Data[END_REF]. Several methods have then been derived [START_REF] Tarboton | A New Method for the Determination of Flow Directions and Contributing Areas in Grid Digital Elevation Models[END_REF].

This grid approach is the most popular one in hydrology because of its intrinsic simplicity. But, as written in [START_REF] Kenny | A methodology for aligning raster flow direction data with photogrammetrically mapped hydrology[END_REF], "a major limitation […] is that watershed boundaries extracted from these data often have significant errors". That the reason why "stream burning offset", "surface reconditioning process" and similar techniques have emerged.

Currently, most of the models suffer from a lack of efficiency because they essentially do not match complex geometries such as urban ones (heterogeneous objects, strong space's partition...). Moreover, most of the grid approaches introduce, in the input data themselves, a great spatial approximation (due to the artificial rasterization process) and also too few degrees of freedom in term of flow directions.

For all those reasons, we have decided to choose a TIN model. Validation will me made using topological and morphological characterization as described in [START_REF] Gironás | Morphologic approach in studying developing urban watersheds[END_REF]. To take benefit from current software development in the French IRSTV institute, we have also decided to develop the model in the context of our OrbisGIS platform.

The OrbisGIS platform and its GDMS layer

2.1 Overall view of GDMS GDMS (Bocher et al., 2007) is related with a similar layer used in the gvSIG project to manage the alphanumeric data access. This layer is called GDBMS [START_REF] Anguix | gvSIG: Open Source Solutions in spatial technologies[END_REF] and one of its main limitation is that it can only be used for alphanumeric purposes.

As noticed in (Bocher et al., 2007), the heterogeneity of source types makes difficult the reuse of algorithms that are tightly coupled with specific format. With an intermediate layer, the work developed by the former will not be coupled with the specificities of each format but with the intermediate layer itself, letting the work to be reused in a much more wide set of scenarios and of course simplify the learning curve for new developers. With GDMS, the idea is to postpone all problems of interoperability across data repositories, but also about the SQL semantics by developing a highly flexible, portable and standards compliant tool to build SQL queries.

With the purpose of simplifying both access and manipulation of data sources, GDMS provides a SQL processor that lets the execution of the common Data Manipulation Language statements against any source mapped by a driver. To avoid introducing a new grammar, GDMS fully preserves the SQL-92 one and adds to this standard, geometric concepts and spatial functions as in OGC simple features SQL specification.

OrbisGIS in action

All the hydrological simulations we have developed, have been written in the context of the OrbisGIS platform. The main objective, with OrbisGIS, is to remain coherent and give to the research community a consistent set of simulations tools and methods for the spatial analysis of the urban environment in the context of a single and cross-operating systems GIS.

Preprocessing phase: the Constrained Delaunay Triangulation (CDT)

Preliminaries

As mentioned in the introduction, the methodological choice we made, lead us first to develop a CDT preprocessing phase. There are a wide amount of algorithms available to build a triangulation for a set of points. To obtain an optimal set of triangles, without adding Steiner points (that is to say vertices of the mesh that are not vertices of the input), the most frequently applied method is the one of Delaunay.

A quite exhaustive review of sequential DT algorithms has been made in [START_REF] Su | A comparison of sequential Delaunay triangulation algorithms[END_REF]. Given a finite set of points in the plane, three points contribute a triangle to the Delaunay Triangulation (DT) if the circumcircle through those three points contains no other point in its strict interior. This is the wellknown Delaunay Property (DP). Such a definition produces an optimal set of well-shaped triangles.

Moreover, the produced triangulation is unique (assuming careful handling of the degeneracies: 4 or more co-circular points) and so repeatable.

In addition to point, linear constraints may also be given as input datas (as in a set of iso-height lines). The objective is to be able to mesh any planar straight-line graph (PSLG). PSLG includes polygons with or without hole(s), dangling edges and isolated vertices. This generalization of the original DT is called the CDT: constraining input segments must become edges of the generated CDT.

The DP itself has to evolve into a weaker version: if the circumcircle of a triangle contains a fourth point, this one must be located behind a constraining edge.

The algorithms for building CDT can be classified in two groups:

the first one considers separately points and line-segments. First process is to mesh given set of points with a classical DT (choosing a strategy: divide and conquer, sweep-line, incremental, high dimensional embedding…) and then add the linear constraints respecting the weaker DP;

the second one processes points and edges simultaneously. [START_REF] Domiter | Sweep-line algorithm for constrained Delaunay triangulation[END_REF] notices that "such algorithms are faster, consume less memory and are less sensitive to data distribution".This is the strategy that is chosen in the algorithm described in [START_REF] Domiter | Sweep-line algorithm for constrained Delaunay triangulation[END_REF].

The Sweep-line CDT algorithm

The CDT algorithm we have implemented uses a sweep-line paradigm combined with Lawson's legalization. As written in [START_REF] Domiter | Sweep-line algorithm for constrained Delaunay triangulation[END_REF], the proposed algorithm consists of three parts:

an initialization phase:

input data points are first sorted. Each point is also associated with the information whether or not it is the upper ending point of one or more edges;

then, two sensible artificial points, an initial triangle and an initial sweep-line (SL) are created,

a sweeping phase: following the sweep-line paradigm, it stops at every point and check its status in case of a point event (isolated or lower edge point):

the point's vertical projection on the SL is computed. Based on this projection point, one or two triangles are generated and legalized. The SL is updated;

to avoid undulation of the advancing front, basin are filled by triangles,

in case of an edge event (upper edge point):

first, each triangle that intersects a constraining edge is removed (taking into account that the SL may be intersected at one or more points) and the SL is updated,

then, what still remains is to re-triangulate the empty area.

a finalization phase: when all points have been swept, two ending tasks still remain removing all the triangles defined by at least one artificial point,

adding the bordering triangles (the union of them all must form the convex hull of input data).

Our own implementation

For coherency reason, we have implemented this CDT algorithm in the OrbisGIS platform, that is based on the robust JTS (Java Topology Suite) library. This library provides spatial indexes that we have included: a MX-CIF quadtree for the triangles and a query-only R-tree created using the STR algorithm for the immutable set of input data. However, the implementation we made is loosely coupled to the GDMS layer. Indeed, the interface part has just to do with inputs (reading the PSLG) and outputs data (storing the resulting triangles).

The asset of the platform is its ability to make the development of spatial tools easier. Thus, it is already possible to verify the validity of the DP for each resulting triangle with a single (GDMS) SQL query: select CheckDelaunayProperty() from mySetOfTriangles. We have made it also possible to compare two sets of triangles (CheckSpatialEquivalence), measure the quality of the produced triangles (QualityMeasuresOfTIN)…

Before implementing the sweep-line algorithm, we have also tried to interface some existing meshing tools. Thus, we have tested, in the platform context, the famous Triangle program developed in C language by Shewchuk, but also a Java implementation made by M. Davis (WaterBuG project), another one written by P. Austin and, at least, a really efficient one written by M. Michaud. As an example we may focus on the Shewchuk's implementation. It is the only one we tested that give the possibility to refine a mesh. As shown in figure 1, the Steiner points that are added increase the minimum angle value (from 11.4 to 21.4 degrees) and therefore produce better shaped triangles.

Hydrological simulations

As written in [START_REF] Mcallister | The Computational Geometry of Hydrology Data in Geographic Information Systems[END_REF], every sample point of a Triangulated Irregular Networks (TINs) must store not only its coordinates and elevation values but also face adjacency information. Starting from a set of triangles, we have to build a topographic graph. The data structure we adopt, is the one described in [START_REF] Van Kreveld | Digital elevation models and tin algorithms[END_REF]. The next step is to build a hydrological graph using this topographic one.

The idea is to determine the contributing area for a river 1 . It presupposes that for each graph node a steepest descent direction is computed (we assume that water always follows the path of steepest descent). The hydrological graph is produced from the topographic one according to the drainage characteristics. All the topographic vertices that contribute to a same common local minimum in the terrain are grouped. These groups are the watersheds. 1 The watershed of a river is the set of terrain points that drain into this river.

Conclusion

We have presented an implementation of a new CDT algorithm and its first application in a TINbased hydrological model. There are several directions for further work. The first one has to do with a CDT refinement algorithm, in order to obtain better shaped triangles. We have also planned to develop a GDMS native format to store both the topographic graph and the hydrological graph.

Figure 1 .

 1 Figure 1. The influence of the Delaunay refinement algorithm on the triangulation. Shewchuk's implementation (on the right side) introduces some Steiner points to obtain a greater minimum angle.

Figure 2 .

 2 Figure 2. Starting from a set of set of iso-height lines and a set of linear constraints, a CDT is build (left side). The watershed is derived from the TIN, through a hydrological graph (right side).

Acknowledgments

The AVuPUR project was funded by the French Agence Nationale de la Recherche (ANR) under contract n°ANR-07-VULN-01.


Nathalie Long

Erwan Bocher

Thomas Leduc

Guillaume Moreau

Sensitivity of spatial indicators for urban terrain characterization

Keywords: urban fabric characterization, indicators, UrbSAT sensitivity. a b c. d

come   L'archive ouverte pluridisciplinaire

INTRODUCTION

At present more than half of the world population is living in urban areas. Therefore urbanization control has become a major issue. Among the social, economical and ecological factors that are taken into account in decision-making, environmental issues now have a greater impact. To address these issues, researchers are faced with new problems among which are tool disparity, knowledge sharing and enhancement. GIS are become unavoidable tools for the development and management of territories and to analyze urbanization impact on environment. An accurate knowledge of the city and its urban fabric is necessary. A physical description like building morphology or land cover allows some characterization of the urban fabric. A specific package for a GIS software has been developed for urban analysis: UrbSAT. It allows extracting knowledge from various data sources and gives some indicators to real applications (sustainable development, management policy, air quality improvement, urban canopy modeling…). This paper proposes a sensitivity analysis of the spatial indicators produced by UrbSAT. Following the introduction, section 2 relates the framework and the methods used in other studies, section 3 presents the UrbSAT software, the databases and the sensitive tests and then, the results of relevant issues end up the paper with future research.

RELATED WORK

The UrbSAT development is based on the DFMap methodology, developed within a French laboratory of Ecole Centrale de Nantes [START_REF] Long | Developement of a software to describe the city morphology and to compute aerodynamic parameters from an urban database[END_REF]. It allows computing not only statistical parameters describing the urban fabric but also aerodynamic parameters like roughness length or displacement height. From the BDTopo database (produced by the French National Geographic Institute, IGN), DFMap converts the initial data into a set of regular squares (vector grid) and computes several functions over each features of the grid. This prototype presents some limitations: among other things, only the BDTopo can be used as data source and GIS facilities (reprojecting a map, rewriting geometric computations or map generation) are not available. To improve this, the OrbisGIS open source GIS software was developed for which urbSAT is a plugin package dedicated to urban analysis.

To characterize the urban fabric, several researchers have developed automated or semi-automated methods to extract urban indicators from database or satellite images. Grimmond & Souch are the first one to use GIS techniques to compute urban indicators for urban climate analysis [START_REF] Grimmond | Surface description for urban climate studies: a GIS based methodology[END_REF]. Then, a small community of American scientists has developed the National Urban Database and Access Portal Tool (NUDAPT) which produces gridded outputs of urban parameters capable of driving urban meteorological models [START_REF] Ching | National urban database and access portal tools (NUDAPT): a project overview[END_REF]. This prototype uses the ESRI ArcGIS software package. With a remote sensing and database approach, Huang et al., propose an analysis of urban form from spatial indicators like compactness, porosity, centrality to compare cities of developing and developed countries [START_REF] Huang | A global comparative analysis of urban form: applying spatial metrics and remote sensing[END_REF]. Then, Longley & Mesev have developed a fine scale data model of population densities to analyze urban sprawl [START_REF] Longley | On the measurement and generalization of urban form[END_REF].

URBSAT SOFTWARE SENSITIVITY

In this section, we present the UrbSAT software and the database used to analyze two French cities: Le Havre and Nantes, as well as the tests realized to evaluate the UrbSAT sensitivity and more precisely, the computation methods.

UrbSAT software

To improve the urban fabric knowledge, the main idea of UrbSAT is to process vector data and to aggregate information on a regular vector grid. It uses a generic processing library called GDMS (Generic Data Source Management System) [START_REF] Bocher | GDMS: An abstraction layer to enhance SpatialData Infrastructures usability[END_REF]. It is a semantic abstract layer between the users and the data sources that provides spatial SQL language-based interface. It is thus possible to extend the semantics from typical SQL with customized functions and queries, developed in Java.

Urban data sources

For this study, the French cities of Nantes and Le Havre are the chosen application cases. The BDTopo database and satellite images are available as information sources. The BDTopo database is produced by the French National Institute (IGN) and describes the elements of the urban fabric in details. It includes an inventory of buildings, vegetation areas, roads, and hydrographic network, among other objects. The outlines of these objects are delineated by polylines including some attributes like building height.

The SPOT satellite images of Nantes have a 2.5m resolution and were taken in May, 2000. A supervised classification is realized on each district selected on the image to distinguish between built, vegetation and road surfaces (treatments realized with ENVI software RSI).

Sensitivity tests

Several tests are realized and presented here to evaluate the quality and accuracy of urban database, first. Then, on the raster process, the grid cell size and its orientation are tested and evaluated, compared to the source information. All these tests are computed with the UrbSAT software.

RESULTS

Quality and accuracy of urban database

To evaluate urban database quality to describe the land cover, a very high-resolution satellite image is used on Nantes. Three districts of 1km² are selected from both data sources: a city-center with high built density, a collective building district and a residential district (figure 1). From a supervised classification, the density of different types of surface are computed from satellite image and compared to the UrbSAT results from the BDTopo.

Table 1 allows evaluating the under-estimation of land cover from the BDTopo. The difference between these two data sources is low in city center but higher in the districts where the built density is lower. In the BDTopo, all the small surfaces (building surface < 50m², vegetation surface < 80m²) are not present in the database. This lack does not allow describing 100% of the land cover. The satellite image analysis gives better results but no 3D information is directly available.

Grid cell size tests

From vector to grid process, as for any sampling, information is modified and the result quality is strongly dependant on the cell size [START_REF] Harris | Giving dimension to point locations: urban density profiling using population surface models[END_REF]. To evaluate this influence, tests are realized on Le Havre with the mean building height variable. In the BDTopo, each building has an height attribute; an average height is then computed for each cell from the buildings that are totally or partially included in the cell. To improve the quality of this urban indicator, the average building height may be weighted by the surface but a. b.

c.

it's not the case in this test. Four districts are selected in Le Havre: a city centre, a collective building district, a residential district and an industrial area (fig 5). For each district, an average building height per cell (ABH) is computed for 100m, 200m, 300m and 500m cell sizes. These results are compared to an average of initial building height (IBH), before the raster process, as a reference.

Fig. 3: Statistical variables of the building height indicator for the city center

The vector to grid process has some consequences on data quality: first, the ABH per cell is higher than the IBH. An over estimate of half a meter to 3.6 m is shown with respect to district type and resolution. The vector grid process influence is lower for the residential district (fig 4) and higher for the city centre (fig 3). The standard deviation (std) shows a homogenization of the urban fabric by the grid process: the higher the cell size the more erased the differences. We must find the best solution to distinguish between the different districts in the city and describe it. For small structures like the residential buildings, a 200m or 300m cell sizes appears to be a satisfactory choice, but for big buildings, a 300m or even 500m for the city center or the industrial area (fig not shown)) cell size appears necessary to describe these districts (fig 2 and3). The results from Table 2 show an improvement of the urban indicator quality when they are computed on an oriented grid, adapted to the main direction of the street network.

Circular grid

The radio-concentric mesh may be interesting for particular applications. For example, for flux simulations over a city, a characterization of building morphology and land cover on the main wind direction is essential to define the flux footprint, (this term describes an upwind area "seen" by the instruments measuring vertical turbulent fluxes). The circular grid allows characterizing this upwind area. The user selects in the grid the cells in a fixed direction, like in figure 6. For the two directions, the built surface density is respectively 0.10 and 0.14 and the vegetation surface density is 0.40 and 0.30, in average. The circular grid appears more appropriate than the rectangular grid for the cell selection.

CONCLUSIONS AND PERSPECTIVES

These different points allow evaluating the UrbSAT sensitivity to define indicators characterizing the urban fabric. The vector to grid process has a strong impact on the data quality. The cell size may be adapted to each urban structure because of the building size variation. In the same way, the grid orientation seems important, and when the grid orientation is the same as the network street, urban indicators are closer to reality.

The vector to grid process is mainly tested here and this study must be supplemented by tests on indicator computation methodology. In other way, the development and the validation of a function defining the main directions of the street network appears essential to complete these first results and automate the process.

Fig. 1 :

 1 Fig. 1: Test areas of Nantes (a.: collective buildings, b.: city centre, c.: residential district, polygons: buildings)

Fig. 2 :

 2 Fig. 2: Statistical variables of the building height indicator for the collective building district

Fig. 4 :

 4 Fig. 4: Statistical variables of the building height indicator for the residential district

Fig. 5 :

 5 Fig. 5: Oriented grid for four districts (a. residential district, b.: industrial area, c.: collective buildings, d.: city centre, polygons: buildings, Le Havre)

Fig. 6 :

 6 Fig. 6: Circular grid with two flux footprints in Nantes (in dark: vegetation area, in grey: buildings, lines: footprints)

Table 1 : Vegetation and built surfaces densities from BDTopo and SPOT image analysis on three districts

 1

		City center	Collective	Residential
			buildings	buildings
	SPOT veg. density	0,23	0,25	0,13
	BDTopo veg. density	0,17	0,08	0,15
	SPOT built density	0,47	0,37	0,33
	BDTop built density	0,38	0,18	0,18

  ).

		Angle	ABH	ABH	IBH
			(oriented	(NS/EW)
			grid)
	Collective 355°	13,63m	14,12m	13,17m
	Industrial	340°	11,38m	11,52m	10,26m
	Residentia	40°	6,28m	6,47m	6,38m
	l
	City	355°	15,14m	15,58m	12,33m
	center

Table 2 : Average building height per cell and for initial building for four districts with different grid orientation

 2