Ghislain Landry

Tsafack Chetsa

Laurent Lefevre
email: laurent.lefevre@ens-lyon.fr

Jean-Marc Pierson
email: pierson@irit.fr

Patricia Stolf
email: stolf@irit.fr

Georges Da Costa
email: dacosta@irit.fr

Beyond CPU Frequency Scaling for a Fine-grained Energy Control of HPC Systems

Keywords: energy optimization, phase identification, hardware performance counters, system adaptation

Modern high performance computing subsystems (HPC) -including processor, network, memory, and IOare provided with power management mechanisms. These include dynamic speed scaling and dynamic resource sleeping. Understanding the behavioral patterns of high performance computing systems at runtime can lead to a multitude of optimization opportunities including controlling and limiting their energy usage.

In this paper, we present a general purpose methodology for optimizing energy performance of HPC systems considering processor, disk and network. We rely on the concept of execution vector along with a partial phase recognition technique for on-the-fly dynamic management without any a priori knowledge of the workload. We demonstrate the effectiveness of our management policy under two real-life workloads. Experimental results show that our management policy in comparison with baseline unmanaged execution saves up to 24% of energy with less than 4% performance overhead for our real-life workloads.

I. INTRODUCTION

In recent years, the number of High Performance Computing (HPC) systems around the world has significantly increased. This can be attributed to the fact that HPC systems are designed to achieve peak performance which in turn enable new levels of innovation and insights for organizations that seek out differentiation with excellence in fields including research and development (R&D), science, etc. For example, in order to bring superior products to market faster, a manufacturer may consider using HPC to build prototypes of its products. While it is tempting to have a HPC system constantly available to tackle big challenges, it is often costly to power and cool such an infrastructure. As a consequence, HPC systems are also ranked based on their power efficiency 1 , i.e., the number of Flops achieved per watt. Beside, the global information and communication technology sector has a considerable share of global carbon dioxide (CO 2 ) emissions which is equivalent to aviation.

Over the past years, the energy consumption problem in HPC systems has widely been investigated. Although nearly all HPC subsystems -including processor, network, memory, and IO -are provided with power saving capabilities, efforts for reducing the power consumption of HPC systems are directed toward the processor to the best of our knowledge. Roughly speaking, existing power saving approaches for high performance computing systems attempt to reduce the processor power consumption ignoring all other subcomponents.

To efficiently address the power consumption issue, HPC systems need to approach things differently than in the past. Processors dominate system power consumption in high performance computing; however, network interconnects, memory and disk altogether consume up to 24% of system power consumption [START_REF] Liu | A survey of the research on power management techniques for high-performance systems[END_REF]. With the processor added, the amount of power consumed by these subsystems exceeds 55% of system power consumption [START_REF] Liu | A survey of the research on power management techniques for high-performance systems[END_REF]. Thus, a fine-grained management of these subsystems can result in significant energy savings.

The majority of HPC workloads alternates during their life cycle between phases of different characteristics. A phase of program behaviour represents an interval of execution during which a measured program or system metric is relatively stable. The most prevalent phases in HPC workloads include memory intensive, compute intensive and communication intensive; they can guide any kind of energy optimization and policy decisions without substantial performance degradation. For example, during a memory intensive phase, the processor's frequency can be shifted from its maximum to a lower value, and the disk sent to sleep. This does not always apply as a typical HPC system is shared by a multitude of applications each having its own characteristics in which case optimizations made for saving energy considering some applications are likely to impact the performance of others. We overcome this limitation by considering phases of the HPC system; instead of those of the workloads.

In this paper, we present an efficient runtime power management architecture that is based on partial phase recognition and system adaptation considering some of the most power hungry HPC subsystems. Partial phase recognition consists of identifying a specific percentage of an ongoing phase with a known phase. Since the implementation of our partial phase recognition is based upon system information, it is independent of any individual application, meaning it does not require any a priori knowledge of applications sharing the infrastructure. Independence from any individual application allows us to offer the same performance to all applications sharing the platform while saving energy.

The rest of the paper is organized as follows: In Section II prior work related to energy reduction in HPC systems is discussed. Partial phase recognition and our global approach can be found in Section III. Section IV presents and analyses experimental results. Finally, Section V concludes our work and presents future directions.

II. RELATED WORK

In the past few years, HPC systems have witnessed the emergence of energy consumption reduction techniques from the hardware level to the software level. At the hardware level, the majority of Information Technology (IT) equipment vendors works either from bottom up, by using the most efficient components in their equipments, and/or by providing their equipments with technologies that can be leveraged to reduce energy consumption of HPC subsystems -such as processor, network, memory, and I/O -during their operation. For example, the majority of modern processors is provided with Dynamic Resource Sleeping (DRS). This technique makes components hibernate to save energy and awakes these on demand. Although major progress has been made, improvements in hardware solutions to reduce energy consumption have been slow, due to the high cost of designing equipments with energy-saving technologies and the increasing demand of raw performance. Our work combines hardware technologies and software solutions to reduce HPC systems' energy consumption.

Unlike hardware approaches, software solutions for reducing HPC systems' energy usage have received extensive attention of researchers over time. In [START_REF] Rountree | Adagio: making DVS practical for complex HPC applications[END_REF], Rountree et al. use node imbalance to reduce the overall energy consumption of a parallel application in a HPC system. They track successive MPI communication calls to divide an application into tasks composed of a communication portion and a computation portion. A slack occurs when a processor is waiting for data during the execution of a task. This leaves the possibility to slow the processor with almost no impact on the overall execution time of the application. Rountree et al. developed Adagio which tracks task execution slacks and computes the appropriate frequency at which it should run. Although the first instance of a task is always run at the highest frequency, further instances of the same task are executed at the frequency that was computed after it is first seen. In [START_REF] Kappiah | Just in time dynamic voltage scaling: exploiting inter-node slack to save energy in mpi programs[END_REF], a tool called Jitter is developed to detect slack moments in performance to perform inter-node imbalance and to use DVFS to adjust the CPU frequency. Our approach differs from that implemented in Adagio in that our fine-grained data collection gives the possibility to differentiate not only computation-intensive and communication-intensive execution portions (which we call phases/regions), but also memory-intensive portions. Memory-intensive phases can be run on a slower core without significant performance penalty [START_REF] Kotla | Characterizing the Impact of Different Memory-Intensity Levels[END_REF]. They also differ in that we consider different green leverages, not only DVFS.

Works done in [START_REF] Choi | Fine-grained dynamic voltage and frequency scaling for precise energy and performance tradeoff based on the ratio of off-chip access to onchip computation times[END_REF], [START_REF] Isci | Live, runtime phase monitoring and prediction on real systems with application to dynamic power management[END_REF] use online techniques to detect applications execution phases, characterize them and set the appropriate CPU frequency accordingly. They rely on hardware monitoring counters to compute runtime statistics such as cache hit/miss ratio, memory access counts, retired instructions counts, etc., which are then used for phase detection and characterization. Policies developed in [START_REF] Choi | Fine-grained dynamic voltage and frequency scaling for precise energy and performance tradeoff based on the ratio of off-chip access to onchip computation times[END_REF], [START_REF] Isci | Live, runtime phase monitoring and prediction on real systems with application to dynamic power management[END_REF] tend to be designed for single task environment. We overcome this limitation by considering each node of the cluster as a black box, which means that we do not focus on any application, but instead on the platform. The flexibility provided by this assumption enables us to track not applications/workloads execution phases, but node's execution phases. Our work also differs from previous works in that we use partial phase recognition instead of phase prediction. Online recognition of communication phases in MPI application was investigated by Lim et al. in [START_REF] Lim | Adaptive, transparent frequency and voltage scaling of communication phases in mpi programs[END_REF]. Once a communication phase is recognized, authors apply CPU DVFS to save energy. They intercept and record the sequence of MPI calls during program execution and consider a segment of program code to be reducible if there are high concentrated MPI calls or if an MPI call is long enough. The CPU is then set to run at the appropriate frequency when the reducible region is recognized again.

While having in common the use of performance counters at runtime, our work differs from those above in two major ways. First our phase detection approach does not rely on a specific HPC subsystem or MPI communication calls. Second, unlike previous research efforts, our model does not focus only on saving the energy consumed by the processor, it also takes advantage of power saving capabilities available on all other HPC subsystems.

