Giuseppe Iurato
email: giuseppe.iurato@unipa.it

A simple phylogenetic remark about certain human blood pressure values

Keywords: evolutionary physiology, humans, life history evolution, phylogenetics

An evolutionary biopsychology hypothesis on some features of the physiological phenomenon of human blood pressure, is put forward.

1. Arterial blood pressure (in short, BP) is the physiological pressure the blood exerts against arterial's walls as it circulates through the human body. It is mainly ruled by a variety of physiological systems, amongst which are neural and hormonal signals from the heart, vasculature, brain, kidneys, and gastrointestinal organs. BP is generally estimated by measuring systolic and diastolic components. Systolic blood pressure (in short, SBP) is the maximal pressure in blood vessels during cardiac systole (heart contraction) and diastolic blood pressure (in short, DBP) is the minimal pressure in blood vessels during cardiac diastole (heart relaxation between contractions). BP is most commonly measured peripherally in the upper arm [START_REF] Piper | Screening for High Blood Pressure in Adults: A Systematic Evidence Review for the U.S. Preventive Services Task Force[END_REF]).

2. The circadian BP variation, although affected by a variety of external factors, represents the influence of internal factors such as ethnicity, gender, autonomic nervous system tone, vasoactive hormones, and hematologic and renal variables. In most individuals, BP presents a morning increase, a small post-prandial valley, and a deeper descent during nocturnal rest [START_REF] Hermida | Circadian variation of blood pressure: The basis for the chronotherapy of hypertension[END_REF].

