Jules Courjault

Baptiste Vrigneau

Matthieu Gautier

Olivier Berder

Accurate LoRa Performance evaluation using Marcum function

Keywords: LP-WAN, LoRa transmission, Binary Error Probability performance analysis, Marcum function

I. INTRODUCTION

During the last decade, the number of small connected nodes grew up exponentially for various applications such as environmental monitoring or connected farms [START_REF] Aoudia | Long-Short Range Communication Network Leveraging LoRa and Wake-up Receiver[END_REF]. Many emerging solutions recently arised to cover the needs for long range and very low energy consumption wireless transmissions. A panel of standards for Low Power Wide Area Network (LP-WAN) was presented in the paper [START_REF] Sanchez-Iborra | State of the art in LP-WAN solutions for industrial IoT services[END_REF]. Among the candidates, LoRa is a long range communication technology promoted by the LoRa Alliance. The LoRa physical layer, which allows long-range, low-power and low data rate communications, is developed by Semtech [START_REF] Seller | Low power long range transmitter[END_REF] [START_REF] Semtech | SX1272/73-860 MHz to 1020 MHz low power long range transceiver[END_REF]. LoRa communications can use frequency bands of 433 MHz, 868 MHz or 915 MHz and its data rate can reach up to 50 kbps. It uses the Chirp Spread Spectrum (CSS) modulation with a linear variation of frequency on time [START_REF] Springer | Spread spectrum communications using chirp signals[END_REF] [START_REF] Reynders | Chirp spread spectrum as a modulation technique for long range communication[END_REF], which allows the reduction of both interference and Doppler effects [START_REF] Le | Energy/reliability trade-off of LoRa communications over fading channels[END_REF]. LoRa modulation can be configured with three parameters: the bandwidth BW , the spreading factor SF , and the coding ratio CR. Although LoRa is described in the patent, rigorous theoretical studies of this technology are still missing. The author in [START_REF] Vangelista | Frequency shift chirp modulation: The lora modulation[END_REF] proposed to rigorously study the modulation and demodulation in a mathematical approach but, as highlighted in [START_REF] Elshabrawy | Closed form approximation of LoRa modulation BER performance[END_REF], the paper still lacks a theoretical analysis of the Bit Error Probability (BEP). As explained later in the present paper, the problem is less to find the theoretical expression, it was already done in [START_REF] Proakis | Digital Communications[END_REF], than an approximation easily computable. All proposed expressions in the literature are based on the complementary error function, while we propose an accurate approximation of the BEP for an Additive White Gaussian Noise (AWGN) channel based on the Marcum function [START_REF] Marcum | Table of Q functions[END_REF] which was already employed in communications theory [START_REF] Simon | Exponential-type bounds on the generalized Marcum Q-function with application to error probability analysis over fading channels[END_REF].

The rest of the paper is organized as follows: Section II presents the required basics and the state of the art of LoRa BER approximations. Section III introduces the novelty based on a Taylor expansion and the main results with Marcum function. Section IV compares proposed results to the state of the art before the conclusion in Section V.

II. STATE OF THE ART

A. LoRa BEP theoretical expression and issues

For conciseness, only the needed equations are introduced and the reader can refer to [START_REF] Vangelista | Frequency shift chirp modulation: The lora modulation[END_REF] and [START_REF] Elshabrawy | Closed form approximation of LoRa modulation BER performance[END_REF] for detailed explanations. We remind that the receiver first applies a down chirp, secondly performs the discrete Fourier transform of the receive signal, and third estimates the symbols with the maximum value criterion. As the symbols have the same probability and the problem is equivalent to a permutation, the study can focus on the probability of correct decision of the symbol m, denoted P m . According to [START_REF] Elshabrawy | Closed form approximation of LoRa modulation BER performance[END_REF], [START_REF] Proakis | Digital Communications[END_REF], this probability can be expressed as follows:

P m = F χ2 | √ N + W p [N -m]| 2 σ 2 N -1 , (1)
where N = 2 SF , SF is the spreading factor and is an integer between 7 to 12, and W p [N -m] is a Gaussian complex noise of variance σ 2 . By using the cumulative density function (cdf) of a χ 2 random with 2 degrees of freedom F χ2 (x) = 1 -exp(-x σ 2 ), eq.( 1) can be written:

P m = 1 -exp | √ N + W p [N -m]| 2 σ 2 N -1 . (2)
By defining the Signal-to-Noise Ratio SNR = 1/σ 2 and the random variable

Z = | √ N + W p [N -m]| 2 σ 2 , (3)
Z follows a non central χ 2 law with 2 degrees of freedom and the noncentrality parameter λ = 2N • SNR. Its probability density function (pdf) is expressed as:

f Z (x) = 1 2 exp - x + λ 2 I 0 √ λx , (4)
with I 0 (.) the modified Bessel function of the second kind of order 0. Finally, the last step is the computation of the average Symbol Error Probability (SEP):

SEP = +∞ 0 1 -1 -e -z/2 N -1 f Z (z)dz, (5) = 1 N N -1 k=1 C k+1 N (-1) k+1 e -N •SNR k (k+1) , (6)
with C k+1 N = N !/(k+1)!(N -k-1)! the binomial coefficient. However, the numerical computation of this result suffers from accuracy issue and the results can not be obtained when SF is greater than 6. This problem was already pointed out in [START_REF] Elshabrawy | Closed form approximation of LoRa modulation BER performance[END_REF], [START_REF] Proakis | Digital Communications[END_REF], and justifies the crucial need for finding accurate approximations. Generally, the performance is expressed with the BEP and due to equiprobability of the wrong symbol, the following relation can be used:

BEP = 2 SF -1 2 SF -1 SEP. (7)
B. Existing approximations of the BEP Reynders et al. proposed in [START_REF] Reynders | Range and coexistence analysis of long range unlicensed communication[END_REF] the approximation:

BEP ≈ Q log 12 (SF ) √ 2 SNR , (8)
with Q(x) the Q-function or the tail distribution of the standard normal distribution. Another study [START_REF] Reynders | Chirp spread spectrum as a modulation technique for long range communication[END_REF] gives the following expression that relies on numerical fitting:

BEP ≈ 0.5 • Q 1.28 √ SF.SNR -1.28 √ SF + 0.4 . (9)
Recently, the authors in [START_REF] Elshabrawy | Closed form approximation of LoRa modulation BER performance[END_REF] proposed two new accurate estimations of the BEP for LoRa communications in AWGN channel and a Rayleigh channel. The principle is to approximate a Ricean distribution by a Gaussian one and the expression for AWGN channel is:

BEP ≈ 0.5 • Q     √ N • SNR -H 2 N -1 -π 2 12 1/4 H N -1 -H 2 N -1 -π 2 12 + 0.5     , (10)
where H m is the m th harmonic number which can be approximated by H m ≈ log(m) + 1 2m + 0.57722. The second approximation is obtained from the former with some simplifications and approximations:

BEP ≈ 0.5 • Q √ 2N • SNR - √ 1.386 • SF + 1.154 . (11)

III. CLOSED-FORM USING MARCUM FUNCTION A. Approximation of P m

As introduced above, the previous works focused on approximating Chi squared or Ricean distributions with a Gaussian one. The novelty of this work is to find an approximation of g(z) = 1 -1 -e -z/2 N -1 while keeping the non-central Chi squared distribution. When g(z) is extended in a sum, the term e -kz/2 is less significant when k and z increase. Thus, the sum can be limited to terms. It is equivalent to find the approximation thanks to the Taylor's series for high values of z and the first major expression is obtained:

g(z) ≈ g(z) = k=1 (-1) k+1 C k N -1 e -kz/2 . ( 12
)
Fig. 1 shows the comparison with the exact g(z) and different approximation orders, = 1 to 7. The good point is that binomial coefficients are limited to low values of k, which solves a computational issue. Since N is high, especially for SF = 12, the binomial coefficient C k N -1 will suffer from loss of precision, typically when k is more than 5. Although the approximation is tightly close for high values of z, the divergence is problematic: the approximation tends toward infinity for odd values of and toward zero for even orders. Nevertheless, a simple and efficient approximation based on a piecewise definition is obvious in order to obtain the initial value 1:

g(z) 1 if z ≤ z c g(z) if z > z c . ( 13
)
where z c is a threshold depending on the values of and SF that will be discussed in the next paragraph. The different values of z c are available in the Fig. 1 and a quick observation is that the piecewise approximation should be accurate for low values of (inferior to 7). The next step is to compute the average as in eq.( 5):

SEP = zc 0 f Z (z)dz + ∞ zc k=1 (-1) k+1 C k N -1 e -kz/2 f Z (z)dz. ( 14
)
By using the definition of the pdf in (4) and the variable changes x 0 = √ z and x k = z(k + 1), the SEP can be expressed by this practical result after some mathematical manipulations:

SEP = 1+ +1 k=1 C k N N (-1) k e -N •SNR(k-1) k Q 1 2N k SNR, kz c , (15)
with Q 1 the Marcum function of order 1 defined by:

Q M (a, b) = ∞ b x x a M -1 e -x 2 +a 2 2 I M -1 (ax)dx. (16)
B. How to find z c ?