In [START_REF] Tsafack | DNA-inspired Scheme for Building the Energy Profile of HPC Systems[END_REF] we developed a cross platform methodology for detecting and characterizing phases in HPC systems. In this previous work we modeled the entire runtime of the system as a state graph associated with a transition matrix. The transition matrix helps determine the next state in the graph. The problem with that approach is that it takes time to get a useful transition matrix since the system may enter different configurations from one execution to another. Hence it may not work well with a system that has varying workloads. Moreover, the entire model is built off-line, while the approach proposed here runs online.

III. METHODOLOGY AND MANAGEMENT POLICIES

High performance computing (HPC) systems are designed to offer peak performance; consequently, any management policy that aims to reduce the energy consumption must not degrade performance. To mitigate performance degradation, it is indispensable to understand the behaviour of the system at hand at runtime. Put simply, optimization proposed is closely related to the behaviour of the system. For instance, scaling the CPU frequency down to its minimum when running CPU-bound workloads may cause significant performance degradation, which is unacceptable. Thus, to efficiently optimize a HPC system at runtime, it is necessary to identify the different behaviours known as phases during execution. In this section, we discuss our phase identification approach along with management policies.

A. System's Phases Tracking and Characterizing

Our methodology relies on the concept of execution vector (EV) which is similar to power vectors (PV) in [START_REF] Isci | Identifying program power phase behavior using power vectors[END_REF]. An execution vector is a column vector whose entries are system's metrics including hardware performance counters, network byte sent/received and disk read/write counts. For convenience, we will refer to these system's metrics as "sensors" in the rest of the paper. Each hardware performance counter related sensor represents the access rate to a specific hardware register over a given time interval, whereas network and disk related sensors monitor network and disk activities respectively. We refer to the literature for selecting hardware performance counters related sensors, these include: number of instructions, last level cache accesses and misses, branch misses and predictions, etc. The sampling rate corresponding to the time interval after which each sensor is read depends on the granularity. While a larger sampling rate may hide information regarding the system's behaviour, a smaller sampling rate may incur a non negligible overhead. We made a trade off of one second considering the literature. As a consequence, each execution vector is timestamped with the time at which it was sampled.

In order to identify system's execution phases, we define the resemblance or similarity between execution vectors as the manhattan distance between them. This similarity is used to cluster EVs along the execution time-line as follow: two consecutive EVs along the execution time-line belong to the same group or are similar if the manhattan distance between them is bellow a similarity threshold (denoted as ST in the following). We define the similarity threshold as a percentage of the maximum known distance between all consecutive execution vectors. For example, given a similarity threshold of 10%, two consecutive EVs belong to the same group if the manhattan distance between them is less than 10% of the maximum existing distance between all consecutive execution vectors. Computing the maximum existing distance between consecutive execution vectors at runtime is quite straightforward; the process is summarized in the following algorithm:

• initially (the beginning of the execution time-line), the maximum existing distance between all consecutive EVs referred to as d max is set to zero. Figure 1. Phase identification using the similarity between consecutive execution vector as phase identification metric. This is a zoomed-in view of the traces collected on one node when the system was running Molecular Dynamics Simulation; "distance" represents the variation of the distance between consecutive execution vectors.

Where d i+1,i is the the manhattan distance between EVs sampled at time t i and t i+1 respectively. By definition, the manhattan distance between two points in an n-dimensional space is the distance between them if a grid-like path is followed. The manhattan distance offers the advantage that it does not depend on the translation of the coordinate system with respect to a coordinate axis, i.e., it weighs more heavily differences in each dimension.

Knowing that the behaviour of the system is relatively stable during a phase (a phase is a time period during which the behaviour of the system is stable) and assuming that this stability is translated into execution vectors sampled during that phase, we define a phase as any behaviour delimited by two successive manhattan distances exceeding the similarity threshold. Therefore, considering Figure 1, where the xaxis represents the execution time-line; with a similarity threshold of 50% and assuming the maximum manhattan distance between two consecutive execution vectors along the execution time-line is 0.2; we can observe seven phases separated by six micro-phases or transition phases (0.1 is 50% of 0.2; referring to the definition of a phase, phases are delimited by distances greater than 0.1). We can also notice that these phases correspond to variations reported in the access rate of plotted performance counters.

The rationale behind phase tracking is the use of characteristics of known phases for optimizing similar phases. An effective phase characterization is therefore needed. To this extend, once a phase is detected, we perform the Principal Component Analyses (PCA) on the dataset composed of EVs pertaining to that phase. After PCA, five sensors among those contributing less to the first principal axes of PCA are kept for phase characterization (PCA can be seen as a variable reduction procedure, it helps finding underlying factors that shape the data). This choice is motivated by the assumption that information regarding what the system did not do during a phase is captured by sensors contributing less to the first principal axes of PCA. These 5 sensors characterise the phase. In addition, we associate each newly detected phase with a reference vector and its duration. The reference vector is the closest vector to the centroid of the group of vectors sampled during the phase. This reference vector is later used for partial phase recognition.

B. Partial Phase Recognition and System Adaptation

A phase cannot be detected unless it is finished, in which case any system adaptation or optimization accordingly is no longer worth. The literature recommends phase prediction. Predicting the next phase permits adapting the system accordingly. Although phase prediction is very effective in some cases, it is not relevant in this context, for we do not have any a priori knowledge of applications sharing the platform. To overcome this limitation, we use partial phase recognition.

Partial phase recognition consists of identifying an ongoing phase (the phase has started and is not yet finished) P i with a known phase P j only considering the already executed part of P i . The already executed part of P i expressed as a percentage of the length (duration) of P j is referred to as the recognition threshold RT . Thus, with a RT % recognition threshold, and assuming that the reference vector of P j is EV Pj and that its length is n j ; an ongoing phase P i is identified with P j if the manhattan distance between EV Pj and each EV pertaining to the already executed part of P i (corresponding in length to RT % of n j ) is within the similarity threshold ST . The following algorithm summarizes partial phase recognition.

• lets consider p i a phase, and v i its reference vector, n i its duration • p j is partially identify with p i if p j lasts RT% n i -8v j EV of p j , the distance between v j and v i is within threshold We use partial phase recognition to guide system adaptation considering three HPC subsystems: processor, disk and network interconnect. For the processor, we consider three computational levels according to the characteristics of the workload:

• high or cpu-bound: the cpu-bound computational level corresponds to the maximum available CPU frequency, and is used for CPU-bound workloads; • medium or memory-bound: corresponds to an average in between the maximum and the minimum available frequencies; it is mainly used for memory-bound workloads. • low: the system is in the low computational level when the processor is set to its lowest frequency. For the disk, we define two states: active and sleep; where active includes both the disk's active and standby modes. Finally, for the network interconnect we define two data transfer speeds:

• the communication-intensive speed: corresponds to the highest available transfer rate of the network card; • low-communication speed: where the speed of the network interconnect is set to the lowest speed. To better understand how our partial phase recognition translates phase characteristics at runtime, we trace the CPU frequency settings for some nodes involved in the computation. Figure 2 -where the x-axis represents the execution time-line and the y-axis the normalized frequencies (0.2 and 0.253 corresponding to 2 GHz and 2.53 GHz respectively)shows Molecular Dynamics Simulation (MDS) [START_REF] Binder | Molecular Dynamics Simulations[END_REF] traces. For the sake of clarity, only the frequency settings of three nodes among twenty five are plotted. Figure 2 indicates that, although running the same application, different settings may apply on different nodes at the same time. This can be attributed to the fact that on a node frequencies changes depend on the behaviour of the node rather than that of the application.

As we mentioned earlier in this paper, principal component analysis (PCA) is applied to vectors belonging to any newly created phase for selecting five sensors which are used as phase characteristics. These characteristics are translated to system adaptation as detailed in Table I, where the left column lists sensors selected from PCA and the right, the associated decisions. For example, when sensors selected for characterizing a phase does not include any IO-related sensor, we assume the system is running a compute intensive workload in which case the processor is set to its maximum, the disk sent to sleep and the network interconnect speed scaled down. Note in passing that changing the disk's state from sleep state to active state does not appear in Table I, this is because the disk automatically enters the active state when it is accessed.