3. An experimental datum of human physiology (Carretero & Oparil, 2000a,b;[START_REF] Chugh | Circadian variations in blood pressure in health and disease: implications for patient management[END_REF][START_REF] Cohen | Meta-analysis of the morning excess of acute myocardial infarction and sudden cardiac death[END_REF][START_REF] Elliott | Circadian variation in the timing of stroke onset: a meta-analysis[END_REF][START_REF] Hassler | Circadian variations in blood pressure: implications for chronotherapeutics[END_REF][START_REF] Hermida | Circadian variation of blood pressure: The basis for the chronotherapy of hypertension[END_REF][START_REF] Kario | Morning surge in blood pressure as a predictor of silent and clinical cerebrovascular disease in elderly hypertensives: a prospective study[END_REF][START_REF] Lacourciere | Sustained antihypertensive activity of telmisartan compared with valsartan[END_REF][START_REF] Manfredini | Morning surge in blood pressure as a predictor of silent and clinical cerebrovascular disease in elderly hypertensive: a perspective study[END_REF][START_REF] Mansoor | Circadian Blood Pressure Variation in Hypertensive Patients With Primary Hyperaldosteronism[END_REF][START_REF] Millar-Craig | Circadian Variation of Blood-Pressure[END_REF][START_REF] O'brien | Working Group on Blood Pressure Monitoring of the European Society of Hypertension International Protocol for validation of blood pressure measuring devices in adults[END_REF][START_REF] Piper | Screening for High Blood Pressure in Adults: A Systematic Evidence Review for the U.S. Preventive Services Task Force[END_REF][START_REF] Redon | Prognostic Value of Ambulatory Blood Pressure Monitoring in Refractory Hypertension. A Prospective Study[END_REF][START_REF] Redon | Uncontrolled early morning blood pressure in medicated patients: the ACAMPA study. Analysis of the Control of Blood Pressure using Ambulatory Blood Pressure Monitoring[END_REF][START_REF] Smolensky | Chronobiology and chronotherapeutics in primary care: patient care[END_REF][START_REF] Weber | The 24-hour blood pressure pattern: does it have implications for morbidity and mortality?[END_REF] regarding circadian measurement of blood pressure, is that its highest SBP and DBP values are reached around sunrise (say SB-DB p sr ) and sundown (say SB-DB p sd ). 4. To be precise, blood pressure rises sharply in the early morning in response to the activation of the sympathetic nervous system when one awakes. This early morning surge is associated with other important hemodynamic and neurohormonal changes, including increases in heart rate, vascular tone and blood viscosity, and decreases in vagal activity. The activity of the sympathetic nervous system is quiescent during sleep, whereas awakening selectively increases epinephrine levels. Subsequent increases in BP and heart rate are controlled by direct sympathetic neural input into the heart and vasculature in response to increases in activity and upright posture, rather than by an endogenous surge of plasma catecholamines [START_REF] White | Importance of Blood Pressure Control Over a 24-Hour Period[END_REF].

5. There is no a proper, satisfying causal explanation to this fact. It is a mere empirical fact, related to the activation of the sympathetic nervous system, to neurohumoral and hormonal variations, but without a valid, effective motivation why this latter takes place. Above all, the occurrence of the sunrise peak, which is the most prominent one out of the two, is truly incomprehensible because one might rightly argue, on the basis of further experimental data [START_REF] Millar-Craig | Circadian Variation of Blood-Pressure[END_REF][START_REF] Mansoor | Circadian Blood Pressure Variation in Hypertensive Patients With Primary Hyperaldosteronism[END_REF], that, after a night of sleep relaxation, this blood pressure value should instead be the lowest rather than the highest.

6. One1 of the earliest descriptions of circadian rhythm dates back to 1729 by Jean-Jacques D'Ortous De Mairan, who described diurnal leaf movements. Nearly three centuries later, having discovered so much regarding circadian biology, we are still searching for answers to the crucial question ''What makes the biological clock tick?'', to which circadian phenomena are closely related. From a clinical standpoint, the discovery of circadian influence on heart rate and blood pressure was parallel to the technological advances in blood pressure monitoring. There exist biophysiological factors governing circadian blood pressure rhythm.

7. In mammals, the mechanism controlling rhythmic biophysiological processes is located in the suprachiasmatic nuclei of the hypothalamus, one of the most ancestral cerebral zone strictly correlated with autonomous nervous system. This group of approximately 20,000 neurons is referred to as the central clock. These cells are autonomous, as they rhythmically express gene products regardless of external influence. Synchronicity feature of neurons within the suprachiasmatic nuclei is achieved through various coupling mechanisms, even in the absence of light. Furthermore, the discovery of autonomous oscillators in other regions of the mammalian brain and peripheral tissues, challenges the paradigm that there is only one major clock coordinating biophysiological rhythms. 8. However, [START_REF] Chugh | Circadian variations in blood pressure in health and disease: implications for patient management[END_REF] stresses on the suprachiasmatic nuclei as the central clock, synchronized each day by exposure to light, coordinating rhythmic variations in biophysiological functions. It is for this main reason that the suprachiasmatic nuclei is considered as a kind of, so to speak, 'pacemaker' of circadian biorhythms. Tight control of genetic expression within these pacemaker cells corresponds to the maintenance of circadian schedules. Some transcription factors are known to be key players in the timing mechanism of the suprachiasmatic nuclei. Under the influence of these factors, for instance multiple protein expression changes occur, ultimately affecting neurohormonal discharge. Various studies exploring the effects of core clock gene mutation have clearly shown a strong association with maintenance of a 24-hour cycle. Activation of clock-controlled genes occurs when clock proteins bind to promoters, enhancement regions, and even non-coding sequences of DNA. 9. Studies regarding clock-controlled gene expression have shown that vasopressin concentrations within the suprachiasmatic nuclei and hepatic production of the steroid Cyp2α5, are under clock gene influence. Using gene microarray analysis, several groups have demonstrated that up to 10% of gene expression, is under circadian control. Furthermore, clock-controlled genes expressed as a result of this influence vary across tissue types. It is important to note that the research supporting these findings has been conducted with the use of animal models, this being methodologically licit because we are studying the most atavistic zones of the brain. 10. Correlative studies in human subjects are accordingly challenging, but new methodologies have been described. In this regard, clinically, the ultimate concerns posed by these molecular details are the end-result effects on circadian blood pressure variations. Multiple neuroanatomical studies have shown that suprachiasmatic nuclei output is mainly restricted to the medial hypothalamus, although there is evidence supporting termination of suprachiasmatic nuclei efferents in the rostral area as well as the thalamus of mammalian brains.

11. Ultimate targets of suprachiasmatic nuclei stimulation are endocrine neurons and preautonomic neurons. Moreover, direct connections to hormone-releasing neurons (e.g., corticotrophinreleasing hormone, thyrotropin-releasing hormone, gonadotropin hormone-releasing hormone) have been shown, as well as the vasopressinergic effect of the suprachiasmatic nuclei exerting on the releasing of glucocorticoids via the hypothalamicpituitary-adrenal axis has been identified. While these neural connections may indicate a level of control on multiple biophysiological processes, they do not encompass the magnitude of the effect observed. Signaling of preautonomic neurons influences a wide amount of biophysiological functions.

12. Studies using transneuronal tracers have identified the paraventricular nucleus as the center mediating signals from the suprachiasmatic nuclei to the periphery by means of sympathetic and parasympathetic neurons. Perhaps the best described examples of autonomic control exerted by the suprachiasmatic nuclei are that of melatonin secretion and hepatic control of glucose metabolism. Hence, current conventional thought is that these mechanisms play more than an ancillary role in the changes seen in the cardiovascular system through the course of a day, with particular attention to circadian blood pressure value range which, in any case, is closely related to these ancestral zones of human brain [START_REF] Chugh | Circadian variations in blood pressure in health and disease: implications for patient management[END_REF].

13. In particular, just in regard to the occurrence of the characteristic pair (SB-DB p sr , SB-DB p sd ) of blood pressure value range, beyond their biophysiological roots (as seen above) which might explain them only as the result of a mere biomechanics causal effect, we would like to propose a simply historical epigenetic motivation for explaining this fact, precisely based on a phylogenetic origin of the occurrence of such typical values. That is to say, we would like to surmise, within evolutionary biopsychology context, that such typical highest values might have been due to the great and uncontrollable fear of primitive humans when sunrise and sundown occur, a fear mainly related with the high probability of aggressions by wild beasts or, anyway, to the immeasurable existential anguish that surely distinguished primitive human life. This fear should therefore have been inherited epigenetically, as clearly corroborated by what has been said at the points 6.-12., since circadian blood pressure values are ruled by biological clock which is just placed into those most ancestral cerebral zones of human brain.

14. Furthermore, since who suffers of essential hypertension has the pressure values SB-DB p sr and SB-DB p sd very high and greater than standard fixed values to be assigned for normotensive patients, to support this our hypothesis, we may invoke the statistical datum [START_REF] Chobanian | Seventh report of the Joint National Committee on Prevention, Detection, Evaluation, and Treatment of High Blood Pressure[END_REF][START_REF] Heidenriech | Forecasting the future of cardiovascular disease in the United States: a policy statement from the American Heart Association[END_REF][START_REF] Nwankwo | Hypertension among adults in the US: National Health and Nutrition Examination Survey, 2011-2012[END_REF] according to which who mostly suffers of hypertension among world population, is just African people living from generations in USA.

The points 6.-12. are drawn, almost verbatim, from(Chugh et al.

2011).


A Emilien
email: arnaud.emilien@inria.fr

† P Poulin
email: poulin@iro.umontreal.ca

M.-P Cani

§ U Vimont

Interactive Procedural Modeling of Coherent Waterfall Scenes

Keywords: Descriptors (according to ACM CCS): I, 3, 5 [Computer Graphics]: Computational Geometry and Object Modeling-Physically based modeling

Figure 1: An example of a waterfall scene, the Iron Hole, on the island of La Réunion. Left: Photo of the real site c Serge Gélabert. Center: Our result after 10 minutes of interactive procedural modeling, starting from a similar terrain model. The scene contains 36 elements interconnected by pools and rivers, deforming the terrain and controlling the flow. Right: Visualization of the control elements that we used to create the waterfall network.

Introduction

Procedural modeling is a great paradigm for automatic modeling of complex objects and scenes, and it has been applied to many problems including terrains, trees, buildings, streets, and cities. In addition to its efficiency for generating complex details, its power lies in its ability to ensure that certain submitted to COMPUTER GRAPHICS Forum (12/2014). constraints are respected (for instance from a physical, biological, or architectural viewpoint), making it much easier for the user to create plausible models. However, when one has a very specific goal in mind, controlling the many intricate parameters of these automatic procedures can be very cumbersome. In such situations, it would be great to have access to an interactive system enabling the user to handle the coarse design, while the system automatically ensures consistency and support for the procedural generation of all the details. In this paper, we apply this paradigm of interconnected procedural generation and interactive user-control, to the difficult task of modeling coherent waterfall scenes.

Waterfalls offer some of the most beautiful settings in nature. Despite this, no easy-to-use method for designing waterfall scenes has been developed so far in computer graphics. One solution consists of using physically based simulation of fluids on top of a modified terrain. Unfortunately, modeling a terrain in order to produce specific waterfalls is an extremely daunting task. Fluid simulation depends on slopes, collisions, riverbeds, source influx, water properties, flow speed, etc. All the work of deforming a terrain does not ensure the final appearance of a waterfall, as these mutually dependent constraints and the governing physical laws are complex and highly nonlinear. Another potential solution consists of manually creating a few waterfalls with standard modeling tools, using manifold meshes and/or particles, and then positioning them along the terrain and connecting them by streams. In this long and tedious process, the artist also needs to manually maintain the consistency between the terrain and waterfalls, as well as the self-consistency of the waterfall network.

In this paper, we combine interactive and procedural methods to enable fast and easy design of plausible waterfall scenes. Our solution, based on a new interactive procedural model for flowing-water networks, allows users to easily shape complex waterfall scenes while automatically enforcing the physical consistency of the results, both in terms of hydraulic flow and of plausible embedding into the terrain. Coherence is enforced at several levels, ensuring that the flow always goes downstream, the user-edited water network follows the terrain (or the terrain is adapted to the network), the flow is distributed in the network without loss or gain of volume, and that the type of waterfall and the appearance of water is adapted to the flow. Our main contributions include:

• a slope-flow diagram-based classification of waterfalls;

• three parametric models for designing waterfall elements; • a procedural method for ensuring waterfall network consistency; • automatic methods for locally adapting user-input water trajectories to the terrain and/or the terrain to the flow.

All these contributions combine into a framework allowing an artistic approach to river and waterfall design. The resulting scene could either be used as a synthetic environment for games or films, or as an initial setup for further refinement through physically based fluid simulation when more sophisticated simulation effects are needed.

Previous Work

Many types of procedural methods [START_REF] Smelik R | A survey on procedural modelling for virtual worlds[END_REF] have been designed for specific types of objects, including terrains [GGG * 13], plants [LRBP12], buildings [LWW08], cities [CEW * 08], road networks [GPGB11], and villages [EBP * 12]. To our knowledge, no previous work addresses the procedural modeling of waterfall scenes. We therefore focus our literature review on the modeling of rivers, waterfalls, and terrains, which are key elements for our problem.

Rivers and Waterfalls

Most previous methods for creating running-water elements, from rivers to waterfalls, consist in simulating shallow water flows over an existing terrain [KM90, HW04, KW06, TMFSG07, YNBH09, LH10, BSW10, CM10]. In particular, particle systems are heavily used for generating waterfalls. However, even when optimized with hierarchical or screenspace methods [BSW10], these approaches suffer from the inherent complexity of simulating networks of running water over large environments. Moreover, they only provide indirect control on the nature of waterfalls and on the trajectories of streams, the latter being dictated by the geometry of the underlying terrain.

To increase user control over complex waterfalls, Sakaguchi et al. [SDZ * 07] introduce pre-tuned particle systems designed to capture specific types of waterfalls. These systems are assembled by the user to generate visually complex scenes. Another solution [GCZ * 06] directly creates waterfall scenes by placing polygonal primitives over the terrain, and renders them using animated textures. However, in both cases, the consistency within the flow network and between the terrain and waterfall elements needs to be manually ensured by the user.

As with these two last approaches, we move away from pure simulation methods in order to give more control to the user. However, we do so while maintaining the main benefits of physically based simulation, namely automatically ensuring a coherent result. This includes both that flow strengths are coherent between the elements of the waterfall network and that each flow segment is geometrically consistent with respect to the underlying terrain. The way our user-designed vector elements control flow trajectories is similar in spirit to the sketch-based approach from Bhat et al. [START_REF] Bhat | Flowbased video synthesis and editing[END_REF] for editing videos of 2D flows, although their method is based on video processing only and does not address the problem of plausibly embedding water flows into a 3D terrain.

submitted to COMPUTER GRAPHICS Forum (12/2014).

Although we do not use their simulation methods, we rely on rendering techniques introduced in previous work to display the aspects of waterfall scenes in real time: Yu et al. [START_REF] Neyret | Scalable real-time animation of rivers[END_REF] derive an animated texture from the motion of particles simulated at the surface of a river, which effectively provides the appearance of a complex fluid. These textures can be augmented with their flow skirting around stones and borders [START_REF] Neyret | Feature-based vector simulation of water waves[END_REF]. Instead, Van Hoesel [vH11] tiles flow textures and modulates their application on the water surface polygons according to the flow speed. His method proves to be quite efficient, compact, and effective. Zhu et al. [ZIH * 11] introduce an interactive flow and diffusion editor with a sketching interface. The flow textures we use combine ideas from all of these approaches, and are enhanced by simple forms of particle systems, inspired by work from Holmberg and Wünsche [START_REF] Holmberg | Efficient modeling and rendering of turbulent water over natural terrain[END_REF].

Terrains

Allowing users to design a specific waterfall scene over some input terrain requires automatic ways of adapting the terrain in order to ensure consistency with the flow. We therefore briefly review previous methods for terrain generation and editing.

Most terrain models are represented using a 2.5D heightfield that is well adapted to GPU processing, enabling the handling of very large environments [LH04, BN08]. In contrast, Peytavie et al. [PGGM09] introduce a specific piling model for capturing complex terrain structures, such as overhangs and caves that can be observed at the top of free-falls or near pools at their bottom. In our work, we instead adapt a horizontal displacement method [START_REF] Gamito | Procedural landscapes with overhangs[END_REF] to capture overhangs over a standard heightfield representation.

Many solutions exist for terrain generation and editing [SDKT * 09], from fully procedural methods to those combining examples or textures with sketch-based interaction. We focus on the latter, since they could be adapted to control some local deformations of the terrain.

Hydraulic erosion methods [BTHB06,MDH07,ŠBBK08]

simulate the effect of water flows on the progressive sculpting of terrains, which results in very realistic landscapes. However, here the control is indirect since the flow depends on the current terrain. Therefore, this approach is not suitable for our problem. Closer to our concerns, Génevaux et al. [GGG * 13] present a procedural method for generating terrains based on hydrology: a realistic terrain is fully generated from dense hydraulic graphs. Unfortunately, their approach does not directly address our goal: we are instead looking for plausible ways to locally edit an existing terrain in order to make it consistent with user-designed flow networks.

Re-using existing terrains to match user specifications, such as sketches, is generally accomplished using texturebased approaches, where heightfield patches from the input terrain are combined to match user constraints [START_REF] Zhou | Terrain synthesis from digital elevation models[END_REF][START_REF] Tasse | Enhanced texturebased terrain synthesis on graphics hardware[END_REF]. We are rather interested in local, feature-based editing of the input terrain, to have it match the range of slopes and local geometry consistently with user-designed flow networks. We are therefore inspired by feature-based terrain generation methods [HGA * 10, BMV * 11, GMS09], and we adapt them for the first time to the editing of an existing heightfield.

Overview of our Method

The key goal of our method is to leave coarse design of waterfall scenes in the hands of the user, while providing automatic ways to generate plausible and detailed results.

In the real world, flow networks formed by complex waterfalls, such as the one in Figure 1(left), include free-falls, segments where running water remains in contact with the terrain, as well as pools. Moreover, each of these segments can be of many different types, from rivers to rapids for the ones in contact with the terrain, and from plunges to cataracts for free-falls. Having to manually select plausible types for each segment of a full network would be both tedious and require specialized knowledge from the user.

Below we propose a two-level classification for segments of waterfall networks, enabling us to leave the choice of lowlevel classes to the user, while automatically computing the most appropriate running-water types from quantitative information, such as the slope of the underlying terrain and the intensity of the flow. The processing pipeline that we use for modeling a waterfall scene, based on this analysis, is presented next.

Two-level Classification for Running Water

Waterfall scenes are comprised of three types of elements: running-water segments that remain in contact with the terrain, free-fall segments where water is in the air, and pools that receive water from the free-falls. Our goal is to provide some coarse, intuitive control to the user, and we leave the choice of these three classes, contact, free-fall, and pool, to the user during interactive design.

In contrast, we would like to free the user from the manual and explicit determination of the precise characteristics that each running-water segment should take, because this can be determined in a more plausible way using an automatic procedure. After studying the existing classifications of streams and falls, we designed a new, slope-flow classification that matches our goals, as explained next.

Running water can take on many forms, from rivers to rapids, and from plunges to cataracts. See Figure 2 for an illustration of these classes. Several classifications of waterfalls are proposed by hydrologists, geologists, and artists, in order to capture this variety: • volume-based classifications of waterfalls [START_REF] Beisel | International Waterfall Classification System[END_REF] sort waterfalls into classes by using a logarithmic scale over the volume of water in the air at a given time. Although easy to compute from quantitative information, this classification provides little clue on the visual aspect of the fall. Moreover, it is restricted to free-falls, and therefore does not fully meet our needs; • geometric classifications, such as the one found in the Waterfall Lover's Guides [DD06, Plu05] or the one depicted in Figure 2, analyze the different types of geometries that can be observed in nature. However, they only provide visual information. No quantitative measurement is proposed to automatically compute the class that a runningwater segment belongs to.

In this work, we would like to classify waterfall segments from quantitative information, while getting visual clues enabling us to generate plausible 3D representations for each segment. We therefore decided to augment the geometric classification of Figure 2 with quantitative evaluation of the classes, as in volume-based approaches. However, we need measures applicable to both running water in contact with the terrain and to free-falls, and therefore our solution is different.

By studying the existing geometric classifications and looking at many real cases, we noticed that the geometric type of running water mostly depends on two important quantitative parameters: the flow value (defined as the volume of water per second traveling through a cross-section of the segment), and the local slope of the terrain. Intuitively, when the flow decreases, a river becomes a stream, a cascade becomes a horsetail, and a free-fall cataract becomes a ledge. Meanwhile, if a given water flow is running on terrains of increasing slope, a river tends to become a rapid, and then eventually a block.

To provide a quantitative classification, we define the different classes of running-water segments as regions in a unified slope-flow diagram. This is done as follows: we first used real examples to find a consistent set of seed values for typical elements of each of the visual types depicted in Figure 2. These representative slope and flow values are listed in the left columns of Table 1. The regions associated with each class are then defined using a qualitative Voronoi segmentation in slope-flow space. The resulting classification is depicted in Figure 3. Note that the types of water flows on the left of the division in red belong to contact waterfall seg-ments, while those on the right belong to free-fall waterfall segments.

We also used the visual classifications serving as reference to associate visual parameters with each type of elements in our classes, given in the columns to the right of Table 1. These values will be used for generating the appropriate visual representation for each element of the waterfall scene.

Processing Pipeline

In the remainder of this paper, we define a waterfall network as a network of waterfall segments and of pools. Waterfall segments can either be of type contact (in contact with the terrain) or of type free-fall (in the air).

Our processing pipeline, based on the two-level classification we just defined, develops as follows (see Figure 4):

1. The user starts by creating a waterfall scene, building an oriented vectorial controller network U over an existing terrain. This is done using three vectorial controllers (contact, free-fall, and pool)all based on Cardinal splines. Contact and pool are interactively created by using control points while free-falls are parabolas, automatically parameterized by their start and end points. During interaction, controllers are not constrained to be placed in contact with the terrain (since the terrain will be adapted later according to the user design), but the flow is constrained to go downhill along each segment of the controller network. While the default flow intensities on each segment of the network can be interactively tuned by the user, a submitted to COMPUTER GRAPHICS Forum (12/2014). give the position of the class seed in the slope-flow diagram of Figure 3. Slope is an average inclination in radians, while flow is expressed in "flow units", which gives an informal notion of relative proportions between waterfall types. Foam, rocks, disturbance, and particle are parameters (∈ [0, 1]) used in our procedural generation of geometry and for rendering.

consistent hydraulic graph G (with fully consistent flow values) is automatically computed at the end of the interactive modeling process. (See Section 4) 2. The next step is the generation of the waterfall network W, which uses the coarse water trajectories from the controller network and the flow information from the hydraulic graph to define a more precise representation of the waterfall. In addition to directly using the control curves that they defined, the user has the option of further refining the geometry of the network through an automatic procedure that locally adapts running-water trajectories to the underlying terrain. Each curve of the network is then divided into a number of water segments, whose sub-class is determined by slope-flow classification from Section 3.1. Finally, the last geometric parameters, such as flow width and depth, are computed for each segment of the waterfall network. (See Section 5) 3. Lastly, the 3D representation for the waterfall network, called the integration mesh, is generated and embedded into the scene through appropriate local deformations of the terrain. Although they can include large changes, such as digging a canyon to allow a stream to find its way downhill (in the extreme case where the user designed water segments go through a mountain), constraints applied to the terrain are mainly aimed at automatically adding all the details that make the scene plausible: this includes borders along streams, riverbeds, and the creation of overhangs behind free-falls. Appropriate decorative elements such as trees and rocks are generated at this stage, using the distance to the closest riverbed and the class of the waterfall segment it belongs to as de- Note that we chose to compute a coherent hydraulic graph G (Step 1 of the pipeline above) based on the controller network U, i.e., from the coarse trajectories defined by the user, before refining these trajectories into a waterfall network W. This enables us to use consistent flow values in Step 2, while refining water trajectories. This helps us, for instance, to prevent large rivers from being refined into a series of short twists when adapted to the terrain, although a smaller stream would be allowed to be more winding. The fact that our algorithm interleaves geometric comsubmitted to COMPUTER GRAPHICS Forum (12/2014). p p j,k-1 j,k smin Figure 5: Minimum slope constraints are imposed on each curve of the controller network, to make it consistent with the user-defined flow direction. Left: A controller network. Right: The algorithm verifies that the minimum slope s min is respected between the contact control points p j,k-1 and p j,k . Since this is not the case, p j,k is lowered. All the subsequent controller points will also respect this constraint. In this figure, the pool control points are therefore lowered.

putation with consistency checks leads to more plausible results. The remainder of this paper details the three aforementioned steps of the pipeline.

From User Input to a Coherent Hydraulic Graph

We will discuss how the user creates the controller network U (Section 4.1), how we generate the associated hydraulic graph G, and how we compute its flow (Section 4.2).

Controller Network Creation

The user builds a controller network U (Figure 4) by creating and manipulating different vectorial elements V i and by interconnecting them.

We define controllers using a Cardinal spline with a controller type α i ∈ { free-fall, contact, pool }. The splines are composed of a series of control points ∆ i = {p i,k }, which can be connected to control points from other controllers to create a controller network (Figure 4 Step 1). Depending on α i , controllers are set differently: a controller contact is created by positioning a series of control points p over the terrain. This controller is used to create all the elements that remain in contact with the ground, such as rivers. The pool controller has been introduced to create flat polygons, and is used to create pool contours. Because a controller pool is flat, the user first positions a horizontal plane in space, and then traces a contour over the plane using control points. Finally, the free-fall controller is used to create an element that loses contact with the ground. Consequently, only a start point and an end point can be positionned on the terrain, with the lower point being typically placed inside a pool. The user associates a flow direction to free-fall and contact curves.

During editing, point and curve magnets are used to facilitate the interactive creation of the connections between elements, like in most classical vectorial editing tools. Moreover, the user can insert, delete, and move control points on each curve, and cut or merge controller curves.

Projection. To facilitate the adaptation of the controller network to the underlying terrain, users are provided with a pro-jection tool, which can be used to project control points of the Cardinal splines onto the terrain. If desired, they can also define an offset from the terrain's local height for each control point. When the projection tool is used for pools, which are constrained to be flat, the pool level is automatically set to the average height of all projected contour points.

Minimum slope. For the scene to remain consistent, the system must ensure that all slopes set to the controller network segments allow the water to flow downstream in the user-defined direction. This is accomplished using automatic correction of the user-defined positions, during a traversal of the controller network (see Figure 5): starting with the flow sources of the network, and traversing it in a topologically sorted order, we check if each control point p j,k of curve ∆ j of each controller maintains the minimum slope s min with respect to its predecessor p j,k-1 . If not, point p j,k is lowered down to match the constraint, i.e., for each connection, we check that all the outgoing nodes are lower than the incoming nodes, and update their positions if needed. Pool contour points are lowered, if necessary, according to the full set of free-fall segments coming into the pool.

Hydraulic Graph Generation

The hydraulic graph G (Figure 4 Step 2) is an oriented graph, generated from the controller network U , in which the flow originates from the sources and exits by the wells. It is composed by a set of nodes N j = (β j , γ j ) and a set of arcs

A k = ((N i , N j ), η k , γ k )
, where β j ∈ { source, well, branch, pool } is its type, γ j is the flow going through the arc or the node, and η k ∈ { contact, free-fall } is the arc type.

For each controller V i we create an arc if its type α i is free-fall or contact, and create a node if α i is pool. Nodes for sources and wells are introduced at the extremities of the graph, branches are created at the intersections between controllers.

Flow Computation

It would be time-consuming and non-intuitive in cases of failure to perform a complex physical simulation on the flow propagation in order to deduce the entire flow properties (intensity, speed, etc.) everywhere along the network. We were inspired by the interactive system of Zhu et al. [ZIH * 11], and therefore simplified our hydraulic model by considering it as a "pipe-like" graph. This enables us to deal only, at this stage, with flow exchanges at the nodes of the hydraulic graph. In the remainder of this section we detail how to compute these flow exchanges using simple and intuitive functions, while providing a coherent visual appearance.

Flow consistency. A hydraulic node should be at equilibrium, i.e., its incoming flow should equal its outgoing flow. Expressing this at each node of the hydraulic graph leads to a system of interconnected equations. There is generally submitted to COMPUTER GRAPHICS Forum (12/2014).

an infinite number of solutions to this system. Therefore, instead of using a global solver to find consistent flow values, we solve them in succession for each node of the graph, exploring it in dependency order from sources to wells. This enables us to take into account the user specifications for the relative strength of the input flows, and to generate a solution that best matches the coarse geometric trajectories in terms of branching angles at each node. Our method for doing so is explained next.

Separating incoming flow into outgoing branches. The flow of an outgoing arc at a branch node N i should depend on the angles between the inflow and outflow arcs in order to capture the natural course of water. For instance, we expect that most of the flow should follow its own original direction.

Let X j be an input arc of N i , and u j its incoming direction (Figure 6). We distribute its flow γ j to each output arc Y k with a direction v k according to a normalized weight:

w jk = 1 + (u j • v k ) w jk = w jk ∑ k w jk .
Thus, the outgoing flow γ k for arc Y k is computed as Separating the flow out of pools. Because the shape of a pool could be complex, distributing the flow according only to the angles between inflows and outflows would not be realistic, while computing a full simulation would be computationally expensive and not suitable for an interactive system [ZIH * 11]. Instead, we simply distribute inflows equally to all the outgoing arcs, except when we need to take the relative flow strength pre-set by the user into account for these branches.

γ k = ∑ j γ j w jk . u v 0 2 v 0 v 1

Waterfall Network Generation

In this section, we discuss the generation of a waterfall network W from the controller network U , which contains the waterfalls coarse trajectories, and the hydraulic graph G, which contains the flow information. The waterfall network is composed of waterfall segments S i , interconnected by waterfall nodes B i .

We start the construction process by subdividing the controller trajectories and locally adapting them to the terrain (Section 5.1). Then, we extract from them the waterfall segments and nodes that form the waterfall network (Section 5.2). We compute the type of each segment using the classification from Section 3.1, and set its other parameters while taking flow and slope into account.

Subdivision and Adaptation to Terrain

In order to procedurally improve the curves created by the user, we provide a fractal-like subdivision scheme based on midpoint displacement, which takes the underlying terrain model into account (see Figure 7). Except for very intricate terrains, our subdivision scheme approximates the natural trajectory of a flow running down a slope while preserving the global shape of the user-defined curve. Let A j be an arc of the graph G, and V i its corresponding controller, with ∆ the points of the curve manually created by the user, and γ the flow going through the arc. Consider the segment formed by two consecutive points p i and p i+1 . We subdivide this segment at its middle point m and move it along direction u, perpendicular to this segment, horizontally. This creates two new segments, V and W , represented by normalized vectors v and w. While the original midpoint displacement method applies a random displacement to m with a maximum amplitude τ, we instead search for x = m + λu, λ ∈ [-τ, τ], minimizing the cost function:

C(x) = wgCg(x) + waCa(x) + wrCr(x)
where Cg(x) is the gradient cost, Ca(x) the angle cost, and Cr(x) the random cost. The values wg, wa, and wr are weight coefficients associated to the costs.

Gradient cost. The gradient cost favors paths that follow the slope of the terrain, and is defined as:

Cg(x) = V g(t)) fp(v, g(t))dt + W g(t)) fp(w, g(t))dt p i+1 -p i , fp(v, g(t)) = -2v • g(t) g(t) + 1,
where vector g(t) is the gradient of the terrain elevation at position t. By integrating the elevation gradient along segments V and W , and by using the scalar product between the gradient and the segment vectors (v or w) as a penalty coefficient, the path will follow the slope and avoid obstacles (Figure 7). Moreover, instead of simply using the scalar product value as a penality coefficient, we use function fp to penalize paths that follow the isolines (i.e., when the scalar product with the gradient is zero); climbing a slope is penalized even more severly. Consequently, this cost is low when v and w are aligned in the same direction as the slope of the terrain, and high otherwise.