Almost all parameters are known for computing (15) except the value of z c . The problem is not trivial and we propose a solution for low values of the order . One should note that the results in the next section and the above remark confirm that the order does not need to be greater than 7. First, by applying the variable change ζ = e -z 2 in (12), solving g(z) = 1 is equivalent to find the roots of a -order polynom. Since the approximation with even values of can not reach the value 1, the study is limited to odd value cases. The only known assumption is that z c is the unique positive real root. However, we are able to obtain the exact solutions for = 1 and = 3. By using the notation ζ c ( ) = e -zc ( ) 2 , we have:

ζ c (1) = 1/(N -1) (17) ζ c (3) = τ - N -4 (N -2)(N -3) 2 1 τ + 1 N -3 (18) τ = (N -4)(N -5) (N -1)(N -2)(N -3) 3 + √ 2(N -4) (N -1)(N -2) 1.5 (N -3) 1.5 1 3 (19) τ (1+ √ 2) 1 3 N (20)
The last simplification for τ is obtained when N is larger than 1. A straightforward numerical observation leads to a heuristic polynomial function of order 1 that gives the exact values of ζ c and z c for = 1, 3 and an approximation for higher values:

ζ c ( ) = α 1 + α 0 (21) α 1 = (ζ c (3) -ζ c (1))/2, α 0 = ζ c (1) -α 1 (22)
Concerning the even cases, we propose to define z c as the value maximizing g(z), i.e. z c is the solution to the equation

∂ g(z)
∂z = 0. The result has also a polynomial form and the same study for odd cases can be applied. However, for conciseness, we proposed the heuristic approximation z c ( ) = z c ( -1) when is even, as shown on Fig. 2. Fig. 2 shows the comparison between a numerical solving of z c and the proposed approximation. We can observe that the approximation is still tight to the numerical value for = 5, 7.

IV. RESULTS AND COMPARISONS

A. Impact and choice of the order

The BEP approximation is illustrated by Fig. 3 for different values of the order = 1 to 7 and SF = 7 to 12. It is obvious that the approximation converges very fast: if the order = 1 is still visible, the other curves are mingled. A quick conclusion is that the order = 1 is already a good approximation with the best trade-off between complexity and accuracy, and it is not necessary to increase higher than 3. In order to confirm this, a convergence parameter is defined by:

∆ BEP,i→j = BEP| =i -BEP| =j / BEP| =i . ( 23
)
It allows the evaluation of the BEP variation when the order is increased. Fig. 4 shows three cases ∆ BEP,1→3 , ∆ BEP,3→5 , and ∆ BEP,5→7 as a function of SNR for different values of SF . The variation when the order increases from 1 to 3 is inferior to 10%, the variation is lower than 3% when the order changes from 3 to 5, and less that 1% for the next step 5 to 7. Moreover, Table I shows maximum values of ∆ BEP,i→i+1 for i = 1 to 6. Since the approximation for even is less accurate, their contributions are less significant than the odd ones. It means that odd orders should be used, typically = 3.

SNR(dB)

B. SNR-corrected zero-order approximation

The complexity computation is linked to the number of Marcum functions in the formula, which is equal to + 1. It means that the lowest approximation, i.e. = 1, needs to compute two Marcum functions. In order to limit the complexity, we propose a heuristic zero order approximation, i.e. only one Marcum function, based on the case = 1 associating the same bound z c (1) and a correction γ applied to the SNR:

SEP = 1 -Q 1 2N γSNR, 2 log(N -1) . (24)
A straightforward numerical study based on minimizing the mean square error was done in order to fit the zero-order expression to the 7th-order one as a function of SF . The obtained numerical values are given in Table II and BEP curves are plotted in Fig. 5. The first observation is that the values of γ depend on SF . A numerical approximation was found that gives the tractable form:

γ √ 2SF /20 + 0.681 (25)
The BEP comparison shows that the approximation is still close to the 7th order. It is confirmed with the relative error plotted in the Fig. 6. The error is less than 10% but it increases very fast at high SNR: the diversity order is different and the values diverge. Moreover, since the error is plotted in absolute values, an inflexion point is visible that corresponds to the equality. In conclusion, the approximation is close to the exact value but is available for a given BEP range, i.e. a value of BEP inferior to 10 -6 in our proposal.

C. BEP comparisons

Fig. 7 compares the proposed approximation with simulated BEP and the two results of [START_REF] Elshabrawy | Closed form approximation of LoRa modulation BER performance[END_REF] which are the best of the state of the art. These simulations validate that the Marcum-based approximation is tightly accurate with a low order = 5. At high SNR, our proposition and [9, eq.( 21)] are equivalent and close to simulation but [9, eq.( 23)] is less accurate. However, differences are visible at low SNR: both solutions of [START_REF] Elshabrawy | Closed form approximation of LoRa modulation BER performance[END_REF] are upper bounds while our approximation is still close to simulations. In order to highlight these differences, Fig. 8 proposes to focus on the relative error of approximations with respect to the simulated BEP. The approximations [9, eq.( 21)] and [9, eq.( 23)] are compared to our proposed solution with = 0, 3, 5. The observation of the curves confirms that the Marcum-based expression is more accurate than both [9, eq.( 21)] and [9, eq.( 23)] and increasing the order enhances the accuracy. At SNR=-15dB, the approximation error for a spreading factor of 7 is about 1%, which represents a gain of almost 10 compare to previous works. On the other hand, the heuristic approximation for = 0 achieves an error quasiequivalent to = 3 and might propose the best trade-off complexity/accuracy.

V. CONCLUSION

In this paper, a new approximation of BEP for LoRa transmissions is proposed for AWGN channel. Instead of approximating a Ricean or Chi squared distribution thanks to a Gaussian one, we proposed to keep the noncentral Chi squared distribution, leading to a Marcum-based form. Simulation results confirm that the proposed solution is more accurate than the state-of-the-art for a full SNR range and may be considered as a numerical reference of the LoRa BEP for helping the design of a communication. On the other hand, a heuristic simplified version was derived in order to keep only one Marcum function to compute. Fig. 8: Relative error between the proposed approximation and [START_REF] Elshabrawy | Closed form approximation of LoRa modulation BER performance[END_REF] with respect to the simulation for = 0, 3, 5 and SF = 7, 12.

  ) and g(z)Values of z Fig.1: Comparison of the exact function g(z) and approximations with different orders from 1 to 7 and SF = 7, 12.Values of z c are plotted with cross marker.

12 SF = 7 SF = 11 SF = 10 SF = 9 SF = 8 Fig. 2 :

 1271110982 Fig.2: Approximation of the threshold z c : numerical solutions are plotted with crosses and the approximations with circles.

7 Fig. 3 :

 73 Fig.3: BEP approximation for order = 1 to 7 and SF = 7 to 12.

7 Fig. 4 :

 74 Fig.4: Accuracy loss when increasing approximation order ∆ BEP,1→3 , ∆ BEP,3→5 , and ∆ BEP,5→7 as a function of SNR for SF = 7 to 12.

Fig. 5 :

 5 Fig. 5: BEP comparison of the SNR-corrected zero-order approximation (cross marquer) with the 7th one for SF equal to 7 to 12 (line).

Fig. 6 :

 6 Fig.6: Relative error of the SNR-corrected zero-order approximation with the 7th one for SF equal to 7 to 12.

Fig. 7 :

 7 Fig.7: Comparison of the proposed approximation with simulation and[START_REF] Elshabrawy | Closed form approximation of LoRa modulation BER performance[END_REF], the order is 5 and SF = 7 to 12.

ACKNOWLEDGMENT

  This work was supported by the ALAMO project labelled by the competitiveness Pole "Image and Networks" and founded by Région Bretagne and Lannion Trégor Communauté.

TABLE I :

 I Numerical values for the maximum value of ∆ BEP,i→i+1 when SF = 7 to 12, i = 1 to 6

	i	1	2	3	4	5	6
	max ∆ BEP,i→i+1 11.6% 0.9% 2.4% 0.4% 0.76% 0.17%

TABLE II :

 II Numerical values for SNR correcting coefficient

							-30 10 -5	-25	-20	-15	-10	-5
	SF	7	8	9	10	11	12
	γ	0.868	0.882	0.894	0.905	0.915	0.924


Carole Scheifer
email: carole.scheifer@aphp.fr

Camille Rolland-Debord
email: camille.rolland-debord@aphp.fr

Edgar Badell
email: edgar.badell-ocando@pasteur.fr

Florence Reibel
email: florencereibel@yahoo.fr

Alexandra Aubry
email: alexandra.aubry@upmc.fr

Alice Perignon
email: alice.perignon@aphp.fr

Olivier Patey
email: olivier.patey@chiv.fr

Sylvain Brisse
email: sylvain.brisse@pasteur.fr

Eric Caumes
email: eric.caumes@aphp.fr

Réémergence de Corynebacterium diphtheriae Re-emergence of Corynebacterium diphtheriae

Keywords: Corynebacterium diphtheriae, migrants, voyageurs Corynebacterium diphtheriae, migrants, travelers Highlights

We report two imported cases of Corynebacterium diphtheriae infection in