IV. EXPERIMENTATION AND VALIDATION A. Evaluation Platform

Our evaluation support is a twenty five node cluster set up on the french large scale experimental platform, Grid5000 [START_REF] Bolze | Grid 5000: a large scale and highly reconfigurable experimental grid testbed[END_REF]. Each node is an Intel Xeon X3440 with 4 cores and 16 GB of RAM. Available frequency steps for each core are: 2.53 GHz, 2.40 GHz, 2.27 GHz, 2.13 GHz, 2.00 GHz, 1.87 GHz, 1.73 GHz, 1.60 GHz, 1.47 GHz, 1.33 GHz and 1.20 GHz. During our experiments, low computational level always sets the CPU frequency to the lowest available which is 1.20 GHz, whereas high and medium computational levels Figure 2. Frequency traces resulting from our management approach. To improve the readability, some frequencies were multiplied by a constant.

set the CPU frequency to the highest available (2.53 GHz) and 2.00 GHz respectively. These values are selected based on empirical evidence; prior experiments we ran indicated that the gab between selected frequencies is somehow important. For example, running the system at 2.00 GHz rather than 1.87 GHz may have an impact on performance, without necessary improving energy performance. Each node uses a SATA II hard drive which supports active, sleep and standby states. Infiniband-20G is used for interconnecting nodes. Linux kernel 2.6.35 is installed on each node where perf event is used to read the hardware monitoring counters. MPICH is used as MPI library. For the experiments, we use two real life applications: Molecular Dynamics Simulation (MDS) [START_REF] Binder | Molecular Dynamics Simulations[END_REF] and the Advance Research Weather Research and Forecasting (WRF-ARW) model [START_REF]The Weather Research and Forecasting Model[END_REF], [START_REF] Skamarock | A description of the Advanced Research WRF Version 2[END_REF]. WRF-ARW is a fully compressible conservative-form non-hydrostatic atmospheric model. It uses an explicit timesplitting integration technique to efficiently integrate the Euler equation. The classical Molecular Dynamics solves numerical Newton's equations of motion for the interaction of the many particles system. We monitored each node's power usage on a per second basis using a power distribution unit. To perform CPU frequency scaling and network speed scaling, we use Linux kernel frameworks CPUFreq and Ethtool respectively (for further information about these tools refer to the documentation available in Linux kernel source).

B. Results Analyses and Discussion

For evaluation, we compare our management policy (referred to as managed) with two other classical system's configuration: on-demand and performance. The on-demand configuration refers to the configuration in which the default Linux on-demand governor is enabled on each individual node of the system; whereas under the performance configuration, Linux performance governor is enabled on each node of the system.

We also consider two levels of system adaptation:

• system adaptation level one: It corresponds to the situation in which only processor related optimization are made.

• system adaptation level two: It embraces level one, and additionally considers optimizing the interconnect and the disk. All results presented in this section were obtained considering a similarity threshold ST of 5% for phase detection and a 10% recognition threshold RT for partial phase recognition. These values were chosen based on past experiments. Although Intel Xeon X3440 processors offer support for Turbo boost, that feature was disabled during the course of these experiments.

1) System Adaptation Level One: Processor's Only Optimization: Once partial recognition is successful, meaning an ongoing phase is identified with an existing phase, the characteristics of the existing phase are used for setting the processor's frequency accordingly. Diagrams of Figure 3 and Figure 4 present the average energy consumption and execution time of the two considered applications under different system's configurations. These diagrams indicates that our management policy can save up to 19% of the total energy consumption with less than 4% performance degradation.

The diagrams of Figure 4 shows that on-demand and performance governors nearly achieve the same performance. This is because Linux on-demand governor do not lower the CPU frequency unless the system load decreases below a threshold. Traces of CPU load under WRF-AWR for some nodes of our cluster are shown in Figure 5. The plot indicates that the CPU load remains above 85%, in which case the on-demand and performance governors almost have the same behaviour. With processor's only optimization, our management policy differs from that of the Linux's ondemand governor in that we do not use the system's load as system adaptation metric, which enables us to be better.

2) System Adaptation Level Two: Processor, Disk and Network Optimization: disk and network interconnect improves energy performance up to 24% with the same performance degradation as with processor's only optimization. Roughly speaking the disk and network allow to save an additional 5%, which is an interesting result given that the energy consumed by both the network and the disk represents 8% of the energy consumed by a typical supercomputer [START_REF] Liu | A survey of the research on power management techniques for high-performance systems[END_REF].

V. CONCLUSION In this paper, we demonstrate that the problem of energy consumption in high performance computing can benefit from more than CPU frequency scaling. We introduce an on-line or runtime methodology for improving energy performance of HPC systems considering three HPC subsystems: processor, disk and network interconnect. As proof of concept, we implemented and evaluated a prototype of our management policy. Experimental results prove the effectiveness of our system adaptation policy using partial phase recognition under real-life workloads. Comparison with baseline unmanaged execution shows that our management policy can save up to 19% of energy with less than 4% performance loss only considering the processor and up to 24% considering processor, disk and network interconnect.

We expect our approach to be extended to a large number of energy-aware clusters for optimizing their energy consumption, given that it does not require any a priori knowledge of applications/workloads sharing the platform.

Future work includes investigating memory optimization and integrating workload consolidation and task migration to our management policy. ACKNOWLEDGMENT This work is supported by the INRIA large scale initiative Hemera focused on "developing large scale parallel and distributed experiments". Some experiments of this article were performed on the Grid5000 platform, an initiative from the French Ministry of Research through the ACI GRID incentive action, INRIA, CNRS and RENATER and other contributing partners (http://www.grid5000.fr).

Figure 3 .

 3 Figure 3. Phase tracking and partial recognition guided CPU optimization results: the chart shows average energy consumed by each application under different configurations.

Figure 4 .

 4 Figure 4. Phase tracking and partial recognition guided CPU optimization results: the chart shows average execution time of each application under different configurations.

Figure 6 Figure 5 .

 65 Figure 5. Load traces for one of the nodes: the cluster was running WRF-AWR under the on-demand configuration.

Figure 6 .

 6 Figure 6. Phase tracking and partial recognition guided processor, disk and network interconnect optimization results: the chart shows average energy consumed by each application under different configurations.

Table I TRANSLATION

 I OF PHASE CHARACTERISTICS INTO SYSTEM ADAPTATION (IO RELATED SENSORS INCLUDE NETWORK AND DISK ACTIVITIES).

		Sensors selected from PCA Decisions
		for phase characterization
				cache references &	CPU frequency set to its maximum
				cache misses &	spin down the disk
				IO related sensors	network speed scaled down
				no IO related sensor	CPU frequency set to its lowest
					network speed scaled up
				instructions &	CPU frequency set to its minimum
		last level cache misses (llc)	network speed scaled up
				instructions or llc &	CPU frequency set to its average value
				IO related sensors	network speed scaled down
					spin down the disk
				IO related sensors	CPU frequency set to its maximum
					spin down the disk
					network speed scaled down
			0.3
							node 2
			0.25
			0.2
	normalized frequency	0	0.15	500	1000	1500	node 3	2000
		0		500	1000	1500	node 4	2000
					times (s)
		0		500	1000	1500		2000

http://www.green500.org/.


Gabriel Synnaeve
email: gabriel.synnaeve@gmail.com

Pierre Bessière
email: pierre.bessiere@imag.fr

A Dataset for StarCraft AI & an Example of Armies Clustering

This paper advocates the exploration of the full state of recorded real-time strategy (RTS) games, by human or robotic players, to discover how to reason about tactics and strategy. We present a dataset of StarCraft 1 games encompassing the most of the games' state (not only player's orders). We explain one of the possible usages of this dataset by clustering armies on their compositions. This reduction of armies compositions to mixtures of Gaussian allow for strategic reasoning at the level of the components. We evaluated this clustering method by predicting the outcomes of battles based on armies compositions' mixtures components.

Introduction