Angle cost. The angle cost prevents undesirable sharp angles that could appear between two consecutive segments, because of their independent subdivisions. We take into account the angles a, b, and c (in radians), created at the introduction of new point x (see Figure 7), as:

Ca(x) = a π 2 + b π 2 + c π 2 .
Random cost. Finally, to add fractal-like details on flat terrains, we use the random cost Cr(x), which is negligible when the other costs are high.

The segments are recursively subdivided until their length is inferior to l. Displacement amplitude τ and detail size l are set as:

τ = u /2 l = γ/σ.
The minimum subdivision length l is proportional to the segment flow, where σ is a user parameter. This enables us to get a more detailed trajectory for small flow values, enabling streams to become more winding than rivers. In our prototype we use the following values: σ = 1/2, wg = 1, wa = 0.1, and wr = 0.2 for all our examples.

Waterfall Network Construction

The waterfall network is composed of waterfall segments

S i = (u i , γ i , κ i , δ i , ε i , ζ i )
where u i is the segment vector, γ i is the flow going through the segment, κ i ∈ { stream, horsetail, cascade, rapid, block, river, ribbon, plunge, ledge, cataract } is the waterfall type according to our classification (Figure 3), δ i is the speed of the flow, and ε i and ζ i are respectively the width and depth of the riverbed. These segments are interconnected by waterfall nodes B j = (µ j , γ j , ζ j )