Paris in 2016 These two cases illustrate the epidemiology of diphtheria (Afgan refugee, French traveler), the clinical presentation (cutaneous, respiratory carriage), the potential consequences (toxin positive, toxin negative) The rising number of unvaccinated or incompletely vaccinated individuals, increased lnumber of imported cases, and growing vaccine skepticism could put Europe at risk of diphtheria. Résumé Introduction. La diphtérie est une maladie réémergente en Europe. En 2015, 36 cas ont été rapportés contre 53 cas de 2000 à 2009. Patients. Nous rapportons deux cas d'infection à Corynebacterium diphtheriae survenus en 2016 dans un hôpital français : une infection cutanée avec recherche de toxine négative chez un voyageur et un portage respiratoire avec recherche de toxine positive chez un réfugié afghan atteint d'une tuberculose pulmonaire. L'histoire vaccinale était inconnue pour le patient afghan.

Introduction

Corynebacterium diphtheriae (CD) is re-emerging in Europe. A total of 53 cases were reported between 2000 and 2009 [START_REF] Wagner | Diphtheria in the postepidemic period, Europe, 2000-2009[END_REF] whereas 48 cases were declared between 2010 and 2013 in Europe [START_REF]Annual Epidemiological Report 2016 -Diphtheria Stockholm[END_REF][START_REF]A case of diphtheria in Spain[END_REF]. Twenty-four cases of CD infections had also been reported in 2014 [START_REF]Annual Epidemiological Report 2016 -Diphtheria Stockholm[END_REF] and 36 in 2015 [START_REF]WHO vaccine-preventable diseases: monitoring system 2016 global summary[END_REF] in Europe. In 2016, a fatal case of diphtheria in an unvaccinated child of immigrant parents was reported in Belgium [START_REF]A fatal case of diphtheria in Belgium, 24 March 2016[END_REF]. We report two imported cases of CD infection declared in Paris in 2016: a cutaneous infection in a returning French traveler and an asymptomatic diphtheria carriage in an Afghan refugee.

Case reports

A 38-year-old French male presented to our outpatient clinic in Paris with a three-month history of skin ulcers of the anterior right upper leg after a five-year professional stay in the capital of Madagascar. The physical examination revealed one crusted, large, superficial ulcer surrounded by a few smaller ones. No discharge, foul smell, or clinical aspect of cellulitis was observed (Figure 1). He had no medical history, no recent leg trauma, and no recent history of tonsil inflammation or respiratory symptoms. The first hypothesis was a pyogenic skin infection and outpatient care was started with oral pristinamycin after collection of superficial skin swabs.

Bacterial cultures yielded Gram-positive cocci identified as Staphylococcus aureus, and Grampositive rods identified as CD by MALDI-TOF mass spectrometry (MALDI Biotyper IVD Library 5627, Bruker). It was confirmed as CD by a home-made multiplex Polymerase Chain Reaction (PCR) test. The biovar was S. mitis. The tox gene detection performed by the national reference center by PCR was negative (tox -) [START_REF] Hauser | Polymerase chain reaction assay for diagnosis of potentially toxinogenic Corynebacterium diphtheriae strains: correlation with ADP-ribosylation activity assay[END_REF]. CD pharyngeal carriage investigation performed one week later was negative. Although skin lesions had improved after a one-week treatment, pristinamycin was switched to amoxicillin for an additional week to cover the recommended two-week treatment for CD infection. The outcome was favorable. The patient reported routine vaccination during childhood but missed the recommended booster set for the age of 25 in France [START_REF]Calendrier des vaccinations et recommandations vaccinales[END_REF]. A booster vaccination against diphtheria, tetanus, and pertussis (DTaP) was administered.

On that same month a 21-year-old Afghan male asylum seeker was admitted to the department of pneumology for suspicion of pulmonary tuberculosis. He had left his country eight months earlier and had been living in a refugee center in Paris for the last two months after a four-month stay in a center in Stockholm. His travel path included stays in Turkey, Greece, Austria, and Germany.

He had a three-month history of cough and thoracic pain, but no evidence of tonsil inflammation or respiratory distress. The CT scan performed at admittance showed cavitation of the upper right lobe. Airborne precautions were taken upon arrival. The sputum analysis was positive for Mycobacterium tuberculosis and CD. The bacterial isolate was confirmed as CD by PCR and MALDI-TOF and was determined to be of S. mitis biovar. It was shown to carry the tox gene (tox +) by PCR and to produce the diphtheria toxin by Elek test. The patient was diagnosed with tuberculosis and asymptomatic carriage of tox+ CD for which he received oral amoxicillin for two weeks. As his vaccination history was unknown, he was vaccinated with the DTaP vaccine but no serotherapy was administered. Six refugees from the same center but no social/health care worker were identified as individuals who had close contact with the patient and local authorities were notified. Among them, two persons had left before testing could be performed eight days later. Testing of the four other contacts was performed in our department, as there was no mobile team to visit the refugee center. All pharyngeal swabs were negative. Every contact person was given vaccination with DTaP (no vaccination history was available) and amoxicillin for 10 days.

Discussion

Diphtheria is once again a growing concern, and these two case reports illustrate some characteristics of this disease and the difficulties faced by physicians when managing such patients.

The two patients were a refugee and a traveler. Refugees are known to be at risk for diphtheria because of frequent incomplete vaccination and promiscuity in the refugee centers. More than 23 cases of cutaneous CD infection (including nine tox+ bacterial strains) have been reported among refugees within the last two years in Europe, mainly in Germany and Switzerland, but also in Denmark and Sweden [START_REF] Meinel | Outbreak investigation for toxigenic Corynebacterium diphtheriae wound infections in refugees from Northeast Africa and Syria in Switzerland and Germany by whole genome sequencing[END_REF][START_REF]Cutaneous diphtheria among recently arrived refugees and asylum seekers in the EU, 30 July 2015[END_REF][START_REF] Sane | Respiratory diphtheria in an asylum seeker from Afghanistan arriving to Finland via Sweden[END_REF]. In these studies and in our patient no information on the vaccination status was available. No data on the seropositivity for diphtheria in migrants is available in France but a study recently showed that 76.1% of refugees in Germany had no longterm protection against diphtheria [START_REF] Jablonka | Tetanus and diphtheria immunity in refugees in Europe in 2015[END_REF].

However, cases reported in returning travelers are uncommon. For instance, 16 cases of imported cutaneous CD infections were reported in Europe between 2000 and 2009 [START_REF] Wagner | Diphtheria in the postepidemic period, Europe, 2000-2009[END_REF]. However, only infections caused by tox+ CD are required to be notified, explaining the lack of data on tox-CD.

The use of the MALDI-TOF approach has improved the identification of CD, and may partly account for the recent increase in CD infection notifications [START_REF] Konrad | Matrixassisted laser desorption/ionisation time-of-flight (MALDI-TOF) mass spectrometry as a tool for rapid diagnosis of potentially toxigenic Corynebacterium species in the laboratory management of diphtheria-associated bacteria[END_REF]. Nevertheless, travelers to endemic regions are probably at risk of diphtheria especially elderly patients, because of a lack of booster immunization in elderly people. A European sero-epidemiological study reported that increasing age was related to an increase in seronegativity. The seroprotection against diphtheria dropped from about 75% in individuals aged <15 years to about 30% for those aged >60 years in the Czech Republic where the diphtheria vaccine coverage rate was >95% [START_REF] Di Giovine | Comparative seroepidemiology of diphtheria in six European countries and Israel[END_REF].

Cutaneous CD infection is difficult to diagnose. A chronic ulcer with a grey membrane is the characteristic presentation, but it frequently presents as a non-specific lesion, often appearing on a pre-existing lesion. Pyogenic bacteria may also be found as co-infection or colonization, making the diagnosis even more difficult, as observed in our patient. The patient acts as a reservoir of the bacterium and can cause outbreaks of the cutaneous and respiratory disease [START_REF] Luis | Clinical manifestations, diagnosis, and treatment of diphtheria 2016[END_REF][START_REF] Koopman | The role of cutaneous diphtheria infections in a diphtheria epidemic[END_REF]. In our case report, the partial improvement with pristinamycin treatment may indicate a coinfection with Staphylococcus aureus and not just colonization by the latter. This possible coinfection shows the importance of a routine use of the MALDI-TOF technique to detect a potentially associated CD infection.

The asymptomatic carriage is usually characterized by the incidental finding of a positive pharyngeal swab even though our patient had a positive sputum analysis. This condition is poorly described. A study conducted in 2007 screened 14,000 individuals in six European countries and observed four cases of symptomatic tox+ pharyngeal diphtheria, two cases of asymptomatic pharyngeal carriage of tox+ CD, and 26 cases of asymptomatic carriage of tox-CD [START_REF] Wagner | Screening for Corynebacterium diphtheriae and Corynebacterium ulcerans in patients with upper respiratory tract infections 2007-2008: a multicentre European study[END_REF].

The clinical management of a diphtheria case and contacts is challenging as few exhaustive official guidelines are available in European countries. Recommendations from National Public Health Authorities may also differ between countries about specific points [START_REF]de l'apparition d'un cas de diphtérie France: Haut conseil de la santé publique[END_REF][START_REF]Diphtheria: public health control and management in England and Wales England[END_REF]. For instance, French guidelines recommend treating a tox+ diphtheria carriage with antibiotics and serotherapy but English guidelines only recommend the antibiotic treatment. The World Health Organization (WHO) guidelines were edited in 1994 [START_REF]Manuel for Management and Control of Diphtheria in the European Region Copenhagen[END_REF].

