Julien Brunel
email: julien.brunel@onera.fr

Peter Feiler

Jérôme Hugues
email: jerome.hugues@isae-supaero.fr

Bruce Lewis
email: bruce.a.lewis.civ@mail.mil

Tatiana Prosvirnova
email: tatiana.prosvirnova@irt-saintexupery.com

Christel Seguin
email: christel.seguin@onera.fr

Lutz Wrage
email: lwrage@sei.cmu.edu

Performing Safety Analyses with AADL and AltaRica

Keywords:

AADL and AltaRica languages can be used to support the safety assessments of system architectures. These languages were defined with different concerns and this paper aims at presenting their principles and how they can be related. A translator from AADL to AltaRica is proposed and its prototype is applied to a simplified flight control system of a UAV. The resulting AltaRica model has been analyzed with the AltaRica safety tools and the experimental results are discussed.

Introduction

The interest of industrial community in Model-Based System Engineering (MBSE) and Model-Based Safety Assessment (MBSA) is gradually increasing. In this paper we consider two modeling languages: AADL and AltaRica.

AADL (Architecture Analysis and Design Language) is a multi-concerns modeling language dedicated to distributed real-time embedded systems [START_REF]SAE: Architecture Analysis and Design Language (AADL) AS-5506[END_REF]. It proposes several annexes to describe embedded systems behavior. The AADL Error Model V2 (EMV2) [START_REF] Delange | Architecture fault modeling with the aadl error-model annex[END_REF] is an error annex focused on Safety Analyses. It offers a terminology and an ontology to capture key features of failure/error propagations.

AltaRica [START_REF] Arnold | The altarica language and its semantics[END_REF][START_REF] Prosvirnova | The altarica 3.0 project for model-based safety assessment[END_REF] is a high level modeling language dedicated to Safety Analyses. Its formal semantics allowed the development of a set of efficient assessment tools, such as compilers to Fault Trees [START_REF] Rauzy | Mode automata and their compilation into fault trees[END_REF][START_REF] Prosvirnova | Automated generation of minimal cut sets from altarica 3.0 models[END_REF] and Markov chains, stochastic and stepwise simulators. It is in the core of several commercially distributed integrated modeling and simulation environments and has been successfully used to perform industrial scale experiments [START_REF] Akerlund | Isaac, a framework for integrated safety analysis of functional, geometrical and human aspects[END_REF].

In this article we study the mapping between AADL EMV2 and AltaRica concepts. To illustrate our purpose, we use as a study case a simplified flight control system of a small UAV quadcopter with a particular focus on safety mitigation architectures, also called safety patterns. The transformation of AADL EMV2 models to AltaRica is interesting because it enables us to enlarge the set of safety assessment tools for AADL and to perform cross check verifications. Indeed, to implement the transformation we use the standardized version of Al-taRica Data-Flow, which is a subset of AltaRica 3.0 [START_REF] Prosvirnova | The altarica 3.0 project for model-based safety assessment[END_REF]. AltaRica 3.0 is supported by the OpenAltaRica platform, which is free of use for research and education purposes. This platform already includes a Fault Tree compiler and a stepwise simulator. A stochastic simulator and a sequence generator are currently under development.

The remainder of this article is organized as follows. Section 2 describes a simplified flight control system of a small UAV quadcopter which is used as a running example of this article. Section 3 gives an overview of AADL EMV2 concepts. Section 4 introduces the AltaRica modeling language. Section 5 explains AADL to AltaRica translation principles. Section 6 presents some analysis results of this study. Section 7 summarizes related works. Section 8 concludes this article.

Running example

The case study is inspired by the flight control system of a small UAV quadcopter. We address a simplified software architecture that encompasses different significant safety patterns. High level functional architecture The control is achieved by three main functions. The ComputeCommand function aims at computing the commands to control the helixes in automated or manual mode. The function needs two different flight parameters to be performed correctly in the automated mode whereas it needs at least one flight parameter and the pilot order in the manual mode. The quality of the inputs and outputs of this function is checked by the function ControlFlightParameters. The function ManageFlightMode adapts the piloting mode and the applicable order according to the checked status. Failure modes and failure conditions Each function may fail and its computation may be erroneous or lost. The loss of control is not catastrophic as long as the crash is controlled in an acceptable area. The erroneous control is a worse case because the UAV may fly away and lose the separation with other aircraft or it may crash over populated area. So the system safety assessment is needed to compute the causes and probability of occurrences of such failure conditions.

In this study we consider two failure conditions:

-FC01: "Quadcopter fly away, i.e. erratic behavior, potentially fly and crash in an unauthorized area, leading at worst to fatalities". -FC02: "Uncontrolled crash, i.e. loss of the quadcopter control".

Safety patterns

Several safety patterns are introduced to reduce the failure occurrences or their effects. A COMMON (COMmand and MONitoring) architecture is proposed to detect and mitigate the effect of potential erroneous computations. The pattern, shown Fig. 2, works as follows. We assume that the computation is achieved in parallel by two different channels so that a computation error of one channel may be detected by comparison with the other channel. If no alarm is raised, the order computed by the command lane is applied. Otherwise, the computed order is no more selected and the order is lost.

The mode manager is proposed to adapt as long as possible the control architecture according to the integrity and availability of its basic functions. Initially, the engaged mode is the automated mode. When the check of the flight parameters detects that the automated flight is no longer safe, the manual mode is engaged. A crash mode is engaged when the manual mode is also estimated unsafe.

Modeling and Safety assessment needs

This case study is modeled using the AADL notation and the AltaRica formal language to better understand each notation and the principles of the translation of the AADL model into an AltaRica one. Then the translation is applied to produce an AltaRica model from the AADL one. AltaRica fault tree generator [START_REF] Prosvirnova | Automated generation of minimal cut sets from altarica 3.0 models[END_REF] is used for both the hand made AltaRica model and the generated one. The comparison of the Minimal Cut Sets for FC01 and FC02 contributes to the validation of our translation process.

AADL EMV2 presentation

The SAE "Architecture Analysis and Design Language" (AADL) [START_REF]SAE: Architecture Analysis and Design Language (AADL) AS-5506[END_REF] is a language for model-based engineering of embedded real-time systems. The AADL allows for the description of both software and hardware parts of a system. It focuses on the definition of clear block interfaces, and separates the implementations from these interfaces. From the separate description of these blocks, one can build an assembly of blocks that represents the full system. To take into account the multiple ways to connect components, the AADL defines different connection patterns: subcomponent, connection, binding.

An AADL description is made of components. Each component category describes well-identified elements of the actual architecture, using the same vocabulary of system or software engineering. The AADL standard defines software components (data, thread, subprogram, process), execution platform components (memory, bus, processor, device, . . . ) and composite components (system, abstract). Besides, the language defines precise legality rules for component assemblies, and both its static and execution semantics.

The AADL defines the notion of properties. They model non-functional properties that can be attached to model elements (components, connections, features, instances, etc.). Properties are typed attributes that specify constraints or characteristics that apply to the elements of the architecture such as clock frequency of a processor, execution time of a thread, bandwidth of a bus, implementation of the functional part.

From this core of elements, AADL allows the designer to attach annex elements that further refine one dimension of the system. The Error Modeling Annex V2 (EMV2) addresses safety modeling concerns. See [START_REF] Delange | Architecture fault modeling with the aadl error-model annex[END_REF] for more details.

EMV2 supports architecture fault modeling at three levels of abstraction:

-Focus on error propagation between system components and with the environment: modeling of fault sources, along with their impact through propagation. It allows for safety analysis in the form of hazard identification, fault impact analysis, and stochastic fault analysis. -Focus on component faults, failure modes, and fault handling: fault occurrences within a component, resulting fault behavior in terms of failure modes, effects on other components, the effect of incoming propagations on the component, and the ability of the component to recover or be repaired.

It allows for modeling of system degradation and fail-stop behavior, specification of redundancy and recovery strategies providing an abstract error behavior specification of a system without requiring the presence of subsystem specifications. -Focus on compositional abstraction of system error behavior in terms of its subsystems.

In addition, EMV2 introduces the concept of error type to characterize faults, failures and propagations. It includes a set of predefined error types as starting point for systematic identification of different types of error propagations providing an error propagation ontology. Users can adapt and extend this ontology to specific domains.

As an illustration consider an AADL EMV2 model given below. It represents a function with an input and an output. It may be in three states: s Ok, s Erroneous and s Lost, representing respectively the nominal behavior, the erroneous behavior and the loss of the function. If the function is in the state s Ok and it receives an error on its input, it propagates it on the output. If the function is in the state s Erroneous it propagates an error on its output and so on. EMV2 follows regular convention for the description of state transition and error propagation, initial state -[trigger] → error event , that reads as follows: when the system is in state initial state , it propagates the corresponding error event . This propagation may be controlled by the trigger.

abstract

AltaRica presentation

AltaRica is a high level formal modeling language dedicated to Safety Analyses [START_REF] Arnold | The altarica language and its semantics[END_REF]. Its Data-Flow version has been created to handle industrial scale models [START_REF] Boiteau | The altarica data-flow language in use: Assessment of production availability of a multistates system[END_REF]. A number of assessment tools have been developed ( [START_REF] Rauzy | Mode automata and their compilation into fault trees[END_REF], [START_REF] Teichteil-Knigsbuch | Lazy forward-chaining methods for probabilistic model-checking[END_REF]). AltaRica Data-Flow is at the core of several Modeling and Simulation tools and has been successfully used for industrial applications [START_REF] Akerlund | Isaac, a framework for integrated safety analysis of functional, geometrical and human aspects[END_REF]. In 2011, an initiative was launched to standardize the syntax of AltaRica Data-Flow.

To implement the transformation of AADL to AltaRica we use the standardized version of AltaRica Data-Flow -a subset of AltaRica 3.0 [START_REF] Prosvirnova | The altarica 3.0 project for model-based safety assessment[END_REF], supported by the OpenAltaRica platform. This platform is developed by IRT SystemX and is free of use for research and education purposes. It already includes a Fault Tree compiler [START_REF] Prosvirnova | Automated generation of minimal cut sets from altarica 3.0 models[END_REF] and a stepwise simulator. Other tools are under development.

In this article we only focus on concepts of AltaRica Data-Flow illustrated using the running example. The interested reader can refer to [START_REF] Batteux | AltaRica 3.0 specification. Tech. rep[END_REF] to know more about AltaRica 3.0.

Basic blocks

The following AltaRica code represents the behavior of a basic function.

domain BasicFunctionStatus {OK, LOST, ERRONEOUS} class BasicFunction BasicFunctionStatus status (init = OK); event fail_loss (delay = exponential(0.001)); event fail_error (delay = exponential(0.0005)); transition fail_loss: status == OK -> status := LOST; fail_error: status == OK -> status := ERRONEOUS; end States: The internal state of the function is represented by means of the state variable status, which takes its value in the domain BasicFunctionStatus. So, the function can be in three states: OK representing the nominal behavior, LOST (loss of the function), and ERRONEOUS (erroneous behavior). The initial value of the state variable is specified by the attribute init. Events: The state of the function changes under the occurrence of an event, introduced with the keyword event. In our example, the function has two failure events: fail error (erroneous behavior), and fail loss (loss of the function). A delay is associated with each event by means of the attribute delay. Delays of the events fail loss and fail error are random exponentially distributed variables. Transitions: A transition is a triple e : G → P , where e is an event, G is a Boolean expression, the so-called guard of the transition, P is an instruction, the so-called action of the transition. In the example above if the state of the function is OK, then two transitions are fireable: the transition labeled with the event fail loss and the transition labeled with the event fail error. If the delay drawn for the transition fail loss is the shortest, then this transition is fired and the variable status is switched to LOST.

then input else status; end

There are two flow variables: input and output, taking their values in the domain BasicFunctionStatus. They represent respectively the quality of the data received and sent by the function. The assertion states that if the state of the function is OK, then its output is equal to its input, if its state is LOST then its output is also LOST, otherwise it is ERRONEOUS.

Hierarchical models

In AltaRica Data-Flow components are represented by classes. Classes can be instantiated in other classes in order to create hierarchical models. Their inputs and outputs can be connected via assertions.

The Figure 2

. end

In this model, we define an immediate event GoToManualMode and the associated transition, which represents the reconfiguration: while in the automated mode, if the alarm is received (inputAlarm is true), then the transition is fired immediately and the mode is switched to manual.

AADL EMV2 to AltaRica translation

Models in AADL with EMV2 and AltaRica are structurally similar so that it is possible to translate one notation into the other. On the structural side AADL components correspond to AltaRica classes. AADL subcomponents correspond to AltaRica variables whose type is an AltaRica class. Thus, given an AADL model it is possible to create an AltaRica model that exhibits the same hierarchical containment structure.

For safety analysis we are interested in the occurrence and propagation of faults throughout the analyzed system. The interface of an AADL component with regard to fault propagation is given by its error propagations as defined in the component's EMV2 annex subclause. Each of these EMV2 error propagations corresponds to an AltaRica flow variable whose domain can be derived from the set of faults that are associated with the error propagation.

Several AADL model constructs are involved in the propagation of faults from one AADL component to another. (a) Connections modeling data and control flow can propagate faults related to, for example, data validity and timing. (b) Faults that occur in an execution platform component, e.g. power loss, propagate to software components that are bound to, i.e. propagation via binding properties. (c) Additional propagations without an explicit path in the architectural model are defined in an EMV2 annex subclause, e.g. fault propagation due to heat transfer between hardware components located in close proximity. All these constructs are translated into "external" AltaRica assertions, i.e. assertions connecting flow variables from different classes.

The internal fault behavior of an AADL component is given using the AADL EMV2 error behavior and component error behavior constructs. The translation of these constructs to AltaRica proceeds as follows:

-An EMV2 error event describes the occurrence of an internal fault that happens in the component. Each such event is translated into an AltaRica event. The occurrence probability of an error event is given by the value of property EMV2::OccurrenceDistribution which is translated into a delay attribute for the AltaRica event. -EMV2 error states are defined as an identifier, i.e. the state's name. All error states of an AADL component are translated into a single state variable in AltaRica. The domain of this variable is the set of symbols created from the EMV2 error state names. -EMV2 state transitions are translated to AltaRica transitions. If the the transition is caused by an internal fault, i.e., an error event, it is translated into an AltaRica transition that is enabled by the corresponding AltaRica event and uses the source state as the guard. The action assigns the target state to the state variable. If, on the other hand, the transition is caused by incoming fault propagations, we create an AltaRica transition that is always enabled and translate the EMV2 error condition into the guard condition.

To enable this kind of transition we add an event with the attribute delay equal to Dirac(0.0) (immediate event) to each AltaRica class. -EMV2 out propagation conditions determine the error type produced at an out propagation based on the state and an error condition involved in propagations. Each out propagation condition is translated to an internal AltaRica assertion that sets the value of an outgoing flow variable to the symbol representing the propagated error type.