Figure 8: Triangular cross section of a waterfall segment, with a constant flow.

where µ j ∈ { source, well, branch, pool } is the node type, γ j the total incoming flow, and ζ j the depth.

Each segment vector u i is directly extracted from the subdivided trajectories (Figure 4), and its flow γ i is equal to the flow of its corresponding graph arc. All consecutive segments are connected by a branch node with only one input and one output. The other segments are connected by the waterfall nodes B j constructed from the hydraulic graph nodes N k and their associated controller V l . The remaining parts of this section describe how the other segment properties are computed.

Waterfall Segment Type

For each waterfall segment S i , we know its slope s i and its flow γ i . These values automatically determine the waterfall segment's type κ i by casting its coordinates (s i , γ i ) in the slope-flow graph of Figure 3. We use Voronoi cells, whose coordinates are detailed in Table 1, to determine to which type the coordinates belong.

Note that if the user did not use a plausible waterfall controller, e.g., if he created a free-fall on a flat terrain, or a contact on a very steep terrain, the waterfall segment may be outside of the valid range of values in the slope-flow graph. In this case, the parametric model is still set using the closest type, but the user is notified (i.e., the related segment is drawn in red). He can then either validate the current design (even if it is not fully realistic), or select more realistic controllers.

Waterfall Segment Properties