The refugee patient was staying in a refugee center in Stockholm, Sweden, where cases of diphtheria were reported among refugees in 2015 [START_REF] Sane | Respiratory diphtheria in an asylum seeker from Afghanistan arriving to Finland via Sweden[END_REF]. However, our patient arrived in Sweden in 2016. Therefore either he was a chronic carrier, or he was recently contaminated with an unknown source in France. This stresses the importance of a close collaboration between countries to share information on the medical conditions of refugees, especially in cases like these where the cough due to pulmonary tuberculosis makes the spread of CD even more likely.

Our second case patient did not receive serotherapy although recommended by French guidelines [START_REF]de l'apparition d'un cas de diphtérie France: Haut conseil de la santé publique[END_REF]. Emergency treatment with diphtheria serotherapy seems difficult to deliver, making an international collaboration even more important. Colleagues from Finland had to use expiredalthough verified to be safeantitoxin for a patient presenting with toxigenic diphtheria [START_REF] Sane | Respiratory diphtheria in an asylum seeker from Afghanistan arriving to Finland via Sweden[END_REF].

Colleagues from Spain had no available antitoxin and had to get it delivered from a neighboring country. The management of an imported case patient was thus delayed [START_REF]A case of diphtheria in Spain[END_REF]. In 2009 at least 14 countries of the European Union had no national stock of diphtheria antitoxin [START_REF] Wagner | A review of the international issues surrounding the availability and demand for diphtheria antitoxin for therapeutic use[END_REF].

Life in refugee centers facilitates transmission and leads to outbreaks because of poor hygiene and promiscuity in individuals originally from countries or regions with no vaccination coverage.

Identified case patients need to be appropriately managed and prevention measures correctly implemented. Routine vaccination for every refugee should be standard along with a medical registration system for life-threatening diseases diagnosed in migrant people. The incomplete vaccine-induced immunity in elderly people is also a public health concern, and measures should be taken to ensure regular administrations of booster doses of vaccines. The rising number of unvaccinated or incompletely vaccinated individuals, increased lifespan, and growing vaccine skepticism could put Europe at risk of diphtheria, a vaccine-preventable disease.

Conclusion

CD is re-emerging. Insufficient vaccination coverage, international travel, and increasing migration call for overcoming the lack of knowledge among European physicians about diphtheria. Guidelines need to be updated, completed, and largely disseminated, including the management of asymptomatic CD carriers.

Figure 1 .

 1 Figure 1. Cutaneous diphtheria in a traveler coming back from Madagascar

Contribution of authors

All authors either contributed to the management of patients or to the literature search and writing of the article. All authors reviewed and approved the final version of the article.


Gert-Jan Duine

Thierry Hedde
email: thierry.hedde@cea.fr

Pierre Roubin

Pierre Durand

Marie Lothon

Fabienne Lohou

Patrick Augustin

Marc Fourmentin

Down-valley flows in complex terrain during KASCADE

Keywords: Nocturnal down-valley flows, stable boundary layer, distinct valley scales, field experiment KASCADE, intersecting valleys, complex terrain

The wind patterns in main and tributary valley systems undergo diurnal variation that are very dependent on the valleys dimensions. Characterising the down-valley winds is important for pollutant dispersion studies and a challenge for numerical simulations. This paper uses the KASCADE (KAtabatic winds and Stability over CAdarache for Dispersion of Effluents) observations to describe down-valley flow development and characteristics in two intersecting, shallow valleys of different size located in the prealpine region of southeastern France. The Durance Valley and its tributary valley, the Cadarache Valley, are different in depth (200 vs. 100 m depth, respectively), slope (0.2 vs. 1.2 ), width (5 to 8 vs. 1 to 2 km) and length (>50 vs. 6 km). Both down-valley flows were dominant features for the period of measurements. The distinct valley scales lead to different characteristics: the smaller Cadarache down-valley wind is primarily thermally driven and attains wind speeds up to 1 -4 m s -1 with its jet nose at around 30 m height a.g.l. The flow at that height initiates 2 hours after sensible heat flux turnover, a timescale that relates to linear wave solution for along-valley winds. The Durance down-valley wind involves larger scales -for both terrain features and large-scale circulation -, is related to a mountain-to-plain circulation and therefore its onset occurs with a delay of 6 hours after sunset. It reaches speeds of 4 -8 m s -1 and a depth of 175 -225 m a.g.l. Its highest occurrence is just after sunrise. In between the two valley winds, a shear layer is observed whose origin could be a combination of a Cadarache valley flow return current and a Durance down-valley wind stream redirected by orography.

Introduction

When contaminants are released in a complex terrain area, there is a need to document both local and regional flows, especially for stable conditions when vertical dilution is very weak. Characterising the stable boundary layer (SBL) is therefore of extreme importance, and is seen as one of the more serious atmospheric conditions for air quality (Zanetti 1990), especially in complex terrains which modify the flow and thus strongly influence dispersion (Bowen et al. 2000;Triantafyllou and Kassomenos 2002;Salmond and McKendry 2005;Largeron 2010). Cadarache, one of the research centers of the Commissariat à l'Energie Atomique et aux Energies Alternatives (CEA), is situated in the prealpine area of Southern Alps in France, where up-and down-valley flows are frequently formed under weak synoptic forcing and clear skies. For the purpose of impact assessment and risk management of the CEA facilities, the French regulation requires a constant monitoring of the dispersion of their emissions to the public and the environment. This necessitates characterising the local meteorology, which can be particularly challenging for the Cadarache area as it is on the intersection of two shallow valleys.