The following two listings show a simple example in AADL, the comparator component used in the COMMON pattern, and the AltaRica code generated by the translator. It is easy to see the correspondence between the two models. event error_propagation (delay = Dirac (0.0)); // not used domain_4 out0 (reset = noerror); // initially no error domain_5 in2_3; domain_5 in1_4; assertion out0 := switch { // all out propagation conditions // aggregated into one switch statement case in1_4 == err and in2_3 == lost: lost case in1_4 == lost and in2_3 == lost: lost case in1_4 == lost and in2_3 == err: lost case in1_4 == noerror and in2_3 != noerror: lost case in1_4 != noerror and in2_3 == noerror: lost case in1_4 == err and in2_3 == err: err default: noerror // noerror is the default }; end Modeling the comparison of error types as several out propagation conditions results in a somewhat unwieldy AltaRica assertion since the translator does not perform any simplification of the generated code. The assertion is equivalent to out0 := if (in1_4 == in2_3) then in1_4 else lost;

Resolution of mismatches between AADL EMV2 and AltaRica

There are a couple of differences between AADL EMV2 and AltaRica which must be taken into account when translating between the two formalisms.

One difference concerns the way error types are defined in EMV2. In general error types can be thought of as typed tokens that propagate through an architecture. However, EMV2 error types are organized in a generalization hierarchy, and all EMV2 allows use of generalized error types wherever an error type can occur. For example, AboveRange is a subtype of OutOfRange, which is a subtype of DetectableValueError. The most general error in this generalization chain is the ItemValueError. AltaRica does not support a notion of generalization in the definition of constants used to create domains. We solve this mismatch by always replacing generalized error types with the most specific error types. For example, a DetectableValueError in an error propagation or condition is replaced with the most specific error types of which DetectableValueError is a generalization, namely OutOfBounds, BelowRange, and AboveRange.

Another difference is that EMV2 uses sets of error types to define error propagations and conditions, whereas AltaRica does not support a built-in notion of sets. This leads to difficulties when assigning a domain to an incoming flow variable if the flow is generated based on an AADL feature that has multiple incoming connections (fan-in). Such an AADL model is valid if the set of error types F i at each of the connected out propagations o i are contained in the type set at the in propagation: ∀i ∈ {1, ..., n}F i ⊆ E. In general all type sets F i may be different. To resolve this mismatch we generate multiple flow variables from an EMV2 in propagation, one per incoming connection. As the domain we use the domain generated from the error of the corresponding out propagation at the other end of the connection. This way both flow variables have the same domain as their type and can be connected using an AltaRica assertion.

Unfortunately, the 1 to n translation of in propagations to flow variables complicates the translation of error conditions on transitions and out propagation conditions. The error conditions contain atomic terms of the form ip(C), with C = {t 1 , t 2 , ..., t k } a set of error types. Such a term is true if and only if one of the error types is propagated in via the in propagation ip. In the simplest case this term is translated into the following Boolean expression in AltaRica:

ip = t 1 ∨ ip = t 2 ∨ ... ∨ ip = t k .
When the in propagation is split into several flow variables, the resulting AltaRica expression is the disjunction of the expressions generated for each of the new flow variables, which can become long and difficult to read.

Another consequence of the 1 to n relationship between in propagations and flow variables is that the number of generated flow variables is not the same for all instances of an AADL classifier. It depends on the context in which the component is used, in particular, how it is connected to other components in the model. Therefore, it is not possible to generate the AltaRica classes based on the declarative AADL model alone. Instead, we translate each AADL component instance into an AltaRica class, potentially resulting in an AltaRica model with several classes that are identical except for identifiers. This increases the size of the AltaRica model but has not much influence on analyses of the models. As an AltaRica model is "flattened" to transform it into a guarded transition system, it is necessary to insert a full copy of a class for each class instance. This is similar to what happens during the generation of an AADL instance model, thus resulting in essentially the same guarded transition system.

AADL EMV2 also has the concept of a typed error state. This is an error state that has an error type associated with it. The error types allowed for a state are enumerated in a type set given with the error state declaration. Even though we have not included typed error states in the translation, it could be extended to generate multiple AltaRica domain constants for each typed EMV2 error state, one per declared type token.

Implementation of the translator

We have implemented the AADL EMV2 to AltaRica translator as a plugin to OS-ATE, the open source AADL tool environment. The translation is implemented using the Atlas Transformation Language. The source code for our translator is available on github6 , and the plugin can be installed into OSATE from our p2 repository. 7The current version supports a subset of AADL EMV2 sufficient for the analysis of the UAV system used in our case study. We are planning to extend that translation to include more EMV2 concepts as needed, e.g. error detections and stochastic state transitions.

Experiments

The automatically generated AltaRica model can be analyzed with the various tools which are provided within the OpenAltaRica platform. In particular, we can simulate the model in order to observe the effect of different failure scenarios on the overall system. We can also analyze the potential causes of the failure conditions by generating Fault Trees and Minimal Cut Sets (MCS).

Let us consider the failure condition FC02 (uncontrolled crash). The analysis returns 12 MCS presented in Fig. 6.

This result fully complies with the result that we obtained with the handwritten AltaRica model of the quadcopter. As a first remark, all the MCS are of order 2. So, we can deduce that no single failure can lead to FC02. As we can see, the pilot order is of prime importance. This is not surprising since one of the main ideas behind this system is to switch from automatic mode to manual mode in case of a problem (loss of a source for the estimation of flight parameters, loss of the automatic command, . . . ). Once the system switches to manual mode, if the function that acquires the pilot order also fails, then an uncontrolled crash occurs. A similar argument holds for the alarm of the function ControlFlightParameters, which is responsible for detecting a problem and switching to manual mode.

The same analysis has been performed for the failure condition FC01. Again, the result is similar to the result obtained from the hand-written model. For this failure condition, there are 2 MCS of order 1, which correspond to the erroneous Fig. 3. Minimal cut sets for FC2: uncontrolled crash behavior of each of the sources for the estimation of the flight parameters. This was expected since the system does not implement any mitigation for this specific failure (contrary to the loss of one source, or the erroneous behavior of the command computation).

We have obtained the same results for the handwritten and the automatically generated AltaRica models. It is a first validation step for our translator.

Related work

Another translation of AADL to AltaRica, based on the definition of ontologies, has been proposed in [START_REF] Mokos | Towards compositional safety analysis via semantic representation of component failure behaviour[END_REF].

Translations of AADL models into classical safety models exist (see e.g. [START_REF] Delange | Architecture fault modeling with the aadl error-model annex[END_REF]). It is worth noting that each type of translation takes advantage of different features of the AADL model according to the targeted safety models. For instance, the translation into a Fault Tree requires the features composite error behavior and error events. Such features describe static dependencies between failure modes of components and are well suited for the compilation into Boolean formulas. Component error behavior gives a more dynamic view of the failure propagation and is useful for the computation of probability of occurrences of states of Markov chains or stochastic processes. AltaRica has been designed to integrate in one formalism both static and dynamic features of the failure propagation; it supports not only the computation of the probability of occurrences of states but also the computation of sets/sequences of events leading to undesired states. So, the proposed translation should ease the analysis of reconfigurable systems with multiple modes as our flight control system.