The final step in the generation of the waterfall network is to compute the last properties of each segment S i , i.e., its speed δ, its riverbed width ε and depth ζ, from the slope and flow information we have. We propose a resolution that leads to satisfying results while being intuitive and fast to compute. Our hypothesis is to consider the waterfall segment as a closed pipe, with a constant flow and a triangular cross section (Figure 8). In this case, the flow can be expressed by:

γ = δ ε ζ 2 . ( 1
)
Since we have one equation with three unknowns, and complex inter-dependencies, a physically accurate solution submitted to COMPUTER GRAPHICS Forum (12/2014). should rely on strong assumptions, which cannot be justified in our context. Instead, we decided to solve the system by providing simple and intuitive functions to compute these variables.

We first propose to solve for the speed as a simple linear function of slope s: δ = ks, k ∈ R + . Then, we set the depth as a function f d of the width and of the segment type as: ζ = f d (ε, γ). Note that our slope constraint enforces that s > 0 everywhere for a waterfall segment. We used k = 10 in all our scenes. In our prototype, we use only one profile function, f d (ε, γ) = 1 2 γ. Applying Equation (1), it is now possible to compute the width and to deduce the depth value.

While simple, our model efficiently provides plausible results using intuitive parameters. In Figure 9, the procedural component of our modeling system correctly handles a reduction of the input flow: when the upstream flow is manually reduced, the following free-fall flow reduces accordingly. The subsequent nodes in the graph are then affected. Note how the type of the outgoing element automatically changes, i.e., from ledge to plunge. Note that nothing prevents our method from being extended to account for more complex river profiles, e.

Terrains Adapting to Waterfalls

In this section, we discuss how to adapt the terrain to the waterfalls. First, we generate the waterfall integration mesh (Section 6.1), then we use it to generate vectorial constraints to deform the terrain (Section 6.2) and to generate additional maps (Section 6.3). The latter are used for the visual integration of the waterfall (e.g., texture changes on the terrain and on the waterfall) and for the generation of procedural decorations.

Integration Mesh Generation

The integration mesh is a 3D mesh defining the waterfall surface. It is composed of the surface meshes for all waterfall elements (Figure 10), which are generated using the waterfall network trajectories and the widths of the elements.

The meshes for the contact elements are computed by extruding the width ε along the segment trajectories. We carefully handle the meshing connections at branch intersections, so no meshes overlap one another. Free-fall meshes are extruded from the borders of their incoming segments (e.g., pool borders) and follow the free-fall element curve. Pool meshes are simply triangulated from their contours. contact pool free-fall branch

Constraint-based Terrain Deformation

The terrain is stored as a heightfield, on which a deformation map (displacement map) is applied to adapt the terrain to waterfalls. This deformation map is computed using a Poisson solution applied to vectorial constraints, as inspired by the terrain modeling method of Hnaidi et al. [HGA * 10]. However, our solution differs from theirs in that it propagates deformations (i.e., height differences) instead of heights, and it generates three different types of constraints (border, riverbed, and overhang) using our procedural model. As a result of our approach, small details on the original terrain will remain on the terrain after deformation.

Border constraints. These constraints force the water to naturally follow waterfall trajectories, by defining height and gradient constraints on the mesh contour (Figure 11(top left)).

Let Ω be the contour of the integration mesh (Figure 12). For each point p = (x, y, z) ∈ Ω, we define d as the difference between the terrain height h and the height at the point, i.e., d = zh(x, y) (Figure 11(top center)). The difference d is the value that will be diffused within our solver of constraints. We also diffuse a gradient constraint to ensure that the neighborhood of the terrain just outside of the waterfall is higher than its border. Please refer to Hnaidi et al. [HGA * 10] for more details about the diffusion algorithm. For a given point x located inside the waterfall border (see Figure 12), we create a height constraint based on its distance to y ∈ Ω, the nearest point of the border, to which a profile function f d is applied. The operation is repeated on a dense sampling basis. This defines a set of additional elevation constraints processed using our solver (see Figure 11(bottom)).

Overhang constraints. In order to create overhangs, we generate a horizontal displacement map based on the same diffusion technique. This map is then used to deform the terrain, as presented by Gamito and Musgrave [START_REF] Gamito | Procedural landscapes with overhangs[END_REF].

Along the border at the top of a free-fall, we set a displacement constraint λu, where u is the free-fall direction and λ a constant defined by the user. In addition, we set a displacement constraint -λu along the border of the receiving pool, under the free-fall.

As shown in Figure 13, these two constraints generate a flipped "S" curve, with the top of the overhang extending out of the terrain in the direction of the flow, and the bottom of the receiving pool extending into the cliff, due to erosion and falling rocks.

Procedural Decoration

To improve the integration of the waterfall into the terrain, we use several procedural decoration maps, generated using the footprint of the integration mesh and the waterfall elements type (Figure 16). A decoration map is computed by rendering the integration mesh viewed from above into a texture, similar to the road footprints of Bruneton and Neyret [START_REF] Bruneton | Real-time rendering and editing of vector-based terrains[END_REF]. During this process, we render each submesh in a grey-scale map, depending on its type and on the map that is being computed.

Terrain decorations. Terrain decorations are generated using the water map, which corresponds to the integration mesh footprint. This map is used to mask the procedural seeding of trees and plants to prevent a generation within the water surface, and also to change the terrain texture to, for instance, a bedrock one.

Water decorations. Table 1 lists a set of parameter values depending on the waterfall type. By using these values as greyscale values during the map computations, we generate a foam map, a rock map, and a disturbance map. The foam map identifies the presence of foam on the water surface and is used to select the water diffuse texture; the rock map indicates the density of rocks to generate; and the disturbance map, the amplitude of the waves on the water surface (Figure 14). The variation of values depending on the waterfall type allows increases of the visual difference between them, and improves the appearance of the scene. Note that some filtering is applied to these maps to reduce visible transitions between different types.

Speed map. The speed map is a texture that represents the 2D speed of all water in contact with the terrain in the scene. It is used for animating the textures of the water surface. Figure 15 shows how the speed of the water is computed depending on the type of the segment. We use three different approaches to compute the surface speed of contacts, pools, and branches. For a contact, the speed at point x is given by

δx = δc 1 - x-c b-c
where c and b are the projections of x along the cross section on the main axis and on the shore respectively. For a pool, a fixed number of 2D fluid simulation steps [START_REF] Stam | Stable fluids[END_REF] are evaluated. For branches, we use 2D interpolation based on a standard technique of weighting by the inverse distance, where each point p i is considered as a velocity constraint δ i . The speed within a branch is given by δx = ∑ i ω i δ i . Interpolation weights ω i are computed using the equation:

ω i = ωi ∑ i ωi ωi = ∏ j 1 - x -p i p j -p i .

Implementation and Results

The system is implemented in C++, using OpenGL and GLSL Compute Shaders. The computations are performed on an NVidia 660GTX GPU and an Intel R Xeon R E5-1650 CPU, running at 3.20 GHz with 16 GB of memory. The system uses two threads: one CPU thread for the interface and computation control, and one CPU/GPU thread for the GPU computations and rendering.

Rendering. The waterfalls in our editor are rendered in real time using the integration mesh and the parameter maps computed earlier (Figure 14). We use the technique of "tiled Evaluation. Figures 17 and18 show an overview of our system under editing, with different stages described in the caption of the figure. An accompanying video gives a much better understanding of our system in action, and illustrates several of the features described in the previous sections.