Real-time strategy (RTS) games AI is not yet at a level high enough to compete with trained/skilled human players. Particularly, adaptation to different strategies (of which army composition) and to tactics (army moves) are strong indicators of human-played games [START_REF] Hagelbäck | A study on human like characteristics in real time strategy games[END_REF]. So, while micro-management (low level units control) has known tremendous improvements in recent years, the broadest high-level strategic reasoning is not yet an exemplary feature neither of commercial games nor of Star-Craft AI competitions' entries. At best, StarCraft bots have an estimation of the available technology of their opponents and use rules encoding players' knowledge to adapt their strategy. We believe that better strategic reasoning is a matter of abstracting and combining the low level states at an expressive higher level of reasoning. Our approach will be to learn unsupervised representations of low-level features.

We worked on StarCraft: Brood War, which is a canonical RTS game, as Chess is to board games. It had been around since 1998, has sold 10 millions licenses and was the best competitive RTS for more than a decade. There are 3 factions (Protoss, Terran and Zerg) that are totally different in terms of units, build trees / tech trees (directed acyclic graphs of the buildings and technologies) and thus gameplay styles. StarCraft and most RTS games provide a tool to record game logs into replays that can be re-simulated by Copyright c 2012, Association for the Advancement of Artificial Intelligence (www.aaai.org). All rights reserved.

1 StarCraft and its expansion StarCraft: Brood War are trademarks of Blizzard Entertainment TM the game engine. That is this trace mechanism that we used to download and simulate games of professional gamers and highly skilled international competitors.

This paper is separated in two parts. The first part explains what is in the dataset of StarCraft games that we put together. The second part showcases army composition reduction to a mixture of Gaussian distributions, and give some evaluation of this clustering.

Related Work

There are several ways to produce strategic abstractions: from using high-level gamers' vocabulary, and the game rules (build/tech trees), to salient low-level (shallow) features. Other ways include combining low-level and higherlevel strategic representation and/or interdependencies between states and sequences.

Case-based reasoning (CBR) approaches often use extensions of build trees as state lattices (and sets of tactics for each state) as for [START_REF] Aha | Learning to win: Case-based plan selection in a real-time strategy game[END_REF][START_REF] Ponsen | Improving Adaptive Game AI with Evolutionary Learning[END_REF] in Wargus. [START_REF] Ontañón | Case-based planning and execution for real-time strategy games[END_REF] base their real-time case-based planning (CBP) system on a plan dependency graph which is learned from human demonstration in Wargus. In [START_REF] Mishra | Situation assessment for plan retrieval in real-time strategy games[END_REF], they use "situation assessment for plan retrieval" from annotated replays, which recognizes distance to behaviors (a goal and a plan), and selected only the low-level features with the higher information gain. [START_REF] Hsieh | Building a player strategy model by analyzing replays of real-time strategy games[END_REF] based their work on [START_REF] Aha | Learning to win: Case-based plan selection in a real-time strategy game[END_REF] and used StarCraft replays to construct states and building sequences. Strategies are choices of building construction order in their model. [START_REF] Schadd | Opponent modeling in real-time strategy games[END_REF] describe opponent modeling through hierarchically structured models of the opponent behavior and they applied their work to the Spring RTS game (Total Annihilation open source clone). [START_REF] Balla | Uct for tactical assault planning in real-time strategy games[END_REF] applied upper confidence bounds on trees (UCT: a Monte-Carlo planning algorithm) to tactical assault planning in Wargus, their tactical abstraction combines units hit points and locations. In (Synnaeve and Bessière 2011b), they predict the build trees of the opponent a few buildings before they are built. Another approach is to use the gamers' vocabulary of strategies (and openings) to abstract even more what strategies represent (a set of states, of sequences and of intentions) as in [START_REF] Weber | A Data Mining Approach to Strategy Prediction[END_REF]Synnaeve and Bessière 2011a). [START_REF] Dereszynski | Learning probabilistic behavior models in real-time strategy games[END_REF] used an hidden Markov model (HMM) whose states are extracted from (unsupervised) maximum likelihood on a Star-Craft dataset. The HMM parameters are learned from unit counts (both buildings and military units) every 30 seconds and "strategies" are the most frequent sequences of the HMM states according to observations. Few models have incorporated army compositions in their strategy abstractions, except sparsely as an aggregate or boolean existence of unit types. Most strategy abstractions are based on build trees (or tech trees), although a given set of buildings can produce different armies. What we will present here is complementary to these strategic abstractions and should help the military situation assessment.

Dataset

We downloaded more than 8000 replays to keep 7649 uncorrupted, 1 vs. 1 replays from professional gamers leagues and international tournaments of StarCraft, from specialized websites 234 . We then ran them using Brood War API5 and dumped: units' positions, regions' positions, pathfinding distance between regions, resources (every 25 frames), all players' orders, vision events (when units are seen) and attacks (types, positions, outcomes). Basically, we recorded every BWAPI event, plus interesting states and attacks. The dataset is freely available6 , the source code and a documentation are also provided7 [START_REF] Forbus | How Qualitative Spatial Reasoning Can Improve Strategy Game AIs[END_REF] have shown the importance of qualitative spatial reasoning, and it would be too space-consuming to dump the ground distance of every position to any other position. For these reasons, we discretized StarCraft maps in two types of regions:

Regions

• Brood War Terrain Analyzer8 produced regions from a pruned Voronoi diagram on walkable terrain [START_REF] Perkins | Terrain analysis in real-time strategy games: An integrated approach to choke point detection and region decomposition[END_REF]. Chokes are the boundaries of such regions.

• As battles often happens at chokes, we also produced choke-dependent regions (CDR), which are created by doing an additional (distance limited) Voronoi tessellation spawned at chokes. This regions set is

CDR = (regions \ chokes) ∪ chokes

Attacks

We trigger an attack tracking heuristic when one unit dies and there are at least two military units around. We then update this attack until it ends, recording every unit which took part in the fight. We log the position, participating units and fallen units for each player, the attack type and of course the attacker and the defender. Algorithm 1 shows how we detect attacks.

We annotated attacks by four types (but researchers can also produce their own annotations given the state available):

• ground attacks, which may use all types of units (and so form the large majority of attacks).

• air raids, air attacks, which can use only flying units.

• invisible (ground) attacks, which can use only a few specific units in each race (Protoss Dark Templars, Terran Ghosts, Zerg Lurkers).

• drop attacks, which need a transport unit (Protoss Shuttle, Terran Dropship, Zerg Overlord with upgrade).

Algorithm 1 Simplified attack tracking heuristic for extraction from games. The heuristics to determine the attack type and the attack radius and position are not described here.

Information in the dataset

Table 1 shows some metrics about the dataset. Note that the numbers of attacks for a given race have to be divided by (approximatively) two in a given non-mirror match-up. By running the recorded games (replay) through Star-Craft, we were able to recreate the full state of the game. Time is always expressed in game frames (24 frames per second). We recorded three types of files:

• general data ( * .rgd files): records the players' names, the map's name, and all information about events like creation (along with morph), destruction, discovery (for one player), change of ownership (special spell/ability), for each units. It also shows attack events (detected by a heuristic, see below) and dumps the current economical situation every 25 frames: mineral, gas, supply (count and total: maxsupply). • order data ( * .rod files): records all the orders which are given to the units (individually) like move, harvest, attack unit, the orders positions and their issue time. • location data ( * .rld files): records positions of mobile units every 100 frames, and their position in regions and choke-dependent regions if they changed since last measurement. It also stores ground distances (pathfindingwise) matrices between regions and choke-dependent regions in the header.

From this data, one can recreate most of the state of the game: the map key characteristics (or load the map separately), the economy of all players, their tech (all researches and upgrades), all the buildings and units, along with their orders and their positions.

Armies composition

We will consider units engaged in these attacks as armies and will seek a compact description of armies compositions.

Armies clustering

The idea behind armies clustering is to give one "composition" label for each army depending on its composing ratio of the different unit types. Giving a "hard" (unique) label for each army does not work well because armies contain different components of unit types combinations. For instance, a Protoss army can have only a "Zealots+Dragoons" component, but it will often just be one of the components (sometimes the backbone) of the army composition, augmented for instance with "High Templars+Archons".