A translation of an AADL subset into NuSMV formal language is a pioneer work to achieve more powerful safety analysis in the framework of the COMPASS project [START_REF] Bozzano | Safety, dependability and performance analysis of extended aadl models[END_REF]. NuSMV language is close to AltaRica Data-Flow: they have similar expressive power. Several safety tools have been developed for NuSMV. However, AltaRica was designed to support safety analysis whereas NuSMV is a more general language. Specific annotations are introduced in a native NuSMV model to distinguish error concepts that should appear in a safety report. The available tools are effective [START_REF] Akerlund | Isaac, a framework for integrated safety analysis of functional, geometrical and human aspects[END_REF], but the annotation process raises practical questions when mixed with AADL view. First versions of COMPASS built the annotations by asking end users to extend AADL nominal models in the COMPASS framework. It could be clearer for end users to perform the model extension in an AADL framework compliant both with the principle of the AADL Error Annex and the level of rigor requested by NuSMV as experimented in our work.

Conclusions

The AADL is a language that is used to describe the software and hardware architectures of system. The AADL EMV2 enables the extension of the architecture models with artefacts that are relevant to generate safety related models. For instance, Fault Trees can be generated from models extended according the principle of the EMV2. AltaRica is a formal language specifically defined to support advanced safety analysis of complex systems that integrate modes and dynamic reconfiguration on specific conditions. We have presented an example of a simplified version of such a kind of dynamic systems and we have modeled it both in AltaRica and in the AADL as starting point. This first activity has enabled a better understanding of the philosophy of both approaches. Both languages have to deal with failure propagations inside a system architecture. In the AADL, the failure extension of a component is flexible and it enables the connexion of components that are subject to heterogeneous failure modes. As a counterpart, end users need to be quite explicit on the various cases of composition of the failure modes. In AltaRica, the connexion of component shall be compliant with stronger typing rules. So end users can less easily propagate heterogeneous failure modes. As a counterpart, the connexion is straight forward and the propagation between connected components results directly from the language semantics. Then we have developed a translator from the AADL EMV2 into AltaRica and we have generated an AltaRica model for our case study. The manual Al-taRica model and the generated one have been simulated and analyzed. They have produced the same results. This is a first encouraging validation step for our translator. Moreover, the analysis exhibits relevant sequences of causes of the undesired events for a multi-mode system. Thus the AADL model can be analyzed with latest generation of safety assessment tools. Further works will aim at short term to more widely validate the translator and the benefit of the coupled approach. At longer term, the AADL EMV2 could be updated to ease modeling in the spirit of some interesting AltaRica features. Conversely, the flexibility of the failure type definition in AltaRica could be extended to better account for extension approach like in the AADL EMV2.

Fig. 1 .

 1 Fig. 1. Overview of the flight control software architecture

Fig. 2 .

 2 Fig. 2. Command and Monitoring pattern

  Flow propagations: In AltaRica the propagation of errors/failures and nominal values is done in the same way: via flow variables and assertions. The value of flow variables are recalculated after each transition firing by means of assertions. Assertions are instructions as are actions of transitions. The difference is that actions of transitions assign only state variables, while assertions assign only flow variables. Consider, for example, the following AltaRica code: class BasicInOutFunction extends BasicFunction; BasicFunctionStatus input, output (reset = LOST); assertion output := if (status==OK)

  is a graphical view of the class COMMONPattern. This class contains two instances of the class BasicInOutFunction, named COM and MON, one instance of the class Selector named selector and one instance of the class Comparator named comparator. The Comparator (cf code below) and Selector are considered as component connectors that are free of failure modes in our case study. Other direct connections between components are represented by plane lines in the figure and by equality assertions linking input and output of connected components in the class COMMONPattern. class Comparator BasicFunctionStatus input1 (reset=LOST); BasicFunctionStatus input2 (reset=LOST); Boolean alarm (reset=false); assertion alarm:= if (input1==input2) then false else true; end Flight modes and reconfigurations Flight modes can be represented by a state variable mode which takes its value in the domain FlightModeDomain. The reconfigurations are represented by immediate transitions, introduced by the attribute delay equal to Dirac(0). They are fireable as soon as their guards are true. As an illustration consider a part of the AltaRica model representing the function ManageFlightMode of the running example. domain FlightModeDomain{AUTO, MANUAL, CRASH, DANGER} class ManageFlightModeFunction Boolean inputAlarm, inputCrashAlarm (reset = false); FlightModeDomain mode(init = AUTO); event GoToManualMode (delay = Dirac(0)); ... transition GoToManualMode:(mode == AUTO) and inputAlarm and not inputCrashAlarm -> mode := MANUAL; ..

  Note that in EMV2 it is not possible to compare error types with each other. For example, the Comparator class as shown in the code in Section 4 uses the expression input1 == input2. Such comparisons must be modeled in EMV2 by using several out propagations that enumerate all possible combinations of error types occurring at the two propagation points.

		--propagated error if in1 != in2
	all -[in1{noerror} and in2]-> out0{lost};
	all -[in1 and in2{noerror}]-> out0{lost};
	all -[in1{lost} and in2{err}]-> out0{lost};
	all -[in1{err} and in2{lost}]-> out0{lost};
		--default: no error is propagated
	end component;
	**};
	end Comparator;
	domain domain_4 {	// domain names are generated
	noerror, lost, err	// constants generated from error types
	}
	domain domain_5 {
	noerror, lost, err
	}
	class Comparator_6
	abstract Comparator
	features
	in1: in data port;	--ports defined in AADL core
	in2: in data port;	--language
	out0: out data port;
	annex EMV2 {**
	use types FailMode;
	error propagations	--propagated errors defined in
		--EMV2 annex
	in1: in propagation {lost, err};
	in2: in propagation {lost, err};
	out0: out propagation {lost, err};
	end propagations;
	component error behavior --no error states needed
	propagations
		--propagated error if in1 = in2
	all -[in1{lost} and in2{lost}]-> out0{lost};
	all -[in1{err} and in2{err}]-> out0{err};

at http://github.com/osate/aadl2altarica

at http://aadl.info/aadl/osate/experimental

Acknowledgements

This material is based upon work funded and supported by the Department of Defense under Contract No. FA8721-05-C-0003 with Carnegie Mellon University for the operation of the Software Engineering Institute, a federally funded research and development center. [Distribution Statement A] This material has been approved for public release and unlimited distribution. Please see Copyright notice for non-US Government use and distribution. DM-0004294


Auguste Gires

Ioulia Tchiguirinskaia

D Schertzer

Multifractal comparison of the outputs of two optical disdrometers

Keywords: radar relations multifractals, disdrometer, radar 1 1

come

Multifractal comparison of the outputs of two optical disdrometers

1) Introduction

Rainfall measurement with the help of disdrometers is rapidly developing for point measurements. The first ones were impact disdrometers [START_REF] Joss | Ein spektrograph fur nieder chlagstropfen mit automatischer auswertung (A spectrograph for raindrops with automatic interpretation)[END_REF], and now optical ones [START_REF] Loffler-Mang | An Optical Disdrometer for Measuring Size and Velocity of Hydrometeors[END_REF][START_REF] Ellis | New laser technology to determine present weather parameters[END_REF][START_REF] Battaglia | PARSIVEL Snow Observations: A Critical Assessment[END_REF][START_REF] Frasson | Assessment of the Thies optical disdrometer performance[END_REF]) are more commonly used for operational and research purposes. The great advantage of these devices with regards to more conventional rain gauges is that they do not measure only a rain rate (or rainfall depth) but information about size and fall velocity for all the hydrometeors passing through the sampling area whose size is few tens of cm 2 . The Drop Size Distribution (DSD) can then be computed from this raw data. From the DSD it is possible to estimate numerous rain related fields such as the rain rate (R) or the liquid water content () or even quantities measured directly by polarimetric weather radars such as the horizontal reflectivity (Zh) and the specific differential phase (Kdp) (Jaffrain and Berne 2012[START_REF] Leinonen | A Climatology of Disdrometer Measurements of Rainfall in Finland over Five Years with Implications for Global Radar Observations[END_REF], Verrier et al. 2013). This widens a lot the range of applications for disdrometer data.