In Figure 1, we show a photo of a real waterfall network, and the result of a 10-minute session with our modeling system; we started with a terrain resembling the original real submitted to COMPUTER GRAPHICS Forum (12/2014). terrain, but without riverbeds. The figure shows that coherent waterfalls similar to those on the photo can be easily modeled, while guaranteeing their physical plausibility.

We organized a user modeling session with two experienced digital artists. After a 20-minute training period, they were both able to create waterfall scenes such as the waterfall presented in Figure 19, all in under 30 minutes.

We asked to reproduce the scene of Figure 1 with classical modeling tools, such as Autodesk R Maya R , it took them more than two days to reach an equivalent level of detail for both the waterfall and terrain deformation; the longest part being the manual deformation of the heightfield to match the waterfall mesh. Of course, Maya R has not been designed to specifically model waterfalls, but this preliminary experiment shows how specialized tools can be beneficial. The artists were pleased by the ease of use of our system and its efficiency. However, they expressed a desire for finer control over the result. Performance. Our system generates a complex waterfall network over a terrain in a few seconds (see Table 2). The number of elements does not have a huge impact on the computation time. Indeed, most of our algorithm uses a fixedsize grid, so its complexity is independent of the number of waterfall elements. When increasing the number of elements, only the time for mesh generation, the riverbed consubmitted to COMPUTER GRAPHICS Forum (12/2014). straints dense sampling (sub-part of the terrain deformation algorithm), and the internal speed computation of pool vary noticeably.

In order to ensure the consistency of our results, each time an element is modified by the user, we reperform all computations. Many simple optimizations can readily detect what needs to be recomputed, and thus greatly improve the efficiency of our system; however, we felt that any such optimizations were not necessary in the current version of our prototype. Limitations. The first limitation of our method is that the terrain is only adapted locally, and therefore does not preserve any global hydraulic properties. Indeed, a waterfall network can be created at an unplausible location on the terrain, failing to respect natural river paths shaped by the terrain slope. While this may lead to unplausible terrains, it also gives more artistic freedom to the user, which we feel is an important property of our system.

#

The heightfield representation of the terrain is another limitation, as it prevents the creation of caves and underground waterfalls, although horizontal displacement maps [START_REF] Gamito | Procedural landscapes with overhangs[END_REF] enable us to create overhangs. With support for stack-based terrains [PGGM09], our system could handle more complex terrain elements.

In addition, our adaptation method relies on a grid-based algorithm, which limits its application to relatively small terrains. In our examples, we used a 2048 × 2048 grid with a resolution of 10 pixels per meter. We could adapt our method to multi-scale editing and rendering [START_REF] Neyret | Scalable real-time animation of rivers[END_REF], by dividing the terrain into several tiles, each calculated independently [START_REF] Ťava | Interactive terrain modeling using hydraulic erosion[END_REF], but we have not yet done so.

Moreover, the recursive nature of our algorithm limits the adaptation of paths in complex cases. Consequently, a good position selection at a given step does not imply a good position for the final shape. For instance, if there are too many obstacles on the path, the heuristic will select a point to avoid them globally, but this selection may prevent further steps to avoid them. A solution inspired by algorithms for procedural roads [GPMG10], should be applicable to procedural river trajectories.

Finally, even if our algorithm supports interactive flow variations that change waterfall geometry and its adaptation to the terrain, we do not support a flow variation without changing the terrain adaptation. For example a drying-out waterfall or a river flood cannot be modeled at this moment with our system.

Conclusion

We presented the first interactive procedural modeling system for the design of waterfalls. Our system relies on a tight coupling of automatic generation and user interaction, where complex constraints and tedious tasks are handled by the procedural component of the system while enabling global and local user control. This leads to an improved user experience and the possibility for more creative modeling.

A number of aspects of our method are dedicated to the special case of waterfalls, including their categorization based on the slope-flow graph, and adapted tools for their interactive modeling. In fact, we are the first to present a system for the interactive design of coherent waterfalls. However, the methodology for the design of our system, with notably the way we interleave high-level user control with automatic processes to check consistency and add details, could be generalized to the modeling of many other natural sceneries, as well as to even less natural complex models.

Future work could focus on the preservation of the hydraulic properties of the terrain during its deformation, or tackle the development of volumetric algorithms allowing the creation of truly 3D waterfall networks.

Figure 2 :

 2 Figure 2: Artistic drawings illustrating the different types of running water and free-falls that can be found in nature.

Figure 3 :

 3 Figure 3: Our classification of waterfall types.

Figure 4 :

 4 Figure4: Processing pipeline used for creating waterfall scenes. The first step is the creation by the artist of the controller network U . Then, a hydraulic graph G is generated; the width of an arc encodes flow quantity. The waterfall network W is then generated with a subdivision algorithm, and the waterfall types are determined. Finally, the integration mesh M is generated and used to deform the terrain and generate the procedural details.

Figure 6 :

 6 Figure 6: Flow repartition in a branch. Left: Input and output directions. Right: Resulting exchange, drawn as segment thicknesses.

Figure 7 :

 7 Figure 7: One step of our recursive subdivision process for waterfall network segments: the segment formed by p i and p i+1 is subdivided at m, which is moved along perpendicular direction u. The final position x is the point that minimizes our cost function C(x).

  g., as introduced by Génevaux et al. [GGG * 13].

Figure 9 :

 9 Figure 9: Varying the flow within a graph. The waterfall elements are changed automatically, in conformity with the classification, and the visual aspects are immediately adapted to the changes.

Figure 10 :

 10 Figure 10: Integration mesh composed by the individual waterfall element meshes.

  d

Figure 11 :

 11 Figure 11: Border and riverbed constraints. Top left: Border constraints. Center: Riverbed constraints. Right: Final terrain. Bottom left: Original terrain. Right: Deformed terrain with apparent riverbed.

Figure 12 :

 12 Figure 12: Integration mesh borders.

Figure 13 :

 13 Figure 13: Top: Horizontal constraints modeling overhangs. Bottom: Free-fall without and with an overhang.

  submitted to COMPUTER GRAPHICS Forum (12/2014).

Figure 14 :

 14 Figure14: Using the integration mesh and the waterfall types, we generate various maps used to render the waterfalls.

Figure 15 :

 15 Figure 15: Internal speed computation for a contact (left), a pool (center), and a branch (right).

Figure 16 :

 16 Figure 16: Incremental representation of procedural decorations. In usual order: Terrain only, adding rocks, water, foam, speed map, final result with vegetation.

Figure 17 :

 17 Figure 17: Overview of the system. Top: Original terrain (left), creating the waterfall network (right). Middle: Generating the control mesh and deducing the types (left), resolving the speeds (right). Bottom: Adapted terrain (left), the final scene (right).

Figure 18 :

 18 Figure 18: Another example of waterfall scene. Top left: Controller network. Top right: Integration mesh with types. Bottom left: Deformed terrain. Bottom right: Final scene.

Figure 19 :

 19 Figure 19: Waterfalls modeled by one of our digital artists.

Table 1 :

 1 Classification of running-water segments that may appear in a waterfall network. The coordinates (slope, flow)

	Type	Slope Flow Foam Rocks Dist. Part.
	Ribbon	π/2	1	0	0.3	0	0
	Plunge	π/2	2	1	0.4	0	0.2
	Ledge	π/2	5	0	0.5	0	0.5
	Cataract π/2	8	1	1	1	1
	Stream π/16	2	0	0.2	0	0
	River	π/32	5	0	0.2	0.1	0
	Rapid	π/16	5	0.2	0.5	0.5	0
	Cascade π/8	3	0.5	1	1	0
	Horsetail π/4	2	1	0	0	0
	Block	π/4	6	1	0	0.2	0

Table 2 :

 2 Computation time in ms Fig. n G W M Terrain Speed Maps Proc.From left to right, the columns of the table list the figure number and the number of waterfall controllers, followed by computation times for the hydraulic graph generation (G), waterfall network generation (W), mesh generation (M), terrain adaptation using a 2048 × 2048 resolution, speed map generation, details map (foam, disturbance, and rocks) generation, and procedural detail generation.

	1 36 2 2 112	758	258	428	284
	17 17 1 1 177	677	384	404	297
	18 29 1 1 77	871	264	494	324
	19 26 1 2 30	1115	284	482	302

submitted to COMPUTER GRAPHICS Forum (12/2014).

Acknowledgements

We thank the anonymous reviewers for their feeback, helping substantially to improve the paper. Derek Nowrouzezahrai helped with the final writing. We also thank the infographists, Laura Paiardini and Estelle Charleroy, for their help to design our prototype application. This work was partially funded by the advanced grant EXPRES-SIVE from the European Research Council (ERC-2011-ADG_20110209), an Exploradoc grant from Région Rhône-Alpes, and the GRAND NCE in Canada.


H Farhat
email: hannafarhat@hotmail.com

H Kdouh

C Brousseau

G Zaharia

G Grunfelder

G El Zein