Because a hard clustering is not an optimal solution, we used a Gaussian mixture model (GMM), which assumes that an army is a mixture (i.e. weighted sum) of several (Gaussian) components. We present the model in the Bayesian programming framework [START_REF] Diard | A survey of probabilistic models using the bayesian programming methodology as a unifying framework[END_REF]: we first describe the variables, the decomposition (independence assumptions) and the forms of the distribution. Then, we explain how we identified (learned) the parameters and lay out the question that we will ask this model in the following parts.

Variables • C ∈ c 1 . . . c K ,

our army clusters/components (C).

There are K units clusters and K depends on the race (the mixture components are not the same for Protoss/Terran/Zerg). • U ∈ ([0, 1] . . . [0, 1]) (length N ), our N dimensional unit types (U ) proportions, i.e. U ∈ [0, 1] N . N is dependent on the race and is the total number of unit types. For instance, an army with equal numbers of Zealots and Dragoons (and nothing else) is represented as {U Zealot = 0.5, U Dragoon = 0.5, ∀ut = Zealot|Dragoon U ut = 0.0}, i.e. U = (0.5, 0.5, 0, . . . , 0) if Zealots and Dragoons are the first two components of the U vector. So i U i = 1 whatever the composition of the army.

Decomposition: For the M battles, the armies compositions are independent across battles, and the unit types proportions vector (army composition) is generated by a mixture of Gaussian components and thus U i depends on C i .

P(U 1...M , C 1...M ) = M i=1 P(U i |C i )P(C i ) Forms • P(U i |C i ) mixture of Gaussian distributions: P(U i |C i = c) = N (µ c , σ 2 c ) • P(C i ) = Categorical(K, p C ): P(C i = c k ) = p k K k=1 p k = 1
Identification (learning): We learned the Gaussian mixture models (GMM) parameters with the expectationmaximization (EM) algorithm on 5 to 15 mixtures with spherical, tied, diagonal and full co-variance matrices, using scikit-learn [START_REF] Pedregosa | Scikitlearn: Machine Learning in Python[END_REF]. We kept the best scoring models (by varying the number of mixtures) according to the Bayesian information criterion (BIC) [START_REF] Schwarz | Estimating the Dimension of a Model[END_REF]. Let θ = (µ 1:K , σ 2 1:K ), being respectively the K different N -dimensional means (µ 1:K ) and the variances (σ 2 1:K ) of the normal distributions. Initialize θ randomly, and let

L(θ; U ) = P(U |θ) = M i=1 K k=1 P(U i |θ, C i = c k )P(C i = c k )
Iterate until convergence (of θ):

1. E-step: Q(θ|θ (t) ) = E[log L(θ; u, C)] = E log M i=1 K k=1 P(u i |C i = c k , θ)P(C i = c k ) 2. M-step: θ (t+1) = argmax θ Q(θ|θ (t) )
Question: For the ith battle (one army with units u):

P(C i |U i = u) = P(C i )P(U i = u|C i )

Counter compositions

In a battle, there are two armies (one for each players), we can thus apply this clustering to both the armies. If we have K clusters and N unit types, the opponent has K clusters and N unit types. We introduce EU and EC, respectively with the same semantics as U and C but for the enemy. In a given battle, we observe u and eu, respectively our army composition and the enemy's army composition. We can ask P(C|U = u) and P(EC|EU = eu).

As StarCraft unit types have strengths and weaknesses against other types, we can learn which clusters should beat other clusters (at equivalent investment) as a probability table. We use Laplace's law of succession ("add-one smoothing") by counting and weighting according to battles results (c > ec means "c beats ec", i.e. we won against the enemy):

P(C = c|EC = ec) =
1 + P(c)P(ec)count battles (c > ec) K + P(ec)count battles with (ec)

Results

We used the dataset presented in this paper to learn all the parameters and perform the benchmarks (by setting 100 test matches aside and learning on the remaining of the dataset). First, we analyze the posteriors of clustering only one army and then we evaluated the clustering as a mean to predict outcomes of battles.

Posterior analysis: Figure 1 shows a parallel plot of army compositions. We removed the less frequent unit types to keep only the 8 most important unit types of the PvP matchup, and we display a 8 dimensional representation of the army composition, each vertical axis represents one dimension. Each line (trajectory in this 8 dimensional space) represents an army composition (engaged in a battle) and gives the percentage of each of the unit types. These lines (armies) are colored with their most probable mixture component, which are shown in the rightmost axis. We have 8 clusters (Gaussian mixtures components): this is not related to the 8 unit types used as the number of mixtures was chosen by BIC score. Expert StarCraft players will directly recognize the clusters of typical armies, here are some of them:

• Light blue corresponds to the "Reaver Drop" tactical squads, which aims are to transport (with the flying Shuttle) the slow Reaver (zone damage artillery) inside the opponent's base to cause massive economical damages.

• Red corresponds to the "Nony" typical army that is played in PvP (lots of Dragoons, supported by Reaver and Shuttle).

• Green corresponds to a High Templar and Archon-heavy army: the gas invested in such high tech units makes it that there are less Dragoons, completed by more Zealots (which cost no gas).

• Purple corresponds to Dark Templar ("sneaky", as Dark Templars are invisible) special tactics (and opening).