Numerous experiments have been carried out to compare the rain rate output of various types of collocated disdrometers along with rain gauges [START_REF] Miriovsky | An Experimental Study of Small-Scale Variability of Radar Reflectivity Using Disdrometer Observations[END_REF][START_REF] Krajewski | DEVEX-disdrometer evaluation experiment: Basic results and implications for hydrologic studies[END_REF][START_REF] Frasson | Assessment of the Thies optical disdrometer performance[END_REF][START_REF] Thurai | Drop size distribution comparisons between Parsivel and 2-D video disdrometers[END_REF]. The temporal evolution [START_REF] Thurai | Drop size distribution comparisons between Parsivel and 2-D video disdrometers[END_REF] or spatial distribution (Jaffrain and Berne 2011) of total drop concentration (Nt) and mass-weighted diameter (Dm), which are commonly used to fully characterize drop size distribution, have also been analysed. However these comparisons are usually done only at a single resolution, most commonly the maximum one available (i.e. using the time series at the recording time step). In this paper we suggest to carry out the comparison not only at a single scale but across scales. Considering various scales at once enables to obtained results robust over a wide range of scales. To achieve this, a theoretical framework relying on scale invariant properties is implemented. Available data enables to carry out analysis on scales ranging from 30 s to approximately 2 hours for a punctual measure. Moreover not only R, , or Nt and Dm, but also the radar quantities Zh and Kdp will be investigated. The output data of two disdrometers, which have been deployed for few months on the roof of the Ecole des Ponts ParisTech building, will be used for this study. The two disdrometers are a Campbell Scientific PWS100 [START_REF] Ellis | New laser technology to determine present weather parameters[END_REF], Campbell Scientific Ltd 2012), which has seldom been used in such study, and an OTT Parsivel 2 [START_REF] Battaglia | PARSIVEL Snow Observations: A Critical Assessment[END_REF][START_REF] Ott | Operating instructions[END_REF].

The variability across scales will be quantified with the help of stochastic universal multifractals which have been extensively used to analyze, model, and simulate geophysical fields extremely variable over wide range of scales such as rainfall [START_REF] Schertzer | Physical modelling and analysis of rain and clouds by anisotropic scaling and multiplicative processes[END_REF], 1997[START_REF] Marsan | Causal space-time multifractal processes: Predictability and forecasting of rain fields[END_REF][START_REF] Olsson | Multifractal analysis of daily spatial rainfall distributions[END_REF][START_REF] Gires | Improvement of measurement with a refraction disdrometer by better taking into account the drops oblateness[END_REF][START_REF] De Lima | Multifractal analysis of 15-min and daily rainfall from a semi-arid region in Portugal[END_REF][START_REF] Lovejoy | Scale, scaling and multifractals in geophysics: Twenty years on[END_REF][START_REF] Nykanen | Linkages between Orographic Forcing and the Scaling Properties of Convective Rainfall in Mountainous Regions[END_REF], Royer et al. 2008[START_REF] De Lima | Investigating the multifractality of point precipitation in the Madeira archipelago[END_REF][START_REF] Mandapaka | Multiscaling analysis of high resolution space-time lidar-rainfall[END_REF][START_REF] De Montera | The Effect of Rain-No Rain Intermittency on the Estimation of the Universal Multifractals Model Parameters[END_REF], Verrier et al. 2010). This framework has seldom been applied to Zh (see [START_REF] Tessier | Universal Multifractals: theory and observations for rain and clouds[END_REF] for an example in space and not in time as here), and not to Kdp to the knowledge of the authors.

Data and the retrieval of the various studied fields are presented in section 2 along with a brief reminder of the Universal Multifractals (UM) framework and a presentation of the methodology implemented. Results are in section 4 where classical

R Z h  of dp K R 
relations are investigated with this data set, and scaling features and UM parameters estimates of the various fields are discussed.

2) Data and methods

2.1) Description of the disdrometer data

The data used in this paper was collected between 27 September 2013 and 18 January 2014 by two disdrometers installed on the roof of the Ecole des Ponts ParisTech building.