Radio Wave Propagation Characterization Between Adjacent Decks on Board Ships

Keywords: wave propagation, channel sounding, ships

The objective of this paper is to characterize the electromagnetic waves propagation between adjacent decks on board a ship. First, narrowband measurements were conducted at 2.45 GHz to study the radio coverage. Then, double directional channel measurements were performed at 3.5 GHz to determine the propagation mechanisms in this hostile environment. The aim of these measurement campaigns is to find a solution to maintain the connectivity of a wireless sensor network between different decks on a ship in spite of the metallic structure of decks.

I. INTRODUCTION

A possible solution to reduce the kilometers of cables used to connect thousands of sensors to control units on board ships is the use of Wireless Sensor Networks (WSNs). Moreover, ferry companies aim to enhance the passenger satisfaction by equipping their boats with WiFi networks. However, the metallic structure of ships severely limits the propagation of Electromagnetic (EM) waves between two points within a ship, especially between adjacent rooms, and between adjacent decks. Therefore, it is necessary to study the EM waves propagation in this particular environment in order to optimize the radio devices placement. A few works have verified the feasibility of wireless communications on board ships [START_REF] Lynch | A summary review of wireless sensors and sensor networks for structural health monitoring[END_REF]. However, spatio-temporal characterization of EM waves propagation between adjacent decks has not been investigated yet. In this paper, we will study the critical communication scenario, in which the transmitter (Tx) and the receiver (Rx) are located in two adjacent decks. First, narrowband measurements were conducted at 2.45 GHz to study the radio coverage. Then, wideband double directional channel [START_REF] Steinbauer | The double-directional radio channel[END_REF] measurements were performed at 3.5 GHz to determine the main propagation mechanisms in this particular configuration. The results can be used to determine the wireless devices best location in order to maintain the connectivity of a shipboard wireless network covering different decks. This paper is organized as follows. In section II, we present the narrowband measurement results. Then, in section III, the double directional channel measurement results are described. Finally, section IV concludes this paper and draws some perspectives of this work.

II. NARROWBAND CHANNEL MEASUREMENTS

Due to the low data rate of a shipboard WSN, Continuous Wave (CW) measurements are sufficient to characterize the radio coverage of a WSN because the bandwidth of the transmitted signal is much less than the coherence bandwidth of the propagation channel. CW measurements were conducted on the 8th and 9th passenger decks (D8 and D9) of the "Armorique" ferry. The Tx is constituted of a signal generator which delivers 0 dBm CW signal at 2.45 GHz. The Rx is composed of a spectrum analyzer operating in zero-span mode and an antenna positioner. Two omnidirectional conical monopole antennas (vertically polarized) were used at the Tx and Rx. The Tx is located in a stairway of D8 and the Rx is successively located at four different locations on D9. To minimize the multipath fading effect, the received power values were averaged along 20 wavelength circular track using 250 power samples. Fig. 1 shows the layout of D8 and D9, the Tx location at D8 and Rx locations (Rx1 to Rx4) on D9 represented by stars. As shown in Fig. 1, D8 and D9 are connected by two stairways. "SW1" is located on the top right corner of the decks and contains Tx in its lower level (D8) and Rx1 at its higher level (D9). "SW2" is located in the middle. Each stairway has two metallic entrance doors at both D8 and D9. "Door 1" and "Door 2" are the two entrance doors of SW1 at D8 and D9 respectively, and "Door 3" is the entrance door of SW2 at D9. Firstly, for Rx1, the Tx and Rx are in the same vertical plan, but the Line-of-Sight (LoS) path is blocked by the metallic floor of D9. However, the wireless communication has a good quality and the received power is around -48 dBm. Communication is possible thanks to the multipath propagation. Secondly, for Rx2 (located at 4.8 m from Door 2 and 8.8 m from Door 3), several cases have been studied with respect to the status of Door 1 and Door 2 (closed or open). When these two doors are open, the received power is around -54 dBm. When closing only Door 2 (at D9), the received power decreases to -74 dBm, which corresponds to a supplementary path loss of 20 dB with respect to the previous case. However, when closing both Door 1 and Door 2, the received power becomes -68 dBm. Hence, the closure of Door 1 at D8 increases the received power at Rx2. This result proves that the main part of energy achieving Rx2 comes from SW1. Closing Door 1 reduced the radio leakage in D8 and kept the energy in SW1, in which it is then reflected to achieve Rx2.

Thirdly, for Rx3 (located in front of Door 3 at 14 m from Door 2), when Door 1 is open and Door 2 is closed, the received power is -77 dBm. It decreases to -88 dBm when closing Door 3. This result shows that, in contrast with location Rx2, a part of energy achieving Rx3 propagates through SW2. Measurement results at Rx2 and Rx3 show that the received power increases when the Rx is located near stairways which may be sources of radio leakage between the two adjacent decks. This view explains the weakness of the received signal at Rx4 (-88 dBm) located far from stairways (at 32 m from Door 2 and 18 m from Door 3). Therefore, it is necessary to place nodes in stairways to maintain the connectivity of a shipboard WSN on different decks. These nodes may cover an area with a radius up to 15 m on the upper deck even when closing all stairways entrance doors. Opening the doors may increase the radius of radio coverage to 32 m. With this measurement configuration we can characterize 360° azimuthal double directional channel at both link sides. The Tx is located at the 7th deck (D7) which is an open hall containing restaurant and bars, with big glazed walls. This reduces the number of EM reflections and facilitates the separation of the main DoDs. The Rx is located on D8 within a stairway. The Tx is located at 25 m from the entrance door of the stairway that contains the Rx. D7 and D8 are connected by three stairways. The entrance doors of the first and third stairways were open during measurements, whereas the door of the second one was closed. We assume a quasi time-invariant channel during the measurements. Attention was paid that no people were moving in the surrounding area. Several measurements were taken and averaged to reduce the noise effect. The objective is to characterize the double directional channel impulse response, from which we gave a special attention to the DoDs and DoAs which are the most important parameters to determine EM leakage between adjacent decks. Collected data are analyzed using the high resolution algorithm SAGE [START_REF] Fleury | Channel parameter estimation in mobile radio environments using the SAGE algorithm[END_REF]. Fig. 2 shows the power normalized (with respect to the strongest path) DoDs at the Tx side (represented by stars). Some principal propagation directions can be easily identified for the DoDs (at 0°, 180°, 210° and between 120° and 150°). The most important paths are clearly directed to the two first stairways (at 180° and between 120° and 150°). This result shows that EM waves propagate mainly through "Stairway 1" and "Stairway 2" that connect D7 and D8. A few number of paths pass through "Stairway 3". This measurement result confirms the conclusion drawn from the narrowband measurements, which considers that stairways are the main sources of radio leakage between adjacent decks.

IV. CONCLUSION

This paper reports on narrowband and double directional channel measurements conducted on a modern ship to characterize the EM waves propagation between adjacent decks. In spite of the totally metallic structure of decks, the obtained results showed that inter-decks radio communication remains possible. EM waves propagate mainly through stairways. This observation brings interesting perspectives for wireless networks deployment in such a complex environment. It would be necessary to place nodes in the stairways to maintain the connectivity of a shipboard wireless network.

Fig. 1 .

 1 Fig. 1. Layout of Deck 8 and Deck 9.

Fig. 2 .

 2 Fig. 2. SAGE estimated DODs.


J L Vignet

GANIL R Revenko

A BUNCH EXTENSION MONITOR FOR THE SPIRAL2 LINAC

A bunch extension monitor for the SPIRAL2 LINAC

INTRODUCTION

The SPIRAL2 project [1] is based on a multi-beam LINAC driver in order to allow both ISOL and lowenergy in-flight techniques to produce RIB. A superconducting light/heavy-ion LINAC capable of accelerating 5 mA deuterons up to 40 MeV and 1 mA heavy ions up to 14.5 MeV/u is used to bombard both thick and thin targets. These beams could be used for the production of intense RIB by several reaction mechanisms (fusion, fission, transfer, etc.) and technical methods (ISOL, IGISOL, recoil spectrometers, etc.). The production of high intensity RIB of neutron-rich nuclei will be based on fission of uranium target induced by neutrons, obtained from a deuteron beam impinging on a graphite converter (up to 10 14 fissions/s) or by a direct irradiation with a deuteron, 3 He or 4 He beam.

The accelerating RF of the LINAC [START_REF] Uriot | Configuration de base du LINAC supraconducteur de SPIRAL2 V2.0[END_REF] is 88.0525 MHz. It means that time distance between two bunches is 11.26 ns. The extension of the phase for bunch (±2σ) is 60˚ or ~1,6ns for bunch length. The LINAC can operate at continuous or pulsed mode with period of macropulse between 100 µs and 1s.