In valley systems, thermally-driven winds develop in stable conditions due to horizontal temperature gradients (Whiteman 2000). During the night, this is initiated by radiative surface cooling and consequently, the air is forced to flow downslope as a result of negative buoyancy (Manins and Sawford 1979;Mahrt 1982). The converging cold air reaching the valley bottom -also referred to as drainage -often results in a weak lifting motion, creating return currents at approximately the valley height (Oke 1987; Zardi and Whiteman 2013). In narrow valleys, the colder air stagnates at the valley bottom and a cold pool develops below a valley inversion (Clements et al. 2003). Down-valley flows often reveal themselves as low-level jets (LLJs - Kraus et al. 1985), their depth, strength, durations and onset times inherently depend on terrain geometries (McKee and O'Neal 1989; Atkinson 1995; Monti et al. 2002; Zardi and Whiteman 2013), which makes local field campaigns a necessity, but generalisation challenging. Besides their thermally-driven nature, valley flows can have multiple origins, some related to above-valley conditions. In this regard, Whiteman and Doran (1993) classified them into the fore-mentioned fully independent thermally-driven flow and three (quasi-) dependent flow types: downward momentum transport, forced channelling and pressure-driven channelling. Downward momentum transport plays a large role in unstable and neutral stratification in wider valleys (Whiteman and Doran 1993). Narrow valleys during unstable or neutral conditions experience the forced channeled flow (Weber and Kaufmann 1998) and pressure-driven channelling can be of importance in wide and shallow valleys under slightly to moderately stable conditions (Carrera et al. 2009). In the present study we therefore focus on the different scales of two confluent, relatively shallow valleys, to investigate their distinct response to stability, and their flow interaction. Two cross-oriented valleys are under study; the large Durance Valley (DV) and its tributary valley, the Cadarache Valley (CV). Their widths (5 to 8 vs. 1 to 2 km, respectively), depths (200 vs. 100 m), average slope angles (0.2 vs. 1.2 km) and lengths (67 vs. 6 km) are distinctly different, making them to behave differently under stratified conditions. We aim to quantify the down-valley flow characteristics by means of observations for these two connected valleys in terms of timing, height, strength and duration, their flow interactions, and to identify the main mechanisms causing the valley flows to develop. Observations come from the KASCADE (KAtabatic winds and Stability over CAdarache for Dispersion of Effluents) field measurement campaign, conducted in prealpine southeastern France during the winter season of 2013, which was characterised as having strong diurnal patterns during weak synoptic forcing.

The paper is organised as follows: First, the context of this study, a detailed description of the study area and details of the KASCADE-campaign are provided in Sect. 2. For general understanding, the observations of regional weather phenomena and stability, and valley flows in particular are described in Sect. 3.1. Characteristics of the valley flows under investigation are given in Sect. 3.2. Conclusions, including a conceptual framework, are given in Sect. 4.

KASCADE experiment

Background

Cadarache is a research centre of the CEA, located in southeastern France (see Fig. 1) at the intersection of the middle DV and the tributary CV immediately upstream of Clue de Mirabeau. The operation of several facilities there requires to control emissions of pollutants and to assess their impact. Clear skies and calm winds are regularly observed in the region (Wrathall 1985), leading to the onset of locally-driven winds in the complex area in and around Cadarache. Understanding the patterns of valley flows, their interactions, and cold-air pool formation in the tributaries confluent to the DV is therefore essential to evaluate the environmental impact of Cadarache's activities. The KASCADE field measurement campaign was conducted in the winter of 2013 and was designed to characterise the local SBL and thermally-driven winds to inform future numerical simulations of pollutant dispersion for environmental impact studies.

Site characteristics

The orography of the region of interest is presented in Fig. 1. The DV has been formed by the Durance river, a braided river flowing from the Southern Alps to the Rhône Valley. The valley is 5 to 8 km wide with an average depth of 200 m. The section of the DV between Sisteron and the Clue de Mirabeau has a length of 67 km with an along-valley axis mean slope angle of 0.2°. The aspect ratio, defined as the ratio of the valley depth to the valley width, is 0.04. The along-valley direction is about 30° North (Fig. 1a). The valley narrows from 5 km to 200 m at the Clue de Mirabeau. The CV is a SE-NW orientated tributary of the DV. The valley is 1 to 2 km wide, 6 km long, with an along-valley axis slope angle of 1.2°, a depth of ~100 m and an aspect ratio of 0.05. The sidewalls have an average slope of around 6°. The land use in the CV is a mixture of broad-leaved and mixed forests, grasslands, artificial surfaces and buildings. [Figure 1 about here.] Other orographic features of different scales (Fig. 1) may influence the local wind patterns. For instance, the Southern Alps, about 70 km northeast of Cadarache, reach heights of at least 1500 m MSL and attain heights of up to 3000 m at a distance of 140 km. In very calm synoptic conditions mountain systems of this scale develop drainage currents that can modify the flow more than 100 km downstream (Jiménez and Cuxart 2014). A sloping plateau with an angle of ~1° (Plateau de Valensole) lies in between the Southern Alps and Cadarache on the DV left bank. This plateau can induce the formation of a LLJ during stable situations and therefore can contribute to flow-channelling in the DV. Two east-west orientated mountain ridges are at moderate distance from the site: the Luberon and the Sainte Victoire. Both have maximum heights of 1000-1100 m. The narrow Clue de Mirabeau facilitates cold-air pool formation, especially during wintertime. Finally, situated within the Cadarache site boundaries, the Maladroit hill is part of the CV northeastern sidewall with an altitude of 394 m, being 110 m above the valley floor. The Maladroit is hypothesized to play an important role in the interaction of the CV and DV flows.

Various important local meteorological phenomena are typical for the region, and appeared to occasionally interfere with the valley winds during the KASCADE experiment. The Mistral is a cold and dry wind that flows down the Rhône Valley, and can occur throughout the year (Reiter 1971;Guenard et al. 2005). At the Cadarache site the Mistral comes from a northwesterly direction. Not all Mistral events, however, reach Cadarache. In the afternoon, due to solar heating and proximity to the Mediterranean Sea, land-sea breezes can develop (Cros et al. 2004;Mestayer et al. 2005). Sometimes the breezes are coupled with the Mistral wind (Bastin et al. 2005) and reach the Cadarache site where they approach from the west. Lastly, precipitation events usually come from the southeast (Margerit 2006) mainly during autumn and early spring months although they are sometimes observed in winter. They have been intensively investigated during the HYMEX campaigns (Drobinski et al. 2013). The topographically induced meteorological mesoscale phenomena (e.g. Mistral) favour clear sky. Wrathall (1985) reported more than 2500 hours of sunshine per year in the region.

Instrumentation and measurement strategy

To investigate the onset of the valley flows, their occurrence and interactions during stable conditions, the KASCADE campaign was conducted during the winter of 2013. Nocturnal atmospheric stability in the region occurs throughout the year, but the relatively long nights enhance SBL formation. The data were collected between 13 December 2012 and 18 March 2013. 23 special intensive observation periods (IOPs) were carried out between 15 January and 2 March 2013. An IOP was planned when clear skies and weak synoptic forcing were forecasted, but were actually executed when the wind was weak enough for the tethered balloon to be operated, making these periods inherently related to weak wind conditions. The dates of IOPs are presented in Table 1.

[Table 1

about here.]

Measurements were made at different locations (see Fig. 1). A detailed description of the instrumentation is provided in Table 2. Cadarache has two permanent weather stations, La Grande Bastide (GBA) and La Verrerie (VER), which measure temperature at 2 m a.g.l. and wind speed and direction at 110 m (GBA) and 15m (VER) a.g.l. The top of the GBA-tower lies above the CV sidewalls. The M30 site (see Fig. 1) has a 30-m flux tower with sonic anemometers at 2, 10 and 30 m, one fast hygrometer at 30 m, longwave and shortwave radiation sensors at 2 and 20 m and two thermohygrometers at 2 and 30 m. This site is located in a flat open section with a 600-m fetch. The land cover is mostly grass with patches of bare soil. Turbulent fluxes of sensible heat, latent heat and momentum were During the IOPs, a tethered balloon was deployed and radiosondes were launched at the M30 site. The tethered balloon was equipped with up to three Väisälä probes separated by a distance of 50 m, measuring temperature, relative humidity, wind speed, wind direction and pressure. The tethered balloon was deployed to make measurements during the sunset and sunrise transitions, but in the last week of the campaign full-night measurements were also made. In general, profiles reached up to 300 m a.g.l. with an average ascent/descent speed of 0.3 m s -1 , yielding profiles every 15 min. At the surface an extra ground level probe was installed to give a pressure reference for height correction. To document the troposphere above 500 m a.g.l., radiosondes were launched during the IOPs, reaching an average altitude of around 5 km. Normally, the releases were at 12, 18 and 06 UTC and during the last campaign week when full-night measurements were made, radiosondes were also launched at 00 UTC.

All sensors were checked for inconsistencies and recalibrated during a two-month inter-comparison experiment right after KASCADE at Centre de Recherches Atmosphériques, Lannemezan, France. Detailed information on the correction procedures can be found in a KASCADE technical report (Appendix C in Duine 2015). Prior to the flux calculations, an additional despiking procedure was developed for the two lower sonic anemometers. The calculated fluxes were evaluated for both normal and de-spiked signals during the inter-comparison with another eddy-covariance package (Baghi et al. 2012) and no discrepancies were found for the de-spiked signal.

Results and discussion

General overview of the observations

The prevalent meteorological conditions

Data for the last 24 days of the campaign, a period during which 15 IOPs were conducted, are presented in Fig. 2. A diurnal temperature range of 15°C was observed regularly during the full campaign (40.9% out of the total of 93 days), with a peak of 20°C on 17 February (see Fig. 2a). For midlatitude regions in wintertime these values are rather extreme (Wallace and Hobbs 2006). They can be explained by the distance of the site to the Mediterranean Sea, generally clear skies and low soil moisture (Jackson and Forster 2010). Furthermore, because of the complex orography, cold air tends to drain from higher slopes to the valley bottom creating strong inversions and increasing the diurnal temperature range (Geerts 2003). A temperature inversion developed during clear-sky nights, as shown by positive temperature differences (T110m -T2m), and can reach 5°C (see Figs. 2a and2c).

Westerly winds prevailed during the afternoon at both the GBA and M30 sites (see Fig. 2b). Once the sun set, a SBL formed rapidly and at the M30 site (in the CV) the observed wind turned to the SE, i.e. down-valley and opposite to the daytime direction. This weak drainage flow will, hereafter, be named the Cadarache down-valley (CDV) wind, and is primarily thermally-driven (Duine et al. 2016). Whenever strong Mistral events occurred (e.g. 8th February, see Fig. 2b), no CDV wind was observed. During cloudy (and rainy) periods (e.g. the nights from 10 to 11 and 17 to 18 February) SE winds often occurred, as can be seen from the wind measurements at 110 m a.g.l. They coincidentally have the same direction as the CDV wind. At 110 m a.g.l. at the GBA site, i.e. well above the CV, the daytime westerly flow often persisted until the early evening and then turned progressively to the NE, i.e. along the DV. We will call this flow the DDV wind hereafter. We will elaborate on characteristics (timing, occurrence, strength) of both the DDV and CDV winds in Sect. 3.2. Figures 3b, 3c and3d show the probability density functions and the median values (for the whole campaign) of the sensible heat flux H, stability index z=L and friction velocity u*, respectively. On average a stable stratification forms 1.5 hour before sunset (see Figs. 3a, 3b and 3c), and is destroyed shortly after sunrise (not shown). In this shallow, narrow CV the sidewalls are able to create shading effects leading to negative values for H before sunset (see Nadeau et al. 2013). H reaches its minimum value half an hour after sunset, then increases gradually to 0, similar to what is observed during the CASES-99 experiment (Steeneveld et al. 2006), and becomes positive 1 hour after sunrise (not shown). In general, a clear transition from unstable to stable conditions was observed. Also the distribution and median value of u* clearly illustrate the transition from unstable, turbulent to stable conditions (see Fig. 3d). The wind speed decreases and weak turbulent mixing is observed during the night, with u_ rarely exceeding 0.2 m s -1 and being often lower than 0.1 m s - 1 . Turbulent mixing increases again from 1.5 hours after sunrise.

[Figure 2 about here.]

The prevalent winds at the VER site are shown by means of the frequency of the wind direction, and are derived from the Sodar measurements (see Fig. 4). Four main wind patterns can be distinguished: the NNE DDV wind, the SE wind that often brings clouds and precipitation, the Mistral from the N to NW sector, and afternoon and evening WSW winds.Westerly winds observed before sunset and just after sunrise are attributed to a mixture of Mistral and general westerly circulation (see Figs. 4a and4b). Southeasterly winds often associated with cloudy and rainy events are detected before sunset. The DDV flow is fully developed only 6 to 9 hours after sunset (see Fig. 4d), and persists at least 3 to 6 hours after sunrise (see Fig. 4h). We elaborate more on the DDV wind in the next sections.

[Figure 3 about here.]

Relationships between within-and above-valley winds

The observations collected during the KASCADE field campaign were very localised, but nevertheless can be used to infer the mechanisms driving the within-valley winds. In the following, we use the diagnostic framework introduced by Whiteman and Doran (1993), which relates within-valley to above-valley winds. The relationships for four possible forcing mechanisms (thermally-driven, downward momentum transport, forced chanelling, and pressure-driven chanelling) are presented in Fig. 5a-d.Within-and above-valley wind direction data, averaged over 2 hour periods, are shown in Fig. 5e for the CV and Fig. 5f for the DV. To estimate the within-and above-valley wind directions for the CV we used wind direction data at 30 m a.g.l. at the M30 site, and at 110 m a.g.l. at the GBA site, respectively. For the DV, the wind directions at 175 m and 475 m a.g.l. derived from the Sodar measurements at VER were used to estimate within-and above-valley wind direction, respectively.

The signature of thermally-driven up-and down-valley winds is clear for the CV (see Fig. 5e), for any above-valley wind direction. Thermally-driven up-valley winds are more variable in direction. Downward momentum transport does appear to be much less observed within the CV and DV (see Fig. 5f) and is principally associated with Mistral (NW) and cloudy (SE) events. Pressure-driven chanelling is driven by large-scale pressure gradients and does not show a diurnal cycle (Whiteman and Doran 1993), hence it can be excluded as a key driving mechanism in the CV. Besides, this mechanism is mostly found in valleys with moderate and cloudy climates, while the Cadarache region mainly experiences clear skies. From this figure, forced chanelling cannot be fully excluded as a mechanism, but the CV flows are primarily thermally-driven (Duine et al. 2016).

A clear diurnal pattern is also found for the DV wind (see Fig. 5f). Note that relationships between within-and abovevalley winds for this valley are more complex because the VER site is off the axis of the valley (see Fig. 1). Besides, the relationships associated with the forced and pressure-driven chanelling mechanisms overlap (see Fig. 5d). Similarly to the CV, no signature for pressure-driven chanelling is observed within the DV. The DDV winds present a signature of forced chanelling as they are associated with above-valley winds with directions in the range [300°-120°] (clockwise), and occur mostly during the late night and early morning hours. Durance up-valley winds are associated with abovevalley wind directions in the range [210°-360°] (clockwise) during the afternoon and early evening. Considering the position of the Southern Alps northeast of the VER site, above-valley flows can result from a mountain-to-plain circulation (Zardi and Whiteman 2013). Therefore, one of the mechanisms driving a DDV wind is drainage from the large-scale circulation. This was investigated and confirmed by mesoscale modelling covering a larger area [see Chap. 7 in Duine (2015)]. It was also found that, in addition to the drainage associated with the mountain-to-plain circulation, drainage from tributary valleys (Vallée de la Bléone and Vall´ee de l'Asse) plays an important role in the development of the DDV flow. Hence there are several possible drivers of the DDV flow: local valley drainage, chanelling effects from larger-scale synoptic winds, LLJs originating from the Plateau de Valensole (see Sect. 2.2), and drainage from the Southern Alps.

[Figure 4 about here.] Both valleys show features related to down-valley flows, which were dominantly present. By picking up individual cases, we will try to infer under which conditions the down-valley winds develop, and their characteristics.

Valley-flow characteristics

This section focuses on the description of the flow characteristics, by exploring a case study under weak synoptic forcing (Sect. 3.2.1), which leads to timing characteristics (Sect. 3.2.2). Then, we investigate the flows under moderate synoptic forcing (Sect. 3.2.3), which allows us to draw conclusions on depth, strength, and occurrence of the valley flows in Sect. 3.2.4.

Development under weak synoptic forcing: IOP 21

We selected IOP 21 to investigate the characteristics of the flows in both valleys under weak synoptic forcing and clear skies. The IOP started at 12:00 UTC on 27 February 2013. Astronomical sunset and sunrise were at 17:21 and 06:16 UTC, respectively. Four days before the IOP started, a precipitation event covered the region with a 10-cm deep snow layer (see Fig. 2a). The snow had mostly melted by the start of the IOP, as revealed by the values of SW↓ (see Fig. 2c). The synoptic situation was characterised by a high pressure area over the United Kingdom that extended eastward and a low pressure area over southern Spain, driving an easterly flow over southern France.

Data from the radiosonde launched at 12:00 UTC reveals a dry mixed layer extending from the surface to 1400 m a.g.l.

(see Fig. 6a). The time evolution of the mixing ratio profile indicates moisture advection during the IOP, with dew formation in the early morning as a result of the strong stability (not shown). A well-marked westerly jet was observed below 900 m at 18:00 UTC (see Fig. 6b). Observations from the Sodar at VER (not shown), and the GBA and M30 stations (see Fig. 2b) indicate that this westerly wind at Cadarache formed around 17:00 UTC and lasted for several hours. By the end of the night (06:00 UTC), the depth of the SBL reached 400 m, which is beyond the height range of the tethersondes.

[Figure 5 about here.]

We selected 14 among the 66 tethersonde profiles obtained during IOP 21 (see Fig. 7). Potential temperature profiles indicate that at 17:13 UTC (that is shortly before sunset) a stable layer had formed in the first 30 m above the surface. By 18:43 UTC the atmosphere in the first 300mabove the surface had cooled by 2°C, and the near-surface temperature inversion strength had reached ~85 K km -1 , see Fig. 7a. The SBL deepened steadily throughout the night while the wind profile changed dramatically within the first hours. The weak winds aloft between 21:40 and 22:31 UTC allowed the CDV flow to develop throughout the whole valley depth, that is ~100 m (see (Fig. 7b). The CDV wind takes the shape of a jet, with a maximum wind speed of about 2.5 m s -1 . The wind at around 175 m a.g.l. turned to the DDV wind direction. At this height the turbulence collapsed (i.e. σw <0.1 m s -1 ) at 23:15 UTC (not shown). At 23:52 UTC, i.e. 6.5 hours after sunset, the DDV wind had developed, and presented a jet with a maximum speed of 3 m s -1 at 175 m a.g.l. The timing for both the CDV and DDV winds is further discussed in Sect. 3.2.2.

[Figure 6 about here.]

During the second part of the night, the surface continued to cool until a minimum temperature of -3°C was reached before sunrise (06:16 UTC), see Fig. 7c. This lead to a diurnal temperature range of 17.4°C (see Fig. 2). The DDV wind has attained speeds up to 4 m s -1 (see Fig. 7d). After 04:00 UTC a shear layer is observed between 50 and 100 m: the wind rotates clockwise with altitude from SE to NNE, through a 265° rotation instead of 95°. This typical behavior was observed during other IOPs featuring both DDV and CDV winds. The following mechanisms may explain these observations: (i) a return current [or closing branch belonging to the thermally driven circulation (Zardi and Whiteman 2013)] flowing just above the CDV flow, (ii) an intrusion of the DDV flow into the CV, partially redirected by the Maladroit hill, or (iii) a combination of (i) and (ii). At the surface, the CV wind reversed from down-valley to up-valley at 07:50 UTC (not shown). The DDV wind lasted until around 10:00 UTC as revealed by the Sodar measurements (not shown). Later on, the SE-direction of the synoptic wind was re-established.

IOP 21 is one of the many examples of calm synoptic forcing conditions. Out of the total of 23 IOPs, the CDV and DDV winds co-existed 14 times (i.e. 61%). This number is restricted by the appearance of a DDV wind, as the existence of the CDV wind is less susceptible to large-scale disturbances. In the next section we will draw some conclusions about the typical timing of both flows.

Timing of the CDV and DDV flows

Figure 8 shows wind direction during the KASCADE campaign as a frequency distribution on a time scale referenced to sunset at 2, 10 and 30 m a.g.l. for the M30 site and 110 m a.g.l. for the GBA site. The most frequent wind direction at all heights in the CV before sunset is WNW, see also Figs. 1b, 4, and7. It is not necessarily a local CV effect as it is also present at 110 m. The wind turns to SE at 2 and 10 m a.g.l. at sunset, indicating the onset of the CDV wind. One hour after sunset the 30 m wind turns to SE as well, as indicated by the median value on Fig. 8d. The CDV flow persists at all heights until 1 hour after sunrise and reverses to up-valley 2 hours after sunrise (see Fig. 9). The transition to up-valley occurs earlier at 30 m a.g.l. than at 2 m a.g.l. This feature was observed regularly from the tethersonde measurements during the morning transition (e.g. IOP 21 in Sect. 3.2.1), and may be explained by the valley shading effect leading to earlier surface heating at the CV floor in the lower part of the valley than in the upper part of the valley. The CDV wind is likely to be eroded from the top, once the WNW flow aloft is re-established. During nighttime (see Fig. 8), a wind from the NW is sometimes observed at 2, 10 and 30 m a.g.l.: this is the signature of Mistral events. These events were observed episodically during the entire campaign.

[Figure 7 The timing of the CDV wind to reach its steady state can be related to the timescale Tl associated with the linear wave solution of the along-valley wind equations derived by Egger (1990), further investigated by Schmidli and Rotunno (2015), and related to the valley length Lx by Arduini et al. (2016). A typical timescale Tss for the along-valley wind component to reach its (quasi-) steady state is found to be about 3 Tl, whereTl = 2Lx=(hN0), with Lx = 6000 m, h = 100 m the valley depth, and N0 proportional to the Brunt-Väisälä frequency (N0 =sqrt ( (g/θv) ∂θv / ∂z) ). The CDV wind maximum is typically around 30 m a.g.l. (see Sect. 3.2.4). Using a characteristic temperature gradient at sunset between 30 and 2 m (see Fig. 3a), we get N0 = 0.05 s -1 , yielding Tl = 40 min. With the coefficient of proportionality of 3, the CDV wind would reach its steady state in 120 min. This coincides well with observations when we consider the time when H changes sign at 30 m a.g.l. as an initial time for the CDV wind to develop (approximately 1 hr before sunset, see Fig. 3b), indicating that the CDV flow reaches its steady state 1 hr after sunset (see Fig. 8d).

The timing of the onset of the DDV flow can be inferred from the wind measurements at 110 m a.g.l. at GBA (see Fig. 8d). Since the tower is situated above the CV (see Fig. 1) it captures the DV flow. After sunset, the daytimeWNWwind continues steadily into the night. The first appearance of the DDV flow is visible about 4.5 hours after sunset when the wind shifts to the DDV direction (NNE direction). Winds from the WNW and NNE directions coexist as the night progresses. The highest occurrence of the DDV wind is found just after sunrise (see Fig. 9). After sunrise the DDV wind persists for longer than the CDV wind, shifting back to the daytime WNW direction around 4.5 hours after sunrise.

To investigate further the time delay in the onset of the DDV flow, we consider the Sodar measurements at VER (see Fig. 1). Figure 10 shows the onset of a wind oriented in the 10°-50° direction range, i.e. along the DDV direction at a height of 225 m a.g.l. The DDV wind is most noticeable between 100 and 300 m a.g.l. (see Fig. 4). We restrict our analysis to winds lasting for at least 1.5 hours to exclude coincidental occurrences in the same wind direction range during low wind speed conditions. The figure shows that DDV winds can start throughout the night at this height, but the highest probability for the onset is between 6 and 9 hours after sunset, although the onset at lower heights could occur earlier.

[Figure 9 about here.]

The turbulence, as shown by the median of σw at 175 m a.g.l., generally decreases to its lowest level 6 hours after sunset (see Fig. 11a). Considering IOP-days only (see Fig. 11b) leads to a similar result. A possible explanation for the development of the DDV wind when turbulence decays is that the DDV wind is not aligned with the wind above the DV, and as long as turbulence maintains a certain level, vertical mixing of momentum prevents the DDV wind from developing in the valley. Figure 11 highlights the strong diurnal patterns observed for both the full campaign and during IOPs.

[Figure 10 about here.]

We have not attempted to find a typical timing for the DDV wind to develop, as was done for the CDV wind, for multiple reasons: N0 cannot be determined due to the lack of stability information within the DV, the theory of the linear wave solution might be less appropriate for large-scale valleys (Schmidli and Rotunno 2015), and the DV is far from ideal with its depth varying greatly along the 70 km-long fetch. We recall however that the distance from the measurement site to the Southern Alps is around 70 km, and that on these scales (both distance and mountain range height) a mountain-toplain circulation might very well explain the found time delay: Jiménez and Cuxart (2014) reported that the Pyrenees mountain ridge can affect flows on horizontal plains more than 100 km away, with a delay of more than 6 hours. This delay is fairly consistent with the findings of the present study. If we consider the DDV flow to be a drainage flow restricted to the valley itself, we come to a similar conclusion: in such case the DDV flow could form around Sisteron (see Fig. 1a), and with a typical wind speed of around 4 -8 m s -1 (see Sect. 3.2.4), it would take around 5 hours to reach the Cadarache area.

To summarise, the CDV and DDV flows have different timings relative to sunset and sunrise. We also note a strong difference in turbulence decay for the two valleys (e.g. compare Figs. 3 and11). Wind reversals in the smaller CV occur close to sunset and sunrise transitions. Onset occurs when a stable layer forms at the surface and cessation occurs shortly after sunrise when atmospheric stability is reduced. The onset and cessation of the DDV flow occur with a delay related to sunset and sunrise transitions. This flow is affected by larger scales; the DV has a longer fetch and the VER site is on a sidewall of the DV at its lower end so that this transition takes time to reach the site's elevation. Besides, the DV is much wider than the CV and is therefore more exposed to external influences. To investigate the importance of such external influences on the differences in the valley flow developments, we pick up a case study under moderate synoptic forcing in the next section.

Development under moderate synoptic forcing: IOP 5

In this section we examine the non-typical IOP 5, which started on 24 January 2013 at 12:00 UTC. Astronomical sunset was at 16:36 UTC and sunrise at 07:01 UTC.

During IOP 5 a low pressure system was present above Italy which resulted in a northerly wind in the region of interest. This IOP was conducted under moderate Mistral conditions during which turbulence was sustained (the average value of σw between sunset and sunrise was 0.19 m s -1 at 125 m a.g.l.) and no clouds were observed. Note from Fig. 11 that the value for σw is indeed higher than the average for the IOPs. At 12:00 UTC a convective boundary layer extended from the surface to ~1000 m height (see Fig. 12a). The wind on the profile is weak and northerly, except in the lowest 300 m a.g.l., where NW winds were observed (see Fig. 12b). At 18:00 UTC, the layer below 800 m has already cooled. Two inversion layers are observed with a neutral layer in between, the lowest one with its top at around 100 m and the upper one at 800 m. Given its height and the short time (1.5 hr) after sunset, the upper inversion cannot be a consequence of the local surface cooling, but has to be related to cold air advection by the Mistral wind. By the early morning (06:00 UTC), the lowest stable layer has deepened to 250 m a.g.l., with a sustained NW flow; the upper inversion has also strengthened.

[Figure 11 about here.]

The tethersonde profiles confirm the moderate Mistral conditions. Among the 50 tethersonde profiles, 7 in the evening and 7 in the morning transition were selected (see Fig. 13). The wind below 200 m turned from NW to W at 14:00 UTC (see Fig. 13a), as revealed by the Sodar at VER and M30 observations (not shown). A similar shift was observed during IOP 21, and is a typical feature for the region as shown in Sect. 3.1.1. Surface cooling started shortly before sunset, together with a weak and variable near-surface wind (16:48 UTC, see Fig. 12b). The CDV wind was observed with a jet maximum at 30 m a.g.l. ~3 h after sunset (19:31 UTC), that is later than the typical onset timing during weak synoptic forcing (see Sect. 3.2.2). The delay is related to the stronger synoptic forcing, which sustained mixing in the lower layers and so delayed the onset of the CDV wind. By the early morning the SBL depth was about 150 m a.g.l. (see Fig. 13c). The Sodar observations reveal that the northerly wind was observed during the full-night period, with speeds of around 4 -6 m s -1 up to 500 m a.g.l. (not shown). No clear DDV wind (30°) was observed below 200 m (see Fig. 13d), unlike during IOP 21. A steady CDV wind persisted near the surface, although it was thinner (~50 m), weaker (1 -2 m s -1 ) and therefore more variable than during IOP 21. The Mistral above the CDV flow prevented its development by eroding its top. We conclude that the depth of the CDV flow is influenced by the wind above the valley [see also Barr and Orgill (1989)]. In the next section we determine typical depths and wind speeds of the CDV and DDV flows.

[Figure 12 about here.] 3.2.4. CDV and DDV wind speed, depth and occurrence

In this section we focus on the data from the IOPs characterized by weak synoptic winds, pointing out typical flow depths and interactions, as well as influences of local features. We have disregarded IOP 2, since a very strong Mistral event prevented tethersonde measurements, and we restrict our analysis to nighttime periods. Figure 14 shows the wind roses at various heights, derived from the tethersonde wind data. The highest frequency of CDV winds occur at around 33 m a.g.l. At 51 m there is still a high frequency of down-valley winds, along with a lower frequency of up-valley winds. At higher elevations (e.g. 75 m) the down-valley wind is barely present and at 111 m it is absent. Typical values of the wind speed for the down-valley jet are 1 -4 m s -1 .

[Figure 13 The influence of clear skies and relatively weak synoptic forcing on the diurnal wind pattern can be demonstrated by comparing Sodar observations at VER for the full winter period (see Fig. 4) with Sodar observations at the same site for the IOPs only (see Fig. 15). The largest differences are found in the SE and NW sectors. Both sectors are characterised by winds that are principally synoptically driven: we recall that the SE large-scale winds are associated mostly with cloudy and rainy events, and the strong NW winds are associated mostly with Mistral events. These features are not visible for the weak synoptic conditions of the IOPs. However, SE winds above 100 m were observed during some IOPs (primarily IOP 15; see Fig. 2b). It is clear that for weak synoptic forcing during the afternoon, westerly winds are predominant. Figure 16 shows tethersonde data for the IOPs at the M30 site. By comparing Figs. 15a and16a, one observes a difference around sunset between the measurements inside and outside the CV. The Sodar observations are divided into two groups: WSW and NW, whereas the tethersonde gather along a single WNW direction. Just after sunset inside the valley, the CDV flow has developed near the surface (see Fig. 16b) and remains well established throughout the night (see Fig. 16d-g). The onset of the CDV flow is consistent with the tower measurements in Sect. 3.2.2 and the data shown in Fig. 8. Figures 16b-c and 15bc further imply that the westerly flow remains until 6 hours after sunset, when the DV flow shifts towards the DDV direction. This confirms the uncertainty on the onset timing of the DDV flow mentioned in Sect. 3.2.2. During the 6 hours preceding sunrise and the 3 hours after sunrise (Figs. 16d-g and 15d-g), the DDV flow is fully developed between 100 and 300 m a.g.l., with wind speed in the range 4 -8 m s -1 (see Fig. 14). LLJs were frequently observed during the IOPs with a maximum wind speed in the range 175 -225 m a.g.l. (not shown).

At the end of the night during IOP 21 a shear layer was observed between 50 and 100 m (see Sect. 3.2.1). This layer, between the CDV and DDV winds, was observed in other IOPs and therefore appears in Fig. 16f-g. In the shear layer the wind turns clockwise with elevation by 265°. It is typically present in the late night and believed to represent large horizontal motions generated in the height range between the CDV and DDV flows when stratification is at its strongest. It is readily identifiable in the IOP mean wind data and is especially apparent during IOPs when both the CDV and DDV flows co-exist: out of 17 early-morning periods when tethersonde data are available, this situation occurred 11 times, and a shear layer was observed 10 times.

The Maladroit hill (see Fig. 1), which lies north of the tethersonde site, may create a sheltering effect that influences the CDV flow. During nighttime at the GBA site, the northeasterly DDV flow is observed at 110 m, whereas at the M30 site it is much less frequent at this height and completely absent below.

[Figure 16 about here.]

Conclusions

Field measurements from the KASCADE field campaign held during the winter 2012/2013 were analysed to characterize the atmospheric stability and valley winds that develop in the connecting Cadarache valley (CV) and Durance valley (DV). The experiment collected high quality data for three months, including 23 periods of intensive observations focused on anticylonic conditions, which are typical for the region. The representative case studies in relation to the full campaign data implied the importance of valley flow characteristics in such shallow (around 100 to 200 m in depth) valleys, and showed their significance for atmospheric dispersion studies.

The distinct dimensions of these two valleys lead to different flow characteristics, and are conceptualised in Fig. 17. The along-valley wind system in the CV is generally thermally-driven, that is up-valley during the day and down-valley during the night. The Cadarache down-valley (CDV) flow develops shortly after sunset. Its onset is immediate, and reaches a quasi-steady state around 2 hrs after the sensible heat flux changes from positive to negative. This is consistent with linear wave theory. Once established, the CDV wind persists throughout the night and ceases within two hours after sunrise. Due to sheltering effect, the CDV flow can also develop under moderate synoptic forcing. During very weak synoptic wind conditions, it extends throughout the valley depth whereas under stronger synoptic forcing it typically reaches a few tens of metres a.g.l. only, as stronger external winds reduce its depth by erosion from the top. The CDV wind speed is typically in the range 1 -4 m s -1 .

[Figure 16 about here.] The Durance down-valley (DDV) flow experienced at Cadarache develops much later than the CDV flow. Although the onset timing is variable, the average delay of 6 hours with respect to sunset could be related to terrain geometry on the mesoscale. The development of the DDV flow is associated with a reduction in turbulence level and an increase of the wind speed towards sunrise, and the maximum wind speeds can even occur after sunrise. Once established, the DDV wind is typically in the range of 4 -8 m s -1 at 175 -225 m a.g.l. The DDV wind generally ceases 3 hours after sunrise.

The analysis revealed an interaction between the DDV and CDV flows at the Cadarache site. The governing processes of this shear layer are hypothesised as being a return current or an intrusion of the DDV flow into the CV. Such knowledge however is crucial for atmospheric dispersion studies. To further increase our understanding of the wind regimes in the Cadarache region during stable conditions, numerical modelling and dispersion studies are currently in progress. Recently, the Weather Research and Forecasting numerical model was tested and evaluated using the field campaign data (Kalverla et al. 2016), showing the potential to be used for model-based dispersion studies (Duine 2015). 1

  Thermally-driven valley flows were investigated extensively during the VTMX (Doran et al. 2002), ASCOT (Clements et al. 1989), METCRAX (Whiteman et al. 2008), and more recently MATERHORN (Fernando et al. 2015) and TRANSFLEX (Fernando et al. 2013) campaigns. Although the literature on valley flows is vast and the governing processes are well described [e.g. Whiteman (1990)], their interaction with tributary flows has received little attention. Complex interactions between tributary and valley wind systems however exist (Porch et al. 1989), as the tributary flow can considerably contribute momentum to the main valley flow (Coulter et al. 1989), can modify the main valley flow (Porch et al. 1989; Zängl and Vogt 2006), and particularly can contribute to clean polluted air in the main valley (Banta et al. 1997). Coulter et al. (1991) compared three different tributary flows, finding that small differences in topography of the tributaries have meaningful effects on the main valley flow. Sakiyama (1990) related the geometry of tributary and main valley flows to their distinct sizes. Despite the effort of these studies, it is clear that different terrain geometries result in different wind characteristics, especially in complex terrains with tributary valley systems.

  Kalthoff et al. (2005) showed the existence of an along-valley wind in the DV only in the early morning during the ESCOMPTE campaign (Cros et al. 2004; Mestayer et al. 2005) and Quenol and Beltrando (2008) investigated the effect of spring frost on the valley wind system in the lower part of the DV.

  calculated using the open-source flux calculation package EddyPro (LI-COR Biosciences, USA), version 4.1.2. The radiometers at two different heights (the lower one on a nearby mast) provided measurements of radiative flux divergence, a crucial variable controlling SBL formation that is still relatively unexplored (Sun et al. 2003; Savijärvi 2006; Hoch et al. 2007; Steeneveld et al. 2010). At VER a Remtech PA2 Sodar provided averaged wind speed and direction, and turbulence characteristics every 15 min with a vertical resolution of 25 m.

  about here.] [Figure 8 about here.]

  about here.] [Figure 14 about here.] [Figure 15 about here.]

Experimental data is available at the data website http://kascade.sedoo.fr/ or by contacting Dr. Thierry Hedde (Thierry.Hedde@cea.fr)

Acknowledgements

The authors gratefully thank S. Derrien, E. Pique, C. Jarnot, P. Medina, C. Darbieu, G. Bret from the Laboratoire d'Aérologie for their invaluable support on the technical and observational work, as well as F. Jourdain, Y. Margerit, S. Bontemps and E. Cohenny from Laboratoire de Modélisation des Transferts dans l'Environnement for their assiduous operation during the the KASCADE-experiment. H. Delbarre from Laboratoire de Physico-Chimie de l'Atmosph`ere is acknowledged for his effective support in the provision of the Sodar equipment. Gert-Jan Steeneveld is thanked for his useful suggestions on the manuscript. The KASCADE campaign and the PhD grant were funded by the CEA.

List of figures

1 The study area. a) provides a zoom on the area which is indicated by the star in the lower left frame, b) is the enlargement of the green rectangle in a). The red lines indicate the Durance (DV) and the Cadarache (CV) valleys with their downvalley directions marked in degrees from North (in red). Measurement locations are indicated by black dots. The most important regional orographic features and their heights MSL are also indicated. The terrain elevation scale applies on all frames. Source maps: Geoportail.gouv.fr 2 Time series during the 24 last days of KASCADE of a) temperature at three different heights and precipitation, b) wind speed and direction (lines and dots, respectively) at two heights, and c) radiation components (SW stands for shortwave, LW for long-wave, in and out for incoming and outgoing, respectively, and Q-net is net-radiation). Precipitation was measured with a rain gauge at VER, radiation close to M30 at 1.2 m height, other measurement locations and heights are indicated in the corresponding boxes. The gray shaded areas indicate nighttime periods. In a) the IOPs are indicated by their number, the orange dashed lines designating the start and end time of each IOP. In b) on the left axis CDV and DDV indicate the down-valley wind direction for CV and DV, respectively. All quantities shown are 10 min averages, except for precipitation which is a 10 min accumulation and winds from M30 which are 30 min averages. See Fig. 1 for the location of the measurements.

List of tables

1 Summary of KASCADE IOPs with the number of radiosondes (RS) launches and tethered balloon (TB) profiles indicated. IOPs start at 12 UTC and last 24 hours. Full-night measurements were made on IOPs 19-23.

2 Summary of measurements at the sites shown in Fig. 1b and their platforms: geographical coordinates, base elevation above mean sea level (MSL), measurement height above ground level (a.g.l.), instruments used, measured parameters, sampling frequency. Abbreviations include T: temperature, RH: relative humidity, U', V', W': sonic wind components, T': sonic temperature, H2O': humidity, CO2': carbon dioxide concentration, SW↓; SW↑; LW↓; LW↑: incoming and outgoing shortwave & longwave radiation, PP: atmospheric pressure, FF: wind speed, DD: wind direction, Td: dew point temperature, vz: ascending/descending speed, see text for site abbreviations. # denotes the number of measurements at different heights.