They are both made of a transmitter that generates one or several laser sheet(s) and receiver(s). The Campbell Scientific PWS100 computes size and fall velocity from light refracted by the hydrometeors [START_REF] Ellis | New laser technology to determine present weather parameters[END_REF], Campbell Scientific Ltd 2012) whereas the OTT Parsivel 2 relies on occluded light (Loffler-Mang and Joss 2000, for an initial version; [START_REF] Battaglia | PARSIVEL Snow Observations: A Critical Assessment[END_REF][START_REF] Ott | Operating instructions[END_REF]. The 50 recorded events of this period are used in this study.

The criteria defining an event is a rainy period during which more than 1 mm is collected and that is separated by more than 15 min of dry conditions before and after. The main output of the disdrometers is a matrix with the number nij of drop recorded according to classes of equivolumic diameter (index i, and defined by a centre Di and a width i D  expressed in mm) and terminal fall velocity (index j, and defined by a centre vj and a width j v  expressed in m.s -1 ). This matrix is recorded for each 30s time step ( t  ). [START_REF] Gires | Improvement of measurement with a refraction disdrometer by better taking into account the drops oblateness[END_REF] which used the same data set noticed that the oblateness of drop was not properly taken into account in the PWS100 rationale and suggested a correction which is used here. Furthermore, as suggested by various authors (Jaffrain and Berne 2011[START_REF] Kruger | Two-Dimensional Video Disdrometer: A Description[END_REF][START_REF] Thurai | Drop Axis Ratios from a 2D Video Disdrometer[END_REF] all the drops whose velocity was more than 60% different from what was expected by Beard's model [START_REF] Beard | Terminal velocity adjustment for cloud and precipitation aloft[END_REF] according to its size were removed because considered as non meteorological measurements.

The rain rate for each time step is then computed as:

   j i i eff i Pars PWS j i Pars PWS D S D n t R , 3 , / , / ) ( 6  (1)
where

) ( i eff D S
is the sampling area of the device which is slightly modified according to the drop size to take into account potential edge effects for large drops. For the Parsivel 2 we have

        2 ) ( i i eff D W L D S
, where L = 180 mm and W = 30 mm are respectively the length and width of the sensing area (LW = 54 cm 2 ) (OTT 2014). The PWS100 is not subject to this issue and eff S is taken as constant equal to 40 cm 2 (Campbell Scientific Ltd 2012). For the selected events the total rainfall depth collected is equal to 172 mm for the PWS100 and 154 mm for the Parsivel 2 . A tipping bucket rain gauge is also located on the same roof and collected roughly 170 mm which is more in agreement with the PWS100, but it should not be "over-interpreted" given that its accuracy is not very high due to high rain rates and low level of maintenance at the beginning of the recording period. The normalized bias between the two disdrometers (computed for the time steps where R > 1mm/h) is equal to 0.11; the correlation equals 0.96; the Nash-Sutcliffe efficiency coefficient equals 0.89 and the RMSE is 1.83. Although not negligible the differences between the two disdrometers are rather low compared to what is commonly observed [START_REF] Miriovsky | An Experimental Study of Small-Scale Variability of Radar Reflectivity Using Disdrometer Observations[END_REF][START_REF] Krajewski | DEVEX-disdrometer evaluation experiment: Basic results and implications for hydrologic studies[END_REF][START_REF] Frasson | Assessment of the Thies optical disdrometer performance[END_REF][START_REF] Thurai | Drop size distribution comparisons between Parsivel and 2-D video disdrometers[END_REF].

A discrete drop size distribution (DSD) is computed from the available data as:

    j j j i i i eff i v n t D D S D N , ) ( 1 ) ( (2)
The number of drops with a diameter in the class i per unit volume (in m -3 ) is given by N(Di) Di. As it is commonly done, we use the total drop concentration Nt (m -3 ) and massweighted diameter Dm (mm) to characterize the DSD. In this framework the DSD is written as

) ( ) ( m t D f N D N 
; where it appears that Dm characterizes the shape of the DSD and Nt its total intensity. These two parameters are defined as (Jaffrain and Berne 2012a, Leinonen et al. For the practical computation of these quantities (and the other ones), the integral is replaced by a summation over all the classes of diameter, i.e.:

    i i i i D D D D N D f dD D N D f ) ( ) ( ) ( ) ( max min (5)
Lastly the liquid water content ( in g.m -3 ) the horizontal reflectivity (Zh in mm 6 .m -3 ) and the specific differential phase (Kdp in °.km -1 ) are estimated with the help of the DSD from which they are basically various moments. These physical parameters are given by Jaffrain and Berne (2012a) and [START_REF] Leinonen | A Climatology of Disdrometer Measurements of Rainfall in Finland over Five Years with Implications for Global Radar Observations[END_REF]:

        max min 3 3 2 ) ( 3 4 10 1 D D w l dD D D N    (6)
With w,the liquid water density in g.m -3 ,

    max min ) ( ) ( 1 1 ; 2 2 5 4 , D D h B v h dD D D N m m Z    (7)       max min ) ( ) ( ) ( Re 1800 10 3 D D vv hh dp dD D N D S D S K   (8)
Where

h B;
 (in mm 2 ) is the backscattering cross session at horizontal polarization, and Re(Shh/vv) (in mm) is the real part of the forward scattering amplitude at horizontal/vertical polarization,  is the radar wavelength (in mm) and m the complex refractive index of water.

The scattering coefficients were computed with the help of the Python PyTMatrix library [START_REF] Leinonen | High-level interface to T-matrix scattering calculations: architecture, capabilities and limitations[END_REF]) which relies on the T-Matrix code by [START_REF] Mishchenko | T-matrix computations of light scattering by nonspherical particles: A review[END_REF].

Computations were carried out for a radar wave length equal to 53.5 mm corresponding to Cband radars (this wave length was chosen to facilitate comparison with other studies, because it is the most widely used and corresponds to the radar currently mostly used by Western European meteorological services), a temperature of 20°C (m=8.633 +1.289i) and an oblate spheroids model for drop shape with an axis ratio -equivolumic diameter relation corresponding to the one implemented in the Parsivel 2 rationale [START_REF] Battaglia | PARSIVEL Snow Observations: A Critical Assessment[END_REF]. The set up considered for drop orientation was the same as in [START_REF] Leinonen | A Climatology of Disdrometer Measurements of Rainfall in Finland over Five Years with Implications for Global Radar Observations[END_REF]; i.e. drops are partially aligned and a normal distribution (mean and standard deviation respectively equal to 0° and 7°, in agreement with the findings of [START_REF] Bringi | Measurements and inferences of raindrop canting angles[END_REF]) characterizes the angle of the symmetry axis.

2.2) Methodology: UM framework

The key elements of the theoretical framework of UM are presented here and the reader is invited to refer to Schetzer and [START_REF] Schertzer | Multifractals, generalized scale invariance and complexity in geophysics[END_REF] for a recent review. In the following   denotes a conservative field at resolution , defined as the ratio between the observation scale l and the outer scale L (=L/l). Practically   is obtained simply by up-scaling (averaging consecutive time steps) the measured field at the maximum resolution. If   is a multifractal field, then its statistical moment orders scale with resolution as:

) (q K q     (9)
Where K(q) is the moment scaling function that fully characterizes the variability across scales of the field   . The quality of the scaling is investigated with the help of the Trace Moment (TM) analysis which simply consists in plotting equation 9 in log-log, the slope of the obtained straight line being K(q).

Most multiplicative processes converge toward Universal Multifractal which are fully characterized with the help of only two scale invariant parameters C1 and  (this a broad generalization of the central limit theorem, Schertzer andLovejoy 1987, 1997).

C1 is the mean intermittency co-dimension and measures the clustering of the (average) intensity at smaller and smaller scales (C1=0 for a homogeneous field).  is the multifractality index 2 0     and measures the clustering variability with regards to intensity level. In this specific framework which is implemented here K(q) is given by:

  q q C q K      1 ) ( 1 (10)
The UM parameter are estimated in this paper with the help of the Double Trace Moment (DTM) technique [START_REF] Lavallée | Nonlinear variability and landscape topography: analysis and simulation[END_REF].

A common framework to deal with a non-conservative field   (i.e. we have

1   
) is to assume it can written as

H        (11)
where H is the non-conservation parameter (H=0 for conservative fields), and   a conservative field characterized by a moment scaling function Kc(q) depending only on UM parameters C1 and . More physically, this non-conservativeness means that the studied field   (the observations here) exhibits stronger correlations than the ones obtained with the help of a simple multifractal cascade and an additional fractional integration is needed to represent it. The moment scaling function K(q) of   is given by:

Hq q K q K c   ) ( ) ( (12)
H can be estimated with the help of [START_REF] Tessier | Universal Multifractals: theory and observations for rain and clouds[END_REF]

): ) 2 ( 2 1 c K H     (13)
where  is the spectral slope. It is the exponent of the power law that characterizes over large range of wave numbers the power spectrum of a scaling field:

   k k E ) ( (14)
Before going on let us clarify the relations between the notions of stationarity and conservation. For stochastic processes, stationarity refers to the fact that given statistics are invariant with respect to time translations. When no given statistics are mentioned, it usually corresponds to the strongest case of stationarity, i.e. the probability itself is time translation invariant and therefore all the statistical moments are stationnary The classical case of "second order stationarity" is presumably the most cited case and is particularly important for additive processes. It corresponds to time translation invariance of the second order moments of the increments. For multiplicative processes, conservation means that a given statistic is strictly independent from scale, in general this statistic is the mean of the field (H=0 in Eq.

11). The latter implies a given form of stationarity because at all scales the fluctuations remain around this mean. However, these fluctuations are not only easily wilder than for an additive process, but their amplitude generally increases with smaller and smaller scales. Conversely, non-conservation implies a strong non stationarity, starting with the fact that the mean has a scaling behaviour (Eq. 11). For instance, a spectral analysis will detect the non conservation (

0  H
), due to the departure of the spectral slope, which is a second order statistic, from that of "a pink noise" (k -1 ), where the exponent one corresponds here to the dimension of the embedding space (time series are studied here). With multifractal fields there is a further correction with the Kc(2) (Eq. 13) which is not negligible (in the range 0.1-0.4 for the retrieved UM parameters).

The TM and DTM techniques are designed for analysing conservative fields (H=0)

and remain reliable as long as H<0.5. In case of greater H, they should be implemented not on   , but on the underlying conservative field   . A fractional integration of order H (equivalent to a multiplication by k H in the Fourier space) is theoretically required to estimate   from   .

However a common approximation [START_REF] Lavallée | Nonlinear variability and landscape topography: analysis and simulation[END_REF]) which provides reliable results consists in taking   at the maximum resolution simply equal to the renormalized absolute value of the fluctuations of the field i.e.:

) ( ) 1 ( ) ( ) 1 ( i i i i                (15)
and then upscaling this field at other resolutions .

Spectral and multifractal analyses are implemented on ensemble average over various samples, i.e. each sample is considered as a realization of the process. For example it means that each sample is up-scaled and taken to the power q independently before taking the average in equation 9. The studied samples are extracted in the following way from the selected events: for each event (i) a sample size is chosen (necessarily a power of two); (ii) the maximum number of samples for this event is computed; (iii) the portion of the event of length equal to the sample size multiplied by the number of samples found in (ii) with the greatest cumulative depth is extracted; (iv) the extracted series is cut into various samples.

Given that the sample size is a power of two there is obviously some data lost during the selection process. Here the percentage of available data actually used is equal to 36, 65, 82 and 91% for samples of size 512, 256, 128 and 64 times steps of 30s respectively. The chosen sample size should result from a trade off between the width of the available range of scales which should be as great as possible, and the amount of wasted data which reduces the reliability of the estimates. Here we chose to analyse the scaling properties with samples of size 256 (approximately 2h, 46 samples available), and given that no scaling break was identified on the main studied fields, UM parameters were estimated with samples of size of 64 (approximately 30 min, 259 samples available) to benefit from the use of more data.

Finally it should be mentioned that since Dm is not defined when there are no drops, only the samples containing drops at all time steps are used for this field.

3) Results and discussion 3.1) Standard radar relationships 10 Given that it was possible to retrieve both rain rates and radar parameters, it is possible to analyse the validity, for this data set, of the standard power law relations that are commonly assumed between these quantities: 1. It appears that the quality of the fitting is comparable for the two relationships and slightly better for the PWS100 than for the Parsivel 2 (r 2 ~ 0.8 vs. 0.75). The estimates are similar for both devices. 2013 for an analysis of this issue for R Z h  relation). Parameters are computed here taking into account all the events at once, meaning that "climatic" (keeping in mind only 4 months of data are used) estimates are studied. An event based analysis will be discussed in future works. It is indeed not needed for the purpose of this paper, which is to quantify the scaling variability observed by the two disdrometers on various fields.