Correct adjustment of the LINAC is necessary for obtaining maximal intensity and luminosity on the target. Adjusting includes synchronization of phase for each acceleration section. For this reason, information about bunch length distribution is needed and will be obtained using a Bunch Extension Monitor (BEM). These diagnostic detectors will be placed inside the warm sections between superconducting acceleration cavities. Each warm section consists of two quadrupoles and diagnostic box placed in between. Five BEM will be mounted into the first five diagnostic boxes at the LINAC entrance where deuteron beam energy will be between 1.46 and 2.29 MeV.

BEM DESCRIPTION

Bunch extension monitor is a non-destructive beam diagnostic detector for estimation of the length of LINAC bunches. The principle of BEM measurement is based on registration of x-rays emitted from a thin tungsten wire due to the interaction with ions beam. All components of BEM should meet UHV requirements since the LINAC will be operated at 10 -8 mbar. Moreover, all BEM materials must satisfy required purity conditions for preventing cavity pollution.

BEM Working Principle

The photons emitted from the tungsten wire will be produced due to ionization of atoms hit by beam ions. The ions can knock-out electrons from inner shells and produce electron vacancy. Electrons from outer shells fill this vacancy and emit the difference of energy of bound states as characteristic x-ray photons. The energies of emitted characteristic photon are unique for each element and in case of tungsten they are 60 keV and 11 keV for Kand L-shell ionization respectively. A two stage microchannel plates is used to produce photoelectrons and multiplied them. Output signal is transmitted through coaxial cable to the input of a constant fraction discriminator CFD 7174. Signal of accelerator RF comes to another channel of CFD. The two logical signals from CFD output are sent to a time-to-amplitude converter Ortec 566 as "start" and "stop" signals. The output signal from TAC (whose amplitude is proportional to the time difference) is digitized by multichannel analyzer CANBERRA Multiport II. The Schema of BEM electronics operation is presented in Fig. 1.

BEM Design

The BEM consists of two parts which are installed inside the diagnostic box that can be seen on Fig. 2. The first part of this detector is the system for wire insertion which allows inserting the wire into the beam with precision better than 100 µm. Three different tungsten wires of 150 microns diameter are fixed on a stainless steel frame. Each holder has dimensions 50mm × 50mm where tungsten wire is fixed along the diagonal. Wires of the next BEM are fixed at a perpendicular direction compared to the previous BEM wires to minimize beam distortion. Three wires are installed to have some spare in case of wire damaging without dismounting. The rod is moved by a brushless motor through a system of screw-nut mounted on air side of flange CF100.

The system for wire insertion allows also to perform measurements of current on the wire. Frame with tungsten wires is isolated and electrically connected through vacuum feedthrough on the flange. Current from the wire is measured by a card developed at GANIL and also used for measurements of current from Faraday Cups and slits.

The second part of BEM is the x-ray detector used for registration of photons emitted from wire. Two MCPs Hamamatsu F1551 ( 18mm) are used assembled in chevron configuration and provide gain up to 10 7 at applied voltage -2kV. For obtaining short signal front, the output electrons are collected by a fast readout coaxial anode of 50 Ohms impedance connected on the CF100 flange with N-type connector. Entrance of first MCP is covered by an 8 mm thick copper plate with a 4 mm hole. This collimator permits to select x-rays coming from the center of the wire.

To avoid background signals from ions coming from residual gas ionization, a deflecting grid polarized at +30V is positioned on the top of the collimator.

TEST OF PROTOTYPE

Test GANIL on IRRsud Ligne

A test was performed at GANIL with a beam of 36 Ar +10 ions of 0,98MeV/n energy. A 60s acquisition time with 2.4 µA intensity has given the spectra shown Fig. 3 with a good signal over noise ratio. Bunch length measurements were done at various intensity up to 2.5 µA using a continuous wave beam with adapted voltage on MCP (see Fig 4). FWHM values are stable and the difference of the three curves is due to a different background subtraction for the FWHM calculation.

Comparison with Theoretical Estimation

Results of the beam test were compared with analytical calculation of x-ray production from the wire (see Fig 5). To perform these calculations the cross section of tungsten wire was divided into small elemental cells. For each cell the value of ionization cross section was calculated. The cross-section can be implemented as a function of scaled velocity of projectile ion to velocity of electron on the shell V proj /V elect and can be represented through Gryzinsky function (see Eq.1) [START_REF] Gryzinski | Classical Theory of Atomic Collisions. Theory of Inelastic Collisions[END_REF].

) ( ) ( 2 0 2 V G U Z N E ionis       (1).
Where N is a number of electrons on the shell (K-shell=2, L-shell=8), Z is charge of projectile ion, U is binding energy of electrons at the given shell, G(V) is a function of the scaled velocity V=V proj /V elect and To estimate the scaled velocity the value of beam energy loss is calculated in cell. The attenuation of x-rays in the wire is also taken into account.

The total value of emitted photon from the wire can be obtained by summing the x-ray production in each cell [START_REF] Orlic | TTPIXAN -A Package of computer programs for quantitative thick target PIXE analysis[END_REF]. The final expression for total number of photons emitted in solid angle ratio Ω/4 and registered by MCP can be written as:

             ij L i ionis ij l k atoms ij e x E I N Q      ) ( 4 , (2).
where Ω/4 is a solid angle ratio of detector, ε is the photon registration efficiency by MCP, N atoms is the concentration of target atoms per cm 2 ,  k,l are fluorescence yields of photons for K-, L-shell ionizations, I ij is number of incident ions of the beam at ij-cell per second, σ ionis is ionization cross-section for K-,L-shells as function of incident ion energy, exp(-µ•L ij ) is x-ray attenuation.

The calculations (see Eq. 2) were done for both cases of ionization for K-and L-shell. The relative yield of photons due to ionization of L-shell is almost 10 5 higher than for K-shell. For this reason, only the L shell calculation is shown in comparison with experimental results (see Fig. 5).

An experimental count rate comparison to the theoretical calculation has been made for different beam intensities and MCP applied voltage (see Fig. 5). Variation of the experimental count rate measurement is similar to the calculation curve but values are lower of these due to different parameters not take in count in the calculation like MCP opacity and gain.

Test of BEM Close to Cryomodule

To verify if x-rays emitted from the cavity could disturbed the BEM acquisition, a BEM has been mounted in a LINAC diagnostic box and connected to a SPIRAL2 cryomodule type B (with two superconducting cavities) (see Fig. 6). This test was performed at ORSAY IPN. The principle of operation for bunch extension monitor (BEM) is based on the registration of x-ray emitted from thin tungsten wire.. A superconducting cavity inside of cryomodule can produced x-rays due to electrons bremsstrahlung inside the cavity at high applied accelerating field. The emitted x-rays from cavity can be directly registered by BEM or can create x-rays on the wire and produce unwanted background in both cases. Two series of measurements for each cavity were carried out and average count rates were measured. Each series of measurements contains several points at different values of electric field gradient starting from zero and exceeding of nominal value (operating acceleration field value is 6.5 MV/m).

Also for each value of electric field gradient two measurements were done with inserted wire and without. The average count rate measured without insertion of the wire permits to estimate global x-ray background generated by cavity. The ratio of count rates obtained with inserted wire and without gives the number x-rays number emitted from the wire due to x-ray fluorescence. At the maximum value of 10.6 MV/m, the count rate was about 50 Hz on the worst cavity (Fig. 7). This value is lower by about two to three orders magnitude compared to normal conditions in beam.

SUMMARY

A prototype of non-interceptive beam diagnostic detector for bunch length measurements was developed. Test with proton and ions beams were performed in various conditions. The time resolution for detector and electronics is approximately 10 ps using analogical electronics. The results of test were compared with theoretical estimation of x-ray production. The test of BEM installed close to the cryomodule has shown that BEM measurements will not be perturbated by it during beam acceleration. A test with high intensity proton beam at SARAF facility is scheduled in October 2014 to estimate the maximum count rate possible with this detector. The set of the five BEM will be installed on the LINAC at the beginning of 2015.

Figure 1 :

 1 Figure 1: Principal schema of operation of BEM electronics.

Figure 2 :

 2 Figure 2: General view of BEM installed at diagnostic box of LINAC.

Figure 3 :

 3 Figure 3: Measured bunch profile of 36 Ar +18 with 2.4 µA beam intensity (Value of MCP applied voltage of 1.95kV and time acquisition =60 s).

Figure 4 :

 4 Figure 4: Measured FWHM bunch of 36 Ar +18 versus beam intensity.

Figure 5 :

 5 Figure 5: Comparison of count rates for different values intensities and applied voltage. Theoretical estimation of count rate due L-shell ionisation is presented.

Figure 6 :

 6 Figure 6: Installation of BEM set close to cryomodule B.

Figure 7 :

 7 Figure 7: value of count rate emitted from cryomodule and detected by BEM.

Proceedings of IBIC2014, Monterey, CA, USA ISBN 978-3-95450-141-0

Copyright © 2014 CC-BY-3.0 and by the respective authors

ACKNOWLEDGMENT

This work is funded in frame of CRISP WP3T1 project.