Figure 2 showcases the dynamics of clusters components: P(EC t |EC t+1 , for Zerg (vs Protoss) for ∆t of 2 minutes. The diagonal components correspond to those which do not change between t and t + 1 (⇔ t + 2minutes), and so it is normal that they are very high. The other components show the shift between clusters. For instance, the first line seventh column (in (0,6)) square shows a brutal transition from the first component (0) to the seventh ( 6). This may be the production of Mutalisks11 from a previously very low-tech army (Zerglings). A soft rock-paper-scissors: We then used the learned P(C|EC) table to estimate the outcome of the battle. For Figure 1: Parallel plot of a small dataset of Protoss (vs Protoss, i.e. in the PvP match-up) army clusters on most important unit types (for the match-up). Each normalized vertical axis represents the percentage of the units of the given unit type in the army composition (we didn't remove outliers, so most top vertices (tip) represent 100%), except for the rightmost (framed) which links to the most probable GMM component. Note that several traces can (and do) go through the same edge. that, we used battles with limited disparities (the maximum strength ratio of one army over the other) of 1.1 to 1.5. Note that the army which has the superior forces numbers has more than a linear advantage over their opponent (because of focus firing12 ), so a disparity of 1.5 is very high. For information, there is an average of 5 battles per game at a 1.3 disparity threshold, and the numbers of battles (used) per game increase with the disparity threshold.

We also made up a baseline heuristic, which uses the sum of the values of the units to decide which side should win. If we note v(unit) the value of a unit, the heuristic computes unit v(unit) for each army and predicts that the winner is the one with the biggest score. For the value of a unit we used: v(unit) = minerals value + 4 3 gas value + 50supply

Of course, we recall that a random predictor would predict the result of the battle correctly 50% of the time. A summary of the main metrics is shown in Table 2, the first line can be read as: for a forces disparity of 1.1, for Protoss vs Protoss (first column),

• considering only military units the heuristic predicts the outcome of the battle correctly 63% of the time.

the probability of a clusters mixture to win against another (P(C|EC)), without taking the forces sizes into account, predicts the outcome correctly 54% of the time.

the probability of a clusters mixture to win against another, taking also the forces sizes into account (P(C|EC) × unit v(unit)), predicts the outcome correctly 61% of the time.

• considering only all units involved in the battle (military units, plus static defenses and workers): same as above.

Results are given for all match-up (columns) and different forces disparities (lines). The last column sums up the means on all match-ups, with the whole army (military units plus static defenses and workers involved), for the three metrics. Also, without explicitly labeling clusters, one can apply thresholding to special units (Observers, Arbiters, Defilers...) to generate more specific clusters: we did not put these results here (they include too much expertize/tuning) but they sometimes drastically increase prediction scores, as one Observer can change the course of a battle.

We can see that predicting battle outcomes (even with a high disparity) with "just probabilities" of P(C|EC) (without taking the forces into account) gives relevant results as they are always above random predictions. Note that this is a very high level (abstract) view of a battle, we do not consider tactical positions, nor players' attention, actions, etc. Also, it is better (in average) to consider the heuristic with the composition of the army ("prob×heuristic") than to consider the heuristic alone, even for high forces disparity. Our heuristic augmented with the clustering seem to be the best indicator for battle situation assessment. These prediction results with "just prob.", or the fact that heuristic with P(C|EC) tops the heuristic alone, are a proof that the assimilation of armies compositions as Gaussian mixtures of cluster works.

Secondly, and perhaps more importantly, we can view the difference between "just prob." results and random guessing (50%) as the military efficiency improvement that we can (at least) expect from having the right army composition. Indeed, for small forces disparities (up to 1.1 for instance), the prediction based only on army composition ("just prob.": 63.2%) is better than the prediction with the baseline heuristic (61.7%). It means that we can expect to win 63.2% of the time (instead of 50%) with an (almost) equal investment if we have the right composition. Also, when we predict 58.5% Table 2: Winner prediction scores (in %) for the three main metrics. For the left columns ("m"), we considered only military units. For the right columns ("ws") we also considered static defense and workers. The "heuristic" metric is a baseline heuristic for battle winner prediction for comparison using army values, while "just prob." only considers P(C|EC) to predict the winner, and "prob×heuristic" balances the heuristic's predictions with C,EC P(C|EC)P(EC).

of the time the accurate result of a battle with disparity up to 1.5 from "just prob.", this success in prediction is independent of the sizes of the armies. What we predicted is that the player with the better army composition won (and not necessarily the one with more or more expensive units).

Conclusion

We delivered a rich StarCraft dataset which enables the study of tactical and strategic elements of RTS gameplay. Our (successful) previous works on this dataset include learning a tactical model of where and how to attack (both for prediction and decision-making), and the analysis of units movements. We provided the source code of the extracting program (using BWAPI), which can be run on other replays. We proposed and validated an encoding of armies composition which enables efficient situation assessment and strategy adaptation. We believe it can benefit all the current StarCraft AI approaches. Moreover, the probabilistic nature of the model make it deal natively with incomplete information about the opponent's army.

Figure 2 :

 2 Figure 2: Dynamics of clusters: P(EC t |EC t+1 ) for Zerg, with ∆t = 2 minutes

Table 1 :

 1 So, there are 7072 Protoss attacks in PvP and there are not 70,089 attacks by Protoss in PvT but about half that. Detailed numbers about our dataset. XvY means race X vs race Y matches and is an abbreviation of the match-up: PvP stands for Protoss versus Protoss.

	match-up	PvP	PvT	PvZ	TvT	TvZ	ZvZ
	number of games	445	2408	2027	461	2107	199
	number of attacks	7072 70089 40121 16446 42175 2162
	mean attacks/game	15.89 29.11 19.79 35.67 20.02 10.86
	mean time (frames) / game	32342 37772 39137 37717 35740 23898
	mean time (minutes) / game	22.46 26.23 27.18 26.19 24.82 16.60
	actions issued (game engine) / game	24584 33209 31344 26998 29869 21868
	mean regions / game	19.59 19.88 19.69 19.83 20.21 19.31
	mean CDR / game	41.58 41.61 41.57 41.44 42.10 40.70
	mean ground distance 9 region ↔ region 2569	2608	2607	2629	2604	2596
	mean ground distance 10 CDR ↔ CDR	2397	2405	2411	2443	2396	2401

http://www.teamliquid.net

http://www.gosugamers.net

http://www.iccup.com

BWAPI http://code.google.com/p/bwapi/

http://emotion.inrialpes.fr/people/ synnaeve/TLGGICCUP_gosu_data.7z

http://snippyhollow.github.com/bwrepdump/

BWTA http://code.google.com/p/bwta/

Mutalisks are flying units which require to unlock several technologies and thus for which player save up for the production while opening their tech tree.

Efficiently micro-managed, an army 1.5 times superior to their opponents can keep much more than one third of the units alive.


Maxime Meilland

Patrick Rives

Andrew I Comport
email: comport@i3s.unice.fr

Dense RGB-D mapping of large scale

come

Dense RGB-D mapping of large scale environments for real-time localisation and autonomous navigation

Maxime Meilland, Patrick Rives and Andrew Ian Comport

Abstract-This paper presents a method and apparatus for building 3D dense visual maps of large scale environments for real-time localisation and autonomous navigation. We propose a spherical ego-centric representation of the environment which is able to reproduce photo-realistic omnidirectional views of captured environments. This representation is composed of a graph of locally accurate augmented spherical panoramas that allows to generate varying viewpoints through novel view synthesis. The spheres are related by a graph of 6 d.o.f. poses which are estimated through multi-view spherical registration. It is shown that this representation can be used to accurately localise a vehicle navigating within the spherical graph, using only a monocular camera for accurate localisation. To perform this task, an efficient direct image registration technique is employed. This approach directly exploits the advantages of the spherical representation by minimising a photometric error between a current image and a reference sphere. Autonomous navigation results are shown in challenging urban environments, containing pedestrians and other vehicles.

I. INTRODUCTION

Acquiring 3D models of large scale environments is currently a key issue for a wide range of applications ranging from interactive personal guidance devices to autonomous navigation of mobile robots. In these applications it is important, not only for human operators but also for autonomous robots, to maintain a world map that holds a rich set of data including photometric, geometric and saliency information. It will be shown in this paper why it is advantageous to define an ego-centric representation of this information that allows fast model acquisition whilst maintaining optimal realism and accuracy.

An a-priori 3D model simplifies the localisation and navigation task since it allows to decouple the structure and motion estimation problems. Current state of the art approaches mostly rely on global 3D CAD models [START_REF] Hammoudi | Generating raw polygons of street facades from a 2d urban map and terrestrial laser range data[END_REF] that are based on tools and representations that have been developed mainly for texture mapped virtual reality environments. Unfortunately, these representations have difficulties in maintaining true photo-realism and therefore they introduce reconstruction errors and photometric inconsistencies. Furthermore, these models are complicated to acquire and often resort to heavy off-line modelling procedures. Whilst efforts are being made to use sensor acquisition systems that automatically acquire these classical virtual 3D models [START_REF] Craciun | Multi-view scans alignment for 3d spherical mosaicing in large-scale unstructured environments[END_REF], it is suggested in this paper that they are not sufficient to precisely represent realworld data. Alternatively, it is proposed to use an ego-centric model [START_REF] Meilland | A spherical robot-centered representation for urban navigation[END_REF] that represents, as close as possible, real sensor measurements.

A well known ego-centric representation model for camera sensors is the spherical panorama. Multiple cameras systems such as in [START_REF] Baker | A spherical eye from multiple cameras[END_REF] allow construction of high resolution spherical views via image stitching algorithms such as reviewed in [START_REF] Szeliski | Image alignment and stitching: a tutorial[END_REF]. However, contrary to virtual reality models, these tools have been developed mainly for qualitative photoconsistency but they rarely require 3D geometric consistency of the scene. This is mainly due to the fact that, in most cases, it is impossible to obtain 3D structure via triangulation of points when there is no or little baseline between images. Another approach is to use a central catadioptric omnidirectional camera [START_REF] Nayar | Catadioptric omnidirectional camera[END_REF] and warp the image plane onto a unit sphere using the model given in [START_REF] Geyer | A unifying theory for central panoramic systems and practical applications[END_REF]. Unfortunately, that kind of sensor has a poor and varying spatial resolution and therefore is not well adapted to a visual memory of the environment. Furthermore, these approaches assume a unique center of projection, however, manufacturing such a system is still a challenging problem [START_REF] Krishnan | Towards A True Spherical Camera[END_REF].

In order to take advantage of both 3D model based approaches and photometric panoramas it is possible to augment the spherical image with a depth image containing a range for each pixel. The multi-camera sensor proposed in [START_REF] Meilland | Dense visual mapping of large scale environments for real-time localisation[END_REF] allows to acquire high resolution spherical images augmented with depth information. An augmented sphere then allows to perform novel view synthesis in a local domain in all directions [START_REF] Avidan | Novel view synthesis in tensor space[END_REF][6] [START_REF] Meilland | A spherical robot-centered representation for urban navigation[END_REF].

These ego-centric models are, however, local and do not provide a global representation of the environment. This problem can be solved by considering multiple augmented spheres connected by a graph of poses that are positioned optimally in the environment. Simple spherical images positioned in the environment are already found in commercial applications such as Google Street View, and more recently in [START_REF] Kopf | Street slide: Browsing street level imagery[END_REF]. The easiest method for positioning spheres would be via a global positioning system (GPS). However, in urban environments this system fails easily due to satellite occlusion. Alternatively, the robot-centred representation introduced in [17] [START_REF] Meilland | Dense visual mapping of large scale environments for real-time localisation[END_REF] positions augmented views globally within a precise topological graph via accurate spherical visual odometry [START_REF] Comport | Real-time quadrifocal visual odometry[END_REF] and does not require any external sensor. The present paper extends this preliminary work, and demonstrates autonomous navigation results in challenging environments using only a monocular camera for localisation.

II. REAL-TIME EGO-CENTRIC TRACKING

The objective of this work is to perform real-time localisation and autonomous navigation using a known environment model (see Fig. 2). The essential part of this paper is therefore divided into two distinct but inter-related aspects:

• Offline learning -This phase consists in acquiring a 3D model of the environment and representing this information in an optimal manner for "on-line" localisation. It has been chosen to develop a learning approach that is also efficient so that, firstly, environments can be acquired rapidly and secondly, so that the approach may be used for online mapping in the near future. Essentially this involves filming, tracking and mapping the 3D environment (≈1Hz depending on the approach). The local ego-centric 3D model and the global graph learning are illustrated in Section III.

• Online localisation and autonomous navigation -The real-time phase involves estimating the 6 d.o.f. pose of a camera at frame-rate (here 45 Hz), onboard the vehicle. This phase must take into account efficient optimisation techniques that require a maximum amount of computation to be performed "off-line" during the learning phase. An accurate vehicle position should be provided to autonomously control the vehicle, as explained in Sections IV and V.

III. SPHERICAL EGO-CENTRED MODEL

An ego-centric 3D model of the environment is defined by a graph G = {S 1 , . . . , S n ; x 1 , . . . , x m } where S i are augmented spheres that are connected by a minimal parametrisation x of each pose as:

T(x) = e [x]∧ = R t 0 1 ∈ SE(3), (1)
where x ab ∈ R 6 is the 6 d.o.f. twist between the sphere a and b (see Fig. 2) defined as:

x = 1 0 (ω, υ)dt ∈ se(3), (2)
which is the integral of a constant velocity twist which produces a pose T. The operator [.] ∧ is defined as follows:

[x] ∧ = [ω] × υ 0 0 , (3)
where [.] × represents the skew symmetric matrix operator.

A. Augmented visual sphere

Each sphere is defined by the set

S = {I S , P S , Z S , W S }, (4)
where:

-I S is the photometric spherical image. This image is obtained from the custom camera system presented in Section III-B by warping multiple images onto the sphere. -P S = {q 1 , . . . , q n } is a set of evenly spaced points on the unit sphere where q ∈ S 2 . These points have been sampled uniformly on the sphere as in [START_REF] Meilland | A spherical robot-centered representation for urban navigation[END_REF]. -Z S are the depths associated with each pixel which have been obtained from dense stereo matching. The 3D point is subsequently defined in the sphere as P = (q, Z). -W S is a saliency image which contains knowledge of good pixels to use for tracking applications. It is obtained by analysing the Jacobian of the warping function so that the pixels are ordered from best to worst in terms of how they condition the pose estimation problem (the interested reader can see [START_REF] Meilland | A spherical robot-centered representation for urban navigation[END_REF] for more details).

B. Spherical acquisition system

To acquire augmented visual spheres, a custom made multi-camera sensor is employed. This sensor, proposed in [START_REF] Meilland | Dense visual mapping of large scale environments for real-time localisation[END_REF] is composed of six wide angle stereo cameras, placed in a ring configuration so that dense stereo matching [START_REF] Hirschmuller | Stereo processing by semiglobal matching and mutual information[END_REF] can be performed between each camera pair. Such a system allows to acquire high resolution visual spheres, with dense depth information for almost each pixels of the spherical image, as it can be seen on Figure 3.

C. Global sphere positioning

To accurately recover the position of the spheres of the graph with respect to one-another, a 6 d.o.f. multi-camera localisation model is used based on accurate dense localisation [START_REF] Comport | Real-time quadrifocal visual odometry[END_REF] [START_REF] Meilland | A spherical robot-centered representation for urban navigation[END_REF]. Considering S * , an augmented sphere defined in Section III-B, the objective is to compute the pose between a reference sphere S * and the next one S. The localisation problem (also known as visual odometry) is solved using a direct 3D image registration technique that provides accurate pose estimation [START_REF] Meilland | Dense visual mapping of large scale environments for real-time localisation[END_REF]. Since this is a local optimisation approach it is assumed that the camera framerate is high (30Hz) and that interframe displacements are small (≤ 2m), meaning a maximum speed of ∼ 200km/h. It is noted here that dense visual odometry is computationally efficient and locally very accurate [START_REF] Comport | Real-time quadrifocal visual odometry[END_REF] so it has been deemed unnecessary to perform costly bundle adjustment on local visibility windows (although this slightly improves the estimate, it makes timely scene acquisition practically infeasible).

In order to avoid reconstructing each sphere from the learning sequence and to obtained a graph with minimal redundancy, a robust statistical criteria is observed to choose where to reconstruct a sphere and add it to the graph. This allows compressing the original sequence of images to only few spherical images (see [START_REF] Meilland | Dense visual mapping of large scale environments for real-time localisation[END_REF]).

Since a visual odometry approach is employed, drift is integrated over large scale trajectories (e.g. ≤ 1%). To overcome this problem, which might lead to inconsistencies in the graph, a spherical loop closure detection is employed [START_REF] Chapoulie | A spherical representation for efficient visual loop closing[END_REF]. This technique detects loop closures from the appearance of the images, using SIFT descriptors [START_REF] Lowe | Distinctive image features from scale-invariant keypoints[END_REF]. Loop closure detections allow to add new edges to the graph, corresponding to new constraints. The final graph can be globally optimised using a graph optimisation approach (e.g. [START_REF] Grisetti | Efficient estimation of accurate maximum likelihood maps in 3d[END_REF]).

IV. REAL-TIME LOCALISATION

It is considered that during online navigation, a current image I, captured by a generic camera (e.g. monocular, stereo or omnidirectional) and an initial guess T of the current camera position within the graph are available. This initial guess permits the extraction of the closest reference sphere S * from the graph. Contrary to non-spherical approaches, a sphere provides all viewing directions and therefore it is not necessary to consider the rotational distance (to ensure image overlap). The closest sphere is subsequently determined uniquely by translational distance. In particular this avoids choosing a reference sphere that has similar rotation but large translational difference which induces self occlusions of buildings and also differences in image resolution caused by distance (which affects direct registration methods).

Since a sphere provides all local information necessary for 6 d.o.f. localisation (geometric and photometric information), an accurate estimation of the pose is obtained by an efficient direct minimization: e(x) = I w( TT(x); s(Z, q * S )) -I * S s(Z S , q * S ) , [START_REF] Chapoulie | A spherical representation for efficient visual loop closing[END_REF] where x is the unknown 6 d.o.f. pose increment. The warping function w(.) transfers the current image intensities onto the reference sphere pixels q S through novel view synthesis [START_REF] Avidan | Novel view synthesis in tensor space[END_REF], using depth information Z S . The function s(.) selects only informative pixels, w.r.t. the saliency map W S which is already pre-computed on the reference sphere [START_REF] Meilland | A spherical robot-centered representation for urban navigation[END_REF]. This selection speeds up the tracking algorithm without neither degrading observability of 3D motion nor accuracy.

The error function e(x) is minimized using an iterative non-linear optimization (IRLS) detailed in Appendix VIII-A. The estimation is updated at each step by an homogeneous transformation:

T ← TT(x), (6)
where T is the current pose estimate with respect to the closest reference sphere which is determined from the previous iterations up to time t -1.

A maximum amount of pre-computation is performed offline during the construction of the spheres (e.g. Jacobian matrices and saliency maps) allowing the online algorithm to be computationally very efficient: the camera pose can be estimated at high frame rate (e.g. 45 Hz for a current image of 800 × 600 pixels in size).

To farther improve performance, a coarse-to-fine optimization strategy is employed by using multi-resolution spheres (e.g. constructed by Gaussian filtering and sub-sampling [START_REF] Burt | A multiresolution spline with application to image mosaics[END_REF]). The minimization begins at the lowest resolution and the result is used to initialize the next level repeatedly until the highest resolution is reached. This greatly improves the convergence domain/speed and some local minima can be avoided. Finally, to ensure robustness to illumination changes that occur between the reference images and the online camera's images, the method proposed in [START_REF] Meilland | Real-time dense visual tracking under large lighting variations[END_REF] is employed. This method combines both model-based tracking (w.r.t. the graph) and a non classic visual odometry approach, which greatly improve robustness to large illumination changes without necessitating the cost of estimating an illumination model.

V. AUTONOMOUS NAVIGATION

During autonomous navigation, the aim is to follow automatically a reference trajectory U generated locally around the learnt graph. The trajectory U = {u * 1 , u * 2 , . . . , u * n } contains n input vectors such that: where the point o * = {x * , y * } is a desired position, ψ * is the yaw angle, U * is the longitudinal velocity and ψ * is the desired angular velocity.

u * = {x * , y * , ψ * , U * , ψ * }, (7)
The control problem can be formulated as detailed in [START_REF] Benhimane | Vision-based control for car platooning using homography decomposition[END_REF]. In the proposed case, a virtual vehicle is followed and used to generate an error in translation and orientation that can be regulated using state feedback. Longitudinal velocity is controlled using a proportional feedback on the longitudinal error and steering angle depends on yaw and transversal errors. These errors are obtained by projecting the current vehicle position onto the closest reference trajectory's edge (cf. Figure 4). The reference position is then selected by translating the projected position along the trajectory by a distance d.

The translation error between the reference point and the current position is then defined by:

e q = e x e y = R T ψ * (o -o * ) = R T ψ * x -x * y -y * , (8)
where the rotation matrix of ψ * can be written by:

R ψ * cos(ψ * ) -sin(ψ * ) sin(ψ * ) cos(ψ * ) . (9)
The angular error is directly defined by:

e ψ = ψ -ψ * , (10)
and the control law, derived from [START_REF] Benhimane | Vision-based control for car platooning using homography decomposition[END_REF] without the state feedback on the velocities is:

U = U * -k x (|U * | + ǫ)e x ψ = ψ * -k y |U * |e y -k ψ |U * | tan(e ψ ) , (11)
where the gains k x , k y , k ψ and ǫ are positive scalars. An accurate 6 d.o.f. vehicle localisation is computed using only a monocular camera as detailed in Section IV. The obtained 6 d.o.f. pose is then converted into a 3 d.o.f. position (x, y, ψ) in the vehicle plane, and the resulting error is regulated using equation [START_REF] Hirschmuller | Stereo processing by semiglobal matching and mutual information[END_REF].

Since the localisation algorithm provides a 45 Hz pose estimation, no additional sensors or filtering are used for the positioning and the vehicle is able to smoothly follow trajectories. To ensure safe experiments, a SICK laser is placed in the front of the vehicle and it is only employed to detect pedestrians and other vehicles. Longitudinal velocity is modulated with respect to the closest objects detected in the laser trace and the vehicle is stopped under a certain distance threshold (cf. Figure 5(e)).

A. Results

The mapping and re-localisation methods have been experimented in large scale environments. The following result shows autonomous navigation in the city centre of Clermont-Ferrand in France, obtained during the final experiments of the CityVIP project (cf. Section VII). A learning phase has been performed on a 490 meters trajectory, by manually driving a vehicle equipped with the spherical acquisition system. In order to ensure an admissible path for the online trajectory (i.e. without obstacles), the trajectory obtained during the learning phase was used as input for the online navigation: however the localisation method is capable of accurately localising a camera navigating within a different region in the graph, as it is demonstrated in [START_REF] Meilland | Dense visual mapping of large scale environments for real-time localisation[END_REF]. The reference longitudinal velocity was set to 1.2m/s for the whole sequence.

Figure 5 shows the desired trajectory in black, and the trajectory followed autonomously by the vehicle in red. The vehicle starts at position (X=0,Y=0) and begins to move along Y axis. The experiment finishes at position (X=-62,Y=61). The vehicle was able to follow autonomously the whole sequence, using only the monocular camera for localisation. As it can be seen on figures (a),(b) and (c), the accuracy of the localisation method allows to navigate in narrowed corridors, whilst the employ of robust estimators ensure robustness to occlusions like pedestrians. Images (d) and (f) show the vehicle navigating in much larger areas (open place). This kind of environment has shown some limitations of vision based navigation. Since geometric information is far from the camera (building facade), accurate estimations of translations are degraded (infinite points are invariant to translations). These effects can be seen on the red trajectory around the landmark (e). However, this lack of precision could be overcomed using additional sensors, such as GPS and inertial measurements.

VI. CONCLUSIONS

The approach described in this paper propose a fully automated mapping and re-localisation system, for autonomous navigation in urban environments. The mapping method allows to reconstruct dense visual maps of large scale 3D environments, using a spherical graph representation. It has been shown that this representation is capable of reproducing photometrically accurate views locally around a learnt graph. Reconstructed spheres acquired along a trajectory are used as input for a robust dense spherical tracking algorithm which estimates the spheres' positions.

During online navigation, an efficient direct registration technique is employed to accurately localise a monocular camera. The robustness of the localisation method has been validated in challenging urban sequences containing a lot of pedestrians and outliers. The proposed localisation system is only based on a standard monocular camera, no additional sensors are used.

Future work will aim at improving the database construction, by geo-referencing the spherical graph in a GIS (Georeferenced Information System), which may contain higher level information such as free space. This geo-location will allow to use advanced path planning algorithms. To farther improve the autonomous navigation solution, it should be interesting to fuse vision based localisation, with GPS signal and inertial measurements.

VII. ACKNOWLEDGMENTS

This work has been supported by ANR (French National Agency) CityVIP project under grant ANR-07 TSFA-013-01.

VIII. APPENDIX

A. Non-linear optimisation

The error function for the real-time tracking (5) is minimized using an iteratively re-weighted least squared nonlinear minimization:

O(x) = arg min x ρ(e(x)), (12)
by ∇O(x)| x=x = 0, where ∇ is the gradient operator with respect to the unknown x defined in equation ( 2) assuming a global minimum is reached at x = x. An efficient second order minimisation approach is employed [START_REF] Malis | Improving vision-based control using efficient secondorder minimization techniques[END_REF], which allows to pre-compute most of the minimization parts directly on the reference image. In this case the unknown x is iteratively updated using a Gauss-Newton like optimization procedure:

x = -(J T DJ) -1 J T De(x), ( 13
)
where T is the transposition operator, J T DJ is the robust Gauss-Newton Hessian approximation. J is the warping Jacobian matrix of dimension n× 6. D is a diagonal weighting matrix of dimension n × n obtained by M-estimation [START_REF] Huber | Robust Statistics[END_REF] which rejects outliers such as occlusions and local illumination changes.

Fig. 1 .Fig. 2 .

 12 Fig. 1. Local representation: augmented sphere S containing intensities I S , depthmap Z S , salliancy W S and a sampling P S .

  (a) Spherical image. (b) Associated depthmap.

Fig. 3 .

 3 Fig. 3. Example of an augmented spherical panorama obtained using a multi-baseline spherical camera.

Fig. 4 .

 4 Fig. 4. Regulated error for visual servoing. The current vehicle position is projected onto the closest reference trajectory edge. A reference position is selected at a distance d, to generate longitudinal and angular errors {eq, e ψ }.

Fig. 5 .

 5 Fig. 5. Autonomous navigation. Top: A 490 meters trajectory followed autonomously, the reference trajectory is shown in black. The current trajectory is shown in red. (a),(b),(c),(d),(e),(f): images captured during the navigation.