More precisely for the

Finally the influence of using the retrieved "climatic" parameters for computing rain rates from radar parameters is assessed. With this purpose, the rain rates computed either directly from raw data (R) or through the radar relations 16 and 17 (respectively

R Z R  and dp K R R 
) once Zh and Kdp have been estimated with raw data are compared.

3.2) Scaling behaviour

The scaling features of the various fields are studied in this section. Figure 2.a displays the spectral analysis (i.e. equation 14 in log-log plot) of the rain rate for the PWS100. A very good scaling (i.e. a straight line) is observed on the whole range of available scales (30s -2h).

The spectral slope is greater than the embedding dimension of the field (1 for time series) suggesting that R is non-conservative (this will be confirmed in the next session with the estimation of H). A practical consequence is that the TM analysis which assumes a conservative field should not be implemented on the field itself but only on its conservative part which can the approximated by the absolute value of its fluctuations (equation 15). The results are displayed on Fig. 2.b, where the good scaling with a unique regime is confirmed.

The r 2 for q = 1.5 which is taken as an indication of the quality of the scaling is greater than 0.99. The same analysis carried out directly on the field yields non aligned points with a flattening for small scales (r 2 = 0.94 for q = 1.5). Same curves for spectral and TM analyses for the Parsivel 2 data are shown in Fig. 2.c and 2.d respectively. It appears that very similar results are found for the Parsivel 2 data. Similar curves for both devices are also obtained for the other studied quantity (, Kdp, Zh, Nt, Dm), therefore in the following only the curves for PWS100 will be showed and discussed in this sub-section.

Very similar results are found for  and Kdp (Fig. 3) with a very good scaling on the whole range of available scales on both the spectra and the TM analysis, which also has to be conducted on the fluctuations of the field. The good scaling behaviour retrieved for Kdp is not surprising and was actually expected if relation 17 is correct. Indeed a power (Kdp here) of a multifractal field (R here) should also behave as a multifractal field (this is the basic concept behind the DTM technique). To the knowledge of the authors, the multifractal behaviour of Kdp has not yet been studied in time with disdrometer data this way, and it opens new perspectives. A potential one would be to compare these outputs with similar analysis performed in space with data provided by weather data. This would enable to study scaling relations in a spatio-temporal framework with a quantity directly measured by the radar (so far the only device providing "rather" high resolution space-time data of rain related fields) 13 without having to rely on tailored relations that may introduce biases in the scaling behaviour as it is the case for the rain rate.

The situation for Zh is more complex. Indeed the energy spectrum (Fig. 4.a) is not linear and could be interpreted as exhibiting two breaks, one at roughly 6 min -1 and the other one (the minimum on the curve) at roughly 2 min -1 . Authors do not have explanation for this behaviour which is not retrieved on R and Kdp. It was not expected and suggests that the "climatic" relation 16 does not hold very well. Indeed if it was true, a good scaling behaviour would be observed on the spectra. To confirm this, the same analysis were carried out on the quantity b aR , where R is the rain rate studied before and a and b are the "climatic" values estimated in section 3.1. As expected for this analysis a good scaling behaviour is indeed retrieved on both the spectra (Fig. of the field (Eq. 15) for which there are a slightly better scaling than for the field itself.

However given the lack of scaling observed in the spectral analysis, the TM analysis might be not very reliable and should not be over-interpreted.

Finally we analysed the scaling features of Nt and Dm which are used to characterize the shape of the DSD. For Nt we find a good spectral behaviour but limited to the range 2 min -2 h, as a flattening of the spectra is observed for small scales (Fig. 5.a). Similarly to the previously studied fields, a spectral slope greater that one is found meaning that the TM analysis should not be performed on the field directly but on its fluctuations (equation 15). It is displayed in Fig. 5.b and it appears that a good scaling behaviour is retrieved on a range of scales (4 min -2 h) similar to the one observed on the spectra (Fig. 5.a). Considering the whole range of scales would lead to a r 2 for q = 1.5 equal to 0.96 whereas is of 0.99 on the limited range of scales. Similar scaling regimes are observed on Dm's spectra (Fig. 6.a), but the slope is much lower meaning that the TM analysis should be conducted directly on the field. It yields a good scaling behaviour with a r 2 for q = 1.5 equal to 0.99 (it is of 0.87 if the TM analysis is conducted on the fluctuations) (Fig. 6.b). The flattening of the spectra (almost a horizontal slope) for small scales corresponds to what would be observed for a white noise. This would suggest that Nt and Dm exhibit a scaling structure down to few minutes and behave as a random homogeneous variable for smaller scales. It is not possible to confirm this interpretation with the help of this data. More data at higher resolution, which would extend the range of available small scales, would be needed to achieve this. A possible interpretation of this could simply be that the sampling uncertainty is more visible at smaller scales and for these quantities than the other ones.

3.3) UM parameters

Estimates of UM parameters for the various studied fields and the two devices are reported in Table 2. It is timely to mention that the proportion of zeros is quite low (12 %, with a fractal dimension of 0.96), which means that estimates are not biased by the multifractal phase transition for small moment orders which is associated with them (see [START_REF] Gires | Influence of the zerorainfall on the assessment of the multifractal parameters[END_REF] for a in-depth analysis of this effect). Given the low quality of the scaling in the spectral analysis for Zh, the values of  and H are not shown for this field since they are not reliable.

For all the fields we find H values greater than 0.5, except for Dm, which confirms that the studied fields are non-conservative, and that the UM parameters  and C1 should indeed be estimated on their fluctuations as it was done. The estimates of the characteristic scaling parameters are very similar for the two devices. It means that despites a roughly 10% difference in terms of rain rate; they both record the same variability across scales. The estimates for the rain rate, which is the field that has been mostly studied, are comparable with the ones usually found for this range of scales by authors focusing the analysis on the rainy portions [START_REF] De Montera | The Effect of Rain-No Rain Intermittency on the Estimation of the Universal Multifractals Model Parameters[END_REF][START_REF] Mandapaka | Multiscaling analysis of high resolution space-time lidar-rainfall[END_REF], Verrier et al. 2010[START_REF] Gires | Development and analysis of a simple model to represent the zero rainfall in a universal multifractal framework[END_REF], although with slightly smaller values of  which were commonly reported to 1.8. With regards to Dm and Zh the differences between the two devices are more pronounced but it is harder to interpret this fact given that the lower quality of the scaling for these fields implies less reliable estimates.

It is possible to give an insight into radar relations 16 and 17 through the UM parameters estimates. Indeed if a field is multifracal, then as previously mentioned a power of it is also multiractal, and there is furthermore a relation between the  and C1 [START_REF] Tessier | Universal Multifractals: theory and observations for rain and clouds[END_REF], Lovejoy et al. 2008) More precisely, if

d dp cK R 
, we have:

dp K R    (18.a) dp K R C d C , 1 , 1   (18.b)
Here for both devices we find a roughly 0.2 difference between  computed for the two fields.

With regards to C1, for the PWS100 data The difficulty here is that the agreement between the estimates of  is worse than for the ). The differences noticed with regards to the underlying theoretical framework highlight the limitation of the "climatic" relations 16 and 17. The better agreement in equations 18 and 19 for the PWS100 data is consistent with the greater quality of the fitting of relations 16 and 17 for the PWS100 than for the Parsivel 2 (Fig. 1). The fact the  are equal for R,  and Nt also suggests that a power-law relation between these quantities could be investigated, which is not surprising since these quantities correspond to various moments of the DSD (Sempere-Torres et al. 2000, [START_REF] Lee | A General Approach to Double-Moment Normalization of Drop Size Distributions[END_REF], Lovejoy and Schertzer 2008, Verrier et al. 2013), even more directly than Zh and Kdp, but this is beyond the scope of this paper.

4) Conclusions

In this paper the output data provided by two optical disdrometers is analysed; the Campbell Scientific PWS100 based on the analysis of the light refracted by drops and the OTT Parsivel 2 based on the analysis of light occluded by drops. Not only the rainfall rate, but also the DSD parameters Nt and Dm, and the radar parameters Zh and Kdp are studied. Furthermore the analysis is not performed only at the maximum resolution as it is commonly done but across scales with the help of the theoretical framework of Universal Multifractals.

It appears that despite a roughly 10 % differences in terms of rain rate, the multifractal analysis yield very similar results for the two measuring devices, with a slightly worse scaling observed on the Parsivel 2 data. R exhibits a very good scaling on the whole range of available scales (30 s -2h) with H=0.7, C1=0.2 and =1.5 which confirms findings of previous studies.

A very good scaling is also retrieve on Kdp and UM parameters estimates are in rather good agreement with what would be expected if the standard relation

  be used in the following sections to help in the interpretation of the observed scaling features exhibited by each quantity. Various authors (Campos and Zawadzki 2000, Jaffrain and Berne 2012b, Verrier et al. 2013) noticed a strong sensitivity of the estimates of a and b on the method implemented to compute them. Here we performed an orthogonal linear regression, which does not assume any dependent variable, on the logs of Zh, R and Kdp. The regressions are performed only for the time steps for which Kdp >10 -2 because we noticed that small values in which we are not interested had a strong influence on the retrieved parameters. Figures 1.a and 1.b display the regressions for the two relations and both devices, and the retrieved values for a, b, c and d are shown in Table

R

  Z h  relation a is slightly greater and b slightly smaller for the PWS100 that for the Parsivel 2 , meaning that the effects of each other are compensating. For the dp K R  relationship c and d are both slightly greater for the PWS100, meaning a given value of Kdp will systematically yield greater estimates of R with the PWS100 values. The estimates of the exponents are in the range of those commonly observed(Jaffrain and Berne 2012b, Figuras I ventura et al. 2013, Verrier et al. 2013). With regards to a values they are in agreement with those found byVerrier et al. (2013), and in the upper range of those reported byJaffrain and Berne (2012b). The c values are similar to those found byJaffrain and Berne (2012b). Relations 16 and 17 are studied only at the maximum resolution (30 s) and the strong 11 scale dependency of the parameters a, b, c and d is not investigated here (see Verrier et al.

  Figure 1.c and 1.ddisplay, for the PWS100 data, a scatter plot for the curves are obtained for the Parsivel 2 data and not shown here. The scattering around the bisector is not negligible with a 30 s time step, and is more that the use of "climatic" values is less acceptable for the former. It should also be noted that there is a tendency of underestimating large rain rate with the than R for R>50 mm/h. The total rainfall depths between the three techniques are very similar with roughly 3% or less differences.

  Fig. 1 is consistent with the fact that we found a good scaling behaviour on R and Kdp and not Zh. Quite surprisingly the scaling breaks observed on the spectra are not visible on the TM analysis (Fig. 4.b) where a unique regime is observed. It was implemented on the fluctuations



  in the relation) which is in rather good agreement with the 0.23 retrieved on the Kdp field. For the Parsivel 2 0.25 retrieved on the Kdp field. The same computations can be carried out for the

  scaling on Zh is not very good and reliable. With regards to C1 we is quite different from the 0.51 found on Zh. Finally it should be mentioned that the estimates of  and C1 on the field PWS b PWS PWS R a (see Fig. 4.c and 4.d for scaling curves) are respectively 1.55 and 0.28 which is roughly in agreement with equation 19

  " parameters computed for this data set. The scaling of Zh is worse, especially on the spectra, which highlights some limitations of the power law relation only holds on the range few min -2 h for Nt and Dm, and there are some hints at a possible random uniform behaviour for smaller scales possibly associated with sampling uncertainty.These results suggest new ways of comparing the outputs of disdrometers by using other fields than the rain rate, and also scaling analysis. The results are particularly promising for Kdp which is also directly measured by polarimetric weather radars, contrarily to the rain rate for which non-trivial transformations potentially biasing the observed scaling are implemented. Multifractal investigations in a spatio-temporal framework on Kdp radar data should be carried out to improve knowledge about rainfall as a space-time process. More data, including spatial ones, should be analysed to confirm the possibility of characterizing DSD parameters with the help of UM, and hence the possibility of developing coupled multifractal cascades to actually simulate DSD fields.

Figure

  Figure caption list:

Figure 1 :Figure 2 :Figure 3 :Figure 5 :

 1235 Figure 1: Computation of the parameters of the

Table 2 :

 2 Scaling parameters of the various studied fields for the PWS100 and Parsivel 2 data

	Table:
						PWS100	Parsivel 2
	Z h 	R	a	349	313
					b	1.49	1.63
	R 	K	dp	c	17.3	15.3
					d	0.72	0.66
	Table 1: Parameters computed for the	R Z h  relation (Equation 16) and the	R 	K	dp
	relation (Equation 17)
							23

Acknowledgments

Authors thanks Serge Botton (from the "Département Positionnement Terrestre et Spatial" of the ENSG) for facilitating access to the roof where the disdrometers are installed. Authors greatly acknowledge partial financial support form the Chair "Hydrology for Resilient Cities" (sponsored by Veolia) of Ecole des Ponts ParisTech, EU NEW INTERREG IV RainGain Project (www.raingain.eu), and EU Climate KIC Blue Green Dream (www.bgd.org.uk)