Alice Lebois

Raoul Florent

Vincent Auvray

GEOMETRY-CONSTRAINED CORONARY ARTERIES MOTION ESTIMATION FROM 2D ANGIOGRAMS -APPLICATION TO INJECTION SIDE RECOGNITION

Keywords:

This paper deals with the 2D motion estimation of coronary arteries. It exploits the geometry of acquisition to strongly constrain the problem, thereby ensuring smooth and robust motion fields. The main contribution of this paper is to formally associate the C-arm angulation from which the heart is observed to the expected parametric motion field. Once this learning phase performed, every new sequence is associated to a fixed motion model, known up to a scale and an origin. The quality of the extracted velocities is visually evaluated over a large representative database. We also present an application of our approach to the motion-based classification of an angiographic exam into right or left artery tree.

INTRODUCTION

Many applications require or benefit from the motion estimation of the coronary arteries. Let us mention the tracking of perfusion basins along the heart cycle [START_REF] Liénard | Quantitative tool for the assessment of myocardial perfusion during X-ray angiographic procedures[END_REF], or the improvement of reconstruction from a rotational angiography [START_REF] Hansis | Projection-based motion compensation for gated coronary artery reconstruction from rotational X-ray angiograms[END_REF]. Depending on the application, it is either desirable to have a strongly regularized (smooth) or a data-driven (accurate) motion field. This topic has been mainly tackled in 3D [START_REF] Blondel | Reconstruction of coronary arteries from a single rotational X-ray projection sequence[END_REF][START_REF] Schäfer | Motioncompensated and gated cone beam filtered backprojection for 3-d rotational X-ray angiography[END_REF]. However, 3D data is available in a minority of cases only. This paper presents a robust framework to automatically estimate the motions of the coronary arteries from 2D angiographies. Our goal is to deliver smooth motion fields (to the possible cost of fine local accuracy), for instance for classification purposes. We are focusing on developing a particularly robust estimation scheme. It is mandatory for a majority of applications to run without manual interaction. Now, the first step of a motion estimation framework specifically dedicated to the coronary arteries is to extract or enhance the considered objects, yielding vessel maps with many false alarms when no user interaction is allowed. A robust estimation method is needed to extract reliable motion fields out of those noisy maps. To ensure the required robustness, we propose to constrain the problem by imposing a known parametric model (here, quadratic) to the motion field. Depending on the perspective from which the heart is observed, the 12 parameters defining the motion field are determined and fixed up to a scale and an origin. Indeed, the 2D observed motions of the coronary arteries is highly correlated to the geometry of acquisition, defined by the two angles of the C-arm. Our main contribution is to formalize this correlation and to exploit it to perform a fast and robust motion estimation of the considered coronary artery tree. We also present an application of our work to the motion-based classification of an angiographic exam into right or left artery.

GEOMETRY-CONSTRAINED PARAMETRIC MOTION ESTIMATION

The same arteries observed from the same perspective exhibit the same kind of motion. To formalize that assumption, we consider a representative database of 80 angiographies. Fig. 1 show their distribution relative to the two angles (φ 1 , φ 2 ) defining the acquisition geometry. We observe that they naturally gather in different clusters typical of the recommended reference diagnostic exams. This section explains how 2D motions representative for each of those clusters can be computed in an offline step (subsection 2.3), and how they can be adapted to a new sequence in an online step (subsection 2.4). The initial vessel map computation is presented in subsection 2.1. In this paper, the estimation is performed between the end diastole and the end systole cardiac phases. As a result, in the following, any sequence will be limited to two vessel maps only: V diast and V syst . This work can of course be applied to any two given instants of the heart cycle.

Vessel map computation

One can hardly rely on the original images I to estimate the motion of the coronary arteries since they are superimposed to the static background [START_REF] Auvray | Joint motion estimation and layer segmentation in transparent image sequences -application to noise reduction in X-ray image sequences[END_REF]. To separate them, we resort to a multiscale ridge filter R [START_REF] Frangi | Multiscale vessel enhancement filtering[END_REF]. To remove the background linear structures, the minimum value of the previous maps is subtracted from the ridge map:

V (x, y, t) = R I(x, y, t) -min t ′ ∈[0,t-1] R I(x, y, t ′ ) (1)

Problem statement

Our aim being to characterize the visual motion of the coronary arteries, we resort to parametric motion fields since they allow to describe a wide range of motions with a few parameters. We experimentally observed that each cluster C (gathering N C sequences) exhibits a dominant motion that can be described with a quadratic field θ C = (a C 1 ...a C 12 ) (12 parameters -Fig. 2). We mean that each sequence S C i in the cluster follows the considered motion, up to a scale parameter λ i and a translation (defined by the parameter x 0 i and y 0 i ). At pixel (x, y), the velocity along x reads:

v Si x (x, y) = λ i a C 1 + a C 2 (x -x 0 i ) + a C 3 (y -y 0 i ) +a C 4 (x -x 0 i ) 2 + a C 5 (y -y 0 i ) 2 +a C 6 (x -x 0 i )(y -y 0 i ) (2)
The velocity v Si y in y has a similar form.

Cluster parameter estimation

To estimate the parameters θ C typical of the cluster, we introduce the following energy to be minimized wrt θ C and the 3N C parameters λ i , x 0 i and y 0 i :

i∈C (x,y)∈Si V syst (x + v Si x , y + v Si y ) -V diast (x, y) 2 (3)
with v Si x and v Si y following (2). Under the assumption of small motions, the difference of displaced vesselness

V syst (x + v Si x , y + v Si y ) -V diast (x, y) (4)
can be linearized as a function of δθ C and δλ i . Solving (3) then boils down to solving a [N C + 12 x N C + 12] linear system. We can therefore jointly estimate the quadratic motion and the N C scales parameter in an efficient multiresolution framework.

It is not acceptable to similarly linearize (4) with respect to x 0 i and y 0 i , since those parameters often prove to be sizable (possibly > 100 pixels). We have chosen to perform an exhaustive search to estimate them, at each resolution level. Finally, the estimation algorithm is the following: at each resolution level

• Update the quadratic (θ C ) and scale (λ i ) parameters by minimizing a linearized version of (3);

• Adjust the translation parameters x 0 i and y 0 i by an exhaustive local search. Since the scales (λ i ) and relative shifts (x 0 i , y 0 i ) of the sequences in the cluster are estimated in the course of the algorithm, we are able to process a large variety of sequences without imposing constraints to the user (like centering the vasculature, using a given scale...).

Sequence parameter estimation

Once the learning phase performed, every new sequence is associated to a given motion parameter set θ C based on its angulation. We then simply have to estimate 3 parameters: the scale λ and the origin (x o , y o ) of the motion. We perform a multiresolution estimation as in subsection 2.3, alternatively updating λ (analytic expression after linearization) and (x o , y o ) (iterative research). Depending on the application, the estimated field can be directly used (if a smooth motion is necessary), or can initialize a more accurate motion estimation.

MOTION ESTIMATION RESULTS

The geometry-constrained motion estimation described in subsection 2.3 has been applied to 80 representative angiographic sequences, split into 10 clusters (left of Fig. 1).

Since the goal is to capture the visual motion, no explicit metric was involved. Instead, we have performed a qualitative assessment. We observed whether the estimated motion matches badly, correctly (the dominant translation as well as the muscle contraction is captured) or excellently (finer motion pattern are also captured). The overall result is shown on the right of figure 1. Globally, 90% of the sequences exhibit a correct (39%) to excellent match (51%) with the visual observation. Six examples of good matching are shown on Fig. 2 a)-f). They present a large variety of anatomies, scale, and vessel visibility. The coronary motions, though different in shape and amplitude, are well captured (including the rotational and shearing components), even though the framework is based on relatively noisy vessel maps. Motion estimation is incorrect in 10% of the sequences. In particular, the estimation performs worse in the cluster 5R (Fig. 1). It corresponds to right coronary arteries presenting a complex twisting without clear global motion, which cannot be captured by one single quadratic model (the velocity field is correct on Fig. 2 f) (little torsion), but not on Fig. 2 g)).

We could either introduce a more complex motion model, or consider that it is up to a motion refinement step to capture the observed twisting. The other poor estimations correspond to isolated sequences in the other clusters. Those sequences present motions that are actually different from the motions observed on the other sequences of the cluster (and indeed they often correspond to angulations at the limits of the considered cluster). More data would be needed to allow for a finer separation of the angulation space in more accurate clusters.

APPLICATION TO INJECTION SIDE RECOGNITION

The clinician typically observes the right and left coronary arteries successively. Being able to establish which side is being observed has a number of applications. For instance, one could position shutters and wedges according to the ongoing exam, which would reduce the dose. In our database, 2 of the 10 clusters correspond to right coronary arteries (17/80 sequences): the clusters 2R and 5R in Fig. 1. If the former can be clearly delineated based on the angulation, the latter corresponds to angulations also used for left coronary exams (cluster 2L). Since the motions observed in both classes are very different (they correspond to cases a) and b) on Fig. 2), we resort to motion to discriminate them. More precisely, we fit both considered quadratic motion models over each sequence using the method presented in 2.4. The class yielding the smallest residual is selected. This leads to excellent classification performance of 94% in a leave-one-out evaluation. The one erroneous classification result corresponds to nearly static left coronary. Integrating the expected amplitude of the motion is a perspective to further improve the classification performance.

CONCLUSION AND PERSPECTIVES

We propose a novel coronary artery estimation method that exploits the strong correlation between the observed motions and the geometry of acquisition. The C-arm angulation space can be split into regions where the expected (quadratic) motion model is known up to a scale and an origin. We report a good matching with the observations, and present a first application of this synthetic motion characterization for classification purposes. Futhermore, our approach characterizes the observed coronary motion with a small set of parameters. As such, they form excellent features allowing to integrate the global coronary motion in a machine learning framework.

We now need to explore to what extent our approach can improve the robustness and speed of motion estimation methods aiming at a finer accuracy. We intend to measure the accuracy obtained by other methods when the geometric constrained approach reported here is used as initialization.

Fig. 1 .

 1 Fig. 1. Database of angiographic sequences as a function of the two C-arm angles (φ 1 , φ 2 ). Left: gathering into 10 clusters (color-coded), according to the angulation and the observed coronary tree branch (left L or right R). The clusters '2L' and '2R' correspond to right and left coronary exams performed from the same angulations. Right: evaluation of the matching between the estimated velocity fields and the observed motion (light green: excellent; dark green: correct; red: bad).

Fig. 2 .

 2 Fig. 2. Different sequences together with the estimated generic field. From left to right: coronary artery at diastole; noncompensated vesselmap difference (illustrating the shape and the amplitude of the motion), estimated motion field (with the proper center and scale), motion-compensated vesselmap difference (illustrating the quality of the motion estimation). Clusters a) 2L b) 2R c) 4L d) 8L e) 8L f) 5R g)5R (see clusters on Fig.1).


Olivier Callot

FastVelo, a fast and efficient pattern recognition package for the Velo

This note describes the new implementation of the Velo pattern recognition package developed in autumn 2010 to cope with the changed requirements for LHCb tracking: No need for RZ track in the trigger, higher collision rate (so called μ) around 2.5 instead of 0.4 and the need for a faster algorithm for HLT1 to cope with the high rate and large occupancy.

After describing the main idea, this note describes the algorithm in details with the meaning of the control parameters. Some performance is also given.

"FastVelo, a fast and efficient pattern recognition package for the Velo" LHCb-PUB-2011-001

Introduction

The most important LHCb tracking detector is the Velo. This is where the position and direction of the produced particles are measured. The rest of the tracking system is there 'just' to measure the momentum. Reconstructing efficiently the tracks in the Velo and with as few as possible wrong tracks (so called ghosts) is thus a critical item for the overall LHCb performance.

The algorithm in use for the first data taking year at LHC was written by me in 2003 [START_REF] Callot | Velo tracking for the High Level Trigger[END_REF], updated in 2004 [START_REF] Callot | Online Pattern Recognition[END_REF] and 2006 and then maintained by David Hutchcroft with some adaptation to the updated geometry [START_REF] Hutchcroft | VELO Pattern Recognition[END_REF].

During the 2010 run, it became clear that the HLT doesn't need RZ tracks, because the Velo can't close completely and thus the RZ tracks are not as good as foreseen on Monte-Carlo for selecting tracks with high impact parameter. The second change is that, instead of working with a single interaction per crossing as foreseen in the preparation of the experiment, we took data with an average number of interactions per crossing (μ) of 2.5, and we intend to run with a value of μ in the range 1.5 to 2.5 in 2011. The code should then be revisited and tuned for this new context. Last, the HLT is saturating the computing resources of the online farm, partially due to the more complex events selected at high μ and an improvement in speed is needed to cope with the expected 1 MHz rate and will give some margin for other modifications in HLT.

During the autumn 2010, a new package has been written from scratch, of course with some inspiration of my previous implementations. This note describes the new approach, released in a package called Tf/FastVelo.

Basic principles

The Velo is made of two types of sensors, R sensor with strips at constant radius covering 45° in a so called sector, and Phi sensors with strips almost radial, grouped in two crowns. A schematic view of the Velo sensors is given in Figure 1 and the overall Velo is shown in Figure 2.

Reference frame

The real detector is not as perfect as in the design: The sensors are not absolutely perpendicular to the Z axis, their centres are not perfectly aligned, and the closing of the two halves is far from perfect, there are about 300 micrometers between the centres of the opposite side sensors. It is then very important to define properly in which reference frame the tracks are searched for. In the FastVelo implementation, one works always in the LHCb global frame. This is where the strip position is best known by the alignment. This frame is easy to use for the phi sensors, where the strip is a straight line. For the R sensors, where the strips are an arc of circle centred off-axis, one has to make an approximation in the first part of the code, where one works in R-Z projection: For each 45° sector, the 'R' coordinate of a strip is the distance along the bisector between the projection (perpendicular to the bisector) of the beam axis and the intercept of the R strip.

R-Z tracking

Tracks mostly come from the beam axis, or close to it, and they originate from the luminous region. In a first approximation, this means that tracks are also straight in the R-Z projection and that they come from a region in Z around 0 with a sigma around 5 cm. A first step of the pattern recognition is to find these (almost) straight lines in the R-Z projection, i.e. using only R sensors. But we are interested in B-decay daughters; these tracks come from a displaced vertex that is not on the beam axis. Similarly, and with larger effect, tracks from or Λ are not at all coming from the beam axis and thus are not straight lines in the R-Z projection. It should be noted that having spider web geometry, i.e. straight strips in 45° sectors, would have made the pattern simpler as straight lines (tracks) stay straight lines in this projection. But one can't change the geometry now! LHCb-PUB-2011-001 It was noted that most of the tracks we want are crossing at least 4 R sensors. In order to be faster, tracks with 4 hits will be searched first. Once all are found, tracks with 3 hits will be looked at, but using only unused clusters: Tracks made with 3 hits have a larger rate of ghost; asking for unused hits, i.e. for hits not already used on longer tracks, reduces the rate of ghosts. The tracking is done first for forward going tracks, then backward going tracks. The quadruplets are searched starting from the last four sensors, where tracks are most separated. Once a quadruplet is found, the track is extended towards the lower Z region as much as possible, allowing for some inefficiency, i.e. missing sensor. If the track is in one of the 45° sectors close to the vertical boundary between the two Velo halves, one searches for clusters on the opposite side that can be on the same track. If there are enough, they are added to the track which is then assumed to be almost in the vertical plane. Clusters are tagged as used and cannot be reused as seed for another quadruplet. The quadruplet is also allowed to have a missing sensor, so 4 out of 5 sensors are enough. Once all hits of the last sensor have been processed, one continues with the next sensor. The search for quadruplet stops when the 4 sensors cannot contain tracks coming from the luminous region, i.e. even at the maximum allowed angle the track would cross the Z axis before the luminous region.

A second step is to perform the same quadruplet search in the reverse direction, i.e. for tracks going towards negative Z. The code is the same, and one stops also as soon as tracks cannot come from the luminous region.

Then short tracks (triplets) are searched for, first backward then forward, requesting that the three clusters are unused. Tracks are also extended, both towards the luminous region, and if relevant in the overlap. As three points is not a very strong constraint to produce a straight line, triplets are compared and if two triplets share a hit, only the best one is kept. This R-Z tracking is in fact very similar to previous implementations, only the search for quadruplet first, plus tuning and details in the logic that will be described later, are really new.

Space tracking

Starting from a RZ track, one searches Phi measurements that can make a straight line. In previous implementations, the real phi coordinate was obtained using the radius in the phi sensor with some trigonometric functions, which is not very fast. The new approach developed in the FastVelo is to use the distance to the strip, a linear function of x and y. However, to select the hits to be considered, one needs an estimator of the distance of the hit plane of the bisector, as hits from the same track should group for this estimator. We use the sine of the angle to the bisector of the R sectors: The radius is known from the RZ track, and one computes the intersection of the phi strip with a circle (R-strip) centred on the sensor's centre. This implies to compute one square root, but there is no trigonometric function. As the track is in a 45° sector, the sine has to be smaller than about 0.4. Hits are stored per station: the phi sensors of the right side can measure coordinates for a track on the left side and vice versa. They are sorted according to this sine value. One searches for a triplet of hits in three stations having a close enough sine, and then performs a full space fit with the R measurements and these 3 first phi measurements. If OK, the track is extended in space as much as possible. The triplet is first searched in consecutive sensors, and then a missing sensor is allowed. Some cleaning is done to remove identical or very similar solutions.

Non pointing tracks

Some good tracks are not found by this procedure. David Hutchcroft (3) had implemented a recovery for non pointing tracks in the previous generation of code, called PatVeloGeneralTracking. A similar recovery for non pointing tracks has been implemented in FastVeloTracking. One tries to find tracks starting from phi measurements only, and only for forward going tracks: backward tracks are useful only to find the primary vertex, non-pointing tracks won't help. As there are fewer constraints in a phi tracking compared to the RZ tracking, this search is performed only on unused hits. A triplet of unused phi measurement has to be found in consecutive sensors. A rough test of compatibility is made using the centre of the phi strips: they have to be reasonably aligned with large tolerance. Then a first R hit is searched for in the R sensor corresponding to the first phi sensor. With these 4 measurements one gets a parameterisation of a space track. Hits in the two other R sensors are searched for. A good track requires 6 hits (3 R, 3 Phi) in 3 consecutive stations. At least 2 of the 3 R hits should be unused so far. The track is then extended as much as possible, using the space parameterisation. It should be noted that there is almost no pointing constraint in this search for track; the only request is that the first two phi hits are reasonably close.

Debugging

In order to understand the inefficiencies, or the ghost rate, it is important to be able to follow the processing of a selected MCParticle. The code should not have any reference to MC truth, but we have implemented a tool with a simple interface to select if a hit is from a given MCParticle. By job options, one can enable this tool with the option

FastVeloTracking().DebugToolName = "PatVeloDebugTool"

and selectively debug the processing of hits from the selected MCParticle by the option

FastVeloTracking().WantedKey = 2021

and follow what the code does, and why the track is not found. This tool also prints the truth information for a hit, in the form of a list of MCParticle keys.

An algorithm is added in the checking sequence of Brunel, called DebugTrackingLosses, which prints the truth information of missed track, or ghost, for a selected part of the tracking. This will give the event number and MCParticle key of the tracks one wants to debug. The following options are useful. Note that the last 4 are exclusive

DebugTrackingLosses().Velo = True DebugTrackingLosses().MinMomentum = 1000. DebugTrackingLosses().FromBeauty = True DebugTrackingLosses().FromStrange = True DebugTrackingLosses().Ghost = True DebugTrackingLosses().Clone = True

In fact, a large fraction of the code is made of debug statements, conditional to the detection of the use of a hit from a selected MCParticle. To use it, one should enable the decoding of MC information before running the reconstruction; this is done in Brunel by adding the option "Truth" in the RecoSequence before the "VELO" option.

Technical implementation

The package has one algorithm to perform the complete tracking. One tool handles the data; the track object handles also the fit. In this note, the main features are described, for all details on the code the only reliable information is to scrutinize the code itself!

In the following descriptions:

 ThisName is the name of a C++ class  ThisOption is the name of a job option, and is usually followed by the default value.

Hit Manager

As the representation of a measurement is different from the previous implementation, and to optimise memory usage and speed, a clean and simple implementation of the hit manager has

LHCb-PUB-2011-001 been coded. The hit manager has two purposes: Give access to a geometry description of each sensor, and give access to the hits of this sensor.

Geometry description: FastVeloSensor

The FastVeloHitManager owns a list of FastVeloSensor objects that can answer many questions about geometry, but also provides the hits of this sensor. The geometry information cached in FastVeloSensor is extracted from the Velo detector element. One large part of it is the parameterisation for the phi strips: This is 2048 times 5 numbers (3 for the line parameters, 2 for the strip centre) for 42 sensors. Access should be fast, and experience showed that one need to reduce the overall size, storing them as float and not as double. As this is pure geometrical information, it could eventually be included in the Velo detector element. The other parameters are the position and inclinations of the sensor itself, and the offsets for the R sectors, i.e. the projection of the centre of the sensor on its bisector.

Hit description: FastVeloHit

The hits are constructed from the VeloLiteCluster decoded from the raw data. The relevant geometrical parameters (position, weight = 1/ , strip centre, parameterisation of the distance) are stored in the hit, plus some working variables: One is the number of tracks using this hit, another is a working coordinate, typically a distance to something, used to sort the hits and get the best ones. For details, see the source files. Contrary to the previous implementation, the R and Phi hits are represented by the same class.

The tool itself: FastVeloHitManager

This is the tool used in the tracking to provide geometry and hits. It owns a pool of FastVeloHit that is reused from event to event. The geometry is constructed at initialisation, and an update is triggered whenever the geometry information (i.e. position in the LHCb frame) of a sensor is changed, this means when the Velo is open/closed. In fact one checks only on the first sensor on each side. The hits are cleared on the BeginEvent incident, and built from the VeloLiteClusters in a dedicated method.

The pattern algorithm: FastVeloTracking

The pattern recognition is a single algorithm that handles the sequencing of the track reconstruction as described earlier, and then stores all the tracks in the official container of LHCb::Track. Inside the algorithm, tracks are built and stored in a dedicated object called FastVeloTrack.

The internal track: FastVeloTrack

The internal track has basically two lists of hits, R and Phi. It has also a few parameters (zone in R sensor, backward flag, valid flag and missed sensors) and the temporary storage and final results of the fits, R-Z fit and full space fit. Several methods are implemented to access this information, add or remove hits, produce the state vector and the covariance matrix. Three methods are providing the transformation of an RZ track to a space track, with different input according to the needs. A second family of methods is related to the space fit: add or remove cluster, fit a space line, remove the worst iteratively until its χ 2 is below a specified value. A last useful method is to select the best cluster in a list to add to the track if good enough: This is the basic method to extend the track once it has a good seed.

The track fit

The fit of the space track is performed using the linear form giving the distance to the strip. For Phi strips, this linear form is constant; this is a property of the strip, more exactly of the cluster as we interpolate between strips. For R strips, this linear form describes the tangent to the strip at the proper azimuth. This means that this form should be updated whenever the parameterization of the track changes, i.e. when the azimuth changes. In practice, there is no need to iterate as the form changes in an insignificant way if the track moves by a phi strip pitch: As the minimal radius is 8 mm, a displacement of 50 micrometer will change the radius by less than 0.4 micron, totally negligible. But care should be taken to update all R hit parameters after adding Phi hits and before fitting a track.

A track has four parameters: x0, tx, y0, ty. The χ 2 to minimize is: ∑ which is:

∑ . .
. Computing the derivatives with respect to x0, tx, y0, ty gives four equations with four unknown, whose coefficients are sums of terms depending on a, b, c and z. The sums are computed and updated when adding/removing a hit. Several methods were tried to solve the system as fast as possible, including the recent Choleski decomposition method available in Root. The simplest solution by substitution was shown to be the fastest, around 5% of the total time compared to the Choleski decomposition.

RZ tracking

As explained earlier, two methods are used: Finding and extending quadruplet, finding and extending triplets of unused hits.

Finding quadruplets

One loops on the 'first' sensor, which would be the last one hit by the track in fact. For this sensor referred as 's0' one loops on the four 45° sectors, named "zone" in the code. It is assumed that a track doesn't cross a zone boundary, which is true for 99.9% of the tracks when the Velo is closed. For the open Velo, tracks crossing the boundaries are found by the thirst step, the 'unused phi' tracking.

Four cases are then processed one after the other 0. take 4 consecutive sensors, X X X X 1. Allow the fourth to be missing X X X . X 2. Allow the third to be missing X X . X X 3. Allow the second to be missing X . X X X

There is no need to allow the first to be missing; this is like using the next sensor...

For each unused hit in s0, starting from the innermost hit, define a range for the hit in s3: The minimum radius is given by the radius of the current s0 hit scaled with the maximum acceptable slope maxRSlope = 450 mrad to the z position of s3. The maximum radius in s3 is the radius in s0 plus one millimeter, this allows finding beam halo tracks, and also beam gas tracks coming from far. For backward tracks, the maximum radius is defined to find track coming from up to z = zVertexMaxBack = 1200 mm. This is enough to find beam2-gas tracks. The hits in sensor s3 are then scanned, starting from the saved first hit: This is the one that was the first inside the search window for the previous s0 cluster: When the radius increases in s0, the search windows in s3 shifts also, and hits in s3 that were below the window are definitely no longer to be looked at further. This gains some speed in the loop. The loop on s3 is aborted as soon as the hit is at a larger radius than the maximum. It should be noted that for the cases 1-3, i.e. with a missing sensor, the hit in s3 should not be used already.

From these two hits, one can predict the radius in s1 and s2, assuming the track is a straight line in the R-Z projection. The tolerance to find a hit is rMatchTol4 = 1.0 in units of the local pitch.

LHCb-PUB-2011-001

This means that the centre of the cluster should be within ± 1 pitch from the expected position. The cluster width is not taken into account. If hits have been found in s1 and s2, a valid quadruplet is found and a track created.

The track is extended towards its origin until either the extrapolation is outside the sensitive area of the sensor, or until there are more than maxMissed = 1 consecutive sensors without hits.

The tolerance to find a hit is defined as rExtraTol = 4 times the error on the extrapolated position, plus one pitch. After that, one checks that at least 2 hits are unused, to avoid finding again the same tracks. Clones should be killed as early as possible!

If the track is in zone 0 or 3, i.e. can be close to the vertical plane and thus have hits also in the other side of the Velo, one search for hits on the opposite side. The prediction is quite accurate, and the tolerance to collect hits is rOverlapTol = 1.0 in units of local pitch. If we have found more than 2 hits, all these opposite side hits are added to the track.

Figure 3 Location of the gap sensor

If a forward track ends at or just before the "gap sensor" indicated in Figure 3, i.e. in an area where the z distance between sensors is big, one tries to extrapolate the track downstream, i.e. find hits in the last sensors. The experience shows that the lever arm of the extrapolation is big, and two short segments before and after the sensor in the gap are frequently built.

The track is now stored in the local track vector. If the track has more than minToTag = 4 hits, the hits are tagged as used.

Finding triplets

The procedure is very similar to the previous one for quadruplets. The difference is that we have only 3 hits to find, so only 3 cases to consider, all 3 hits should be unused and the tolerance for the central hit is rMatchTol3 = 1.1 times the pitch. If the track has only 3 hits after extension, one checks that the track has not missed too many sensors downstream i.e. doesn't stop too early in the Velo: The density of ghosts in short incomplete tracks is too high. There is also a special handling around the "gap sensor": We frequently find tracks with 3 hits in the last 3 sensors, a missing hit, and a single extra hit in the gap sensor or before. This is because the extrapolation tolerance is proportional to the distance in z of the extrapolation, which becomes quite large in this case. It was found that this extra far hit should be removed as it is almost always a randomly picked up hit.

A special operation is performed to clean up 3-hits tracks: If two tracks of only 3 hits share at least a hit, only the one with the best χ 2 is kept. The accepted tracks are stored in the same local track vector and the hits tagged if the track is long enough.

Space tracking

Each RZ track is processed in turn. In order to build the more constrained tracks first, the RZ tracks are sorted by length, and longest tracks are processed first. Any criterion of unused hit is then robust, as previous tracks are longer and have then very few ghosts.

A first processing step is to define the first and last Phi sensor to use. This is done by checking when the radius of the track (in the local sensor frame) is inside the silicon. Then Phi hits in this selected sensor range are collected, selected and stored per station. The grouping per station is needed to handle properly the overlap of the two sides of the Velo when closed.

Computation of the selection variable.

The key ingredient here is the computation of the selection variable, which is the sine of the angle to the bisector.

Figure 4 Computation of the intersection point

For that, and as shown on Figure 4, we start from the blue line, i.e. the R value at the central azimuth of the R zone. We want the intersection of the circle, passing by this point and centred at the centre of the sensor, with the Phi strip. This is the intersection of a line and a circle, a second degree equation, with two solutions. However, as the Phi strips are almost radial, only one solution is in the proper R zone. A small technical detail: We work in Cartesian coordinates, and to avoid accuracy problems, the resolution is done in the X coordinate for zones close to the vertical, and in Y coordinate for zone near the horizontal.

For RZ tracks in the overlap, i.e. having R hits on both sides, the reference azimuth is the vertical line, and angular tolerances are reduced.

Before resolving the equation, a fast filtering of hits is applied, as we process all the hits of the two Phi sensors of the station: The projection of the strip's centre to the bisector should have a minimal value. This excludes all the opposite side strips plus part of the same side but far from the zone, at the cost of a simple scalar product. A second fast filter is that the distance of the starting point to the strip, measured by the linear form of the strip, is within a reasonable range. This range is defined as the radius times PhiMatchZone = 0.410, except for RZ tracks in the overlap, where the angular tolerance is the radius times PhiCentralZone = 0.040.

The selected hits should have their sine inside ±PhiMatchZone (or PhiCentralZone if in overlap) to be selected. The z position of the hit is adjusted as we now know a good approximation of the coordinate of the track. This sine is used as sorting variable inside each station. Stations with hits are counted, and the search range is restricted by ignoring the first and last stations if they have no hits.

LHCb reference frame

Sensor reference frame R from RZ, at central azimuth of the R zone Phi strip Wanted point

Search for a triplet of hits

We start from the first station with hits, counting always in the opposite direction to the track propagation, i.e. for forward going track this is the last station along the track. Five cases are looked at for the 3 stations s1, s2 and s3: 0. Only unused hits on sensors X X X 1. All hits on sensors X X X 2. All hits on sensors X X . X 3. All hits on sensors X . X X 4. All hits on sensors . X X X When a track with a full length, i.e. as many Phi hits as R hits, is found at the end of a case, the loop on the five cases is aborted. For each hit in the first station s1, the hit with the closest 'distance' in s2 is selected, if below MaxDelta2 = 0.05. The 'distance' is the absolute value of the slope perpendicular to the bisector plane, i.e. the ratio of the change of distance to the bisector to the difference in z. This is a measure of the fact that the track comes more or less from the beam line. With these two points (which are space points) one compute the straight line joining them in space. For all hits in the station s3 the χ 2 of the hit to the track is computed, and the smallest selected. This is the square of the distance, computed with the linear form, multiplied by the weight, inverse of the square of the error. If the smallest χ 2 chi2 is below 4 times MaxChi2ToAdd = 40, a good triplet is found.

The track is updated as space track adding the 3 hits; the Z position of the R hits is adjusted as we know now x and y; the track is fitted, and the highest χ 2 is computed and has to be below MaxChi2PerHit = 12 for the candidate to be kept. This cut is tighter than the previous one, because now we have the complete and exact geometry, while the computation of the distance used in the first selection has approximations if all sensors are not perfectly aligned.

If in overlap, hits from the other Phi sensor of the same stations are added if compatible, i.e. their χ 2 is lower than MaxChi2ToAdd. This is because, due to the dogleg shape of the sensor, there is some overlap between the two Phi sensors of a station, even more as their Z is not the same so the radius at which the track crosses them is different. Having as many tracks as possible with hits on the two sides of the Velo is a must for the alignment of the detector,

Extend the track

Phi hits are searched in the next stations, until two consecutive stations have no hits, or until we exit the selected range of stations. Hits should have a χ 2 contribution less than MaxChi2ToAdd, except if the z distance to the previous hit is larger than 140 mm. This is again to handle properly the extrapolation in the large gap shown on Figure 4. Similarly, if in the overlap, hits from the other sensor of the same station are added. A final fit is made, removing the worst Phi hit until its χ 2 contribution is less than MaxChi2PerHit. If there are not enough Phi hits after, the candidate is discarded. This minimal number of Phi hits is FractionFound = 0.70 times the number of stations. A quality factor is computed, this is the χ 2 per degree of freedom, plus twice the fraction of already used Phi hits. If the quality factor is higher than MaxQFactor = 6, the candidate is discarded. A test is performed to ensure that the track crosses the R sensor in the proper zone, at most one incompatible hit is allowed. A final check is made that there are both R and Phi hits on the same side of the Velo, i.e. that a track is not made with R on one half and Phi only on the other half. At this point, the track is kept and stored as space track. The loop continues for the next hits in s1, then for the next case except if a full length track has been found after case 1. There is no loop on the first sensor, as we know in which sensor the track should start, and case 4 recovers the case where the first sensor is missing.

Missed overlap tracks

If a RZ track in the overlap has not been upgraded as space track, a first recovery is performed by taking all the selected hits of all the stations, and trying to fit a track. For certain track, and with the alternated orientation of the Phi sensors, it may occur that a track has Phi hits only in every other station. The previous sequence will fail as there is no case like that. As the overlap has a small angular window, the hit multiplicity is small. The track is fitted, and hits removed until the highest χ 2 is below MaxChi2PerHit. It should be noted that the method to fit-andremove first removes the worst hit in sensors with several hits, and only when all hits in multihit sensors have a good χ 2 it removes hits which are alone in their sensor. This helps in finding a track.

If this recovery fails also, this indicates that the 'overlap' track is in fact a coincidence of two tracks in the two halves having the same RZ projection. The input RZ track is then split in two RZ tracks, one with hits from the left half, one with hits from the right side, and these two tracks are space-extended as described before, as tracks not in the overlap.

Final cleanup

If several candidates are found for the same input RZ tracks, there is a possibility to have generated clones, even if the code tries to avoid it. A check is made to detect two space tracks sharing more than FractionForMerge = 0.70 of their Phi hits. If this is the case, the shortest track, or the one with the worse quality factor, is discarded.

Tracks with only 3 Phi hits, and where all these hits are already used, are also discarded.

If there are several candidates, a minimal number of Phi hits is computed: This is the largest number of Phi hits in the candidates, decremented by one if this number is higher than 4. The lowest quality factor is then searched within the track having enough phi hits. If the minimal length is higher than 3, then a tolerance of DeltaQuality = 0.5 is added to this minimal quality factor, to get the maximal acceptable quality for a track.

Tracks passing the criteria of length and quality are then filtered again: If the track has only 3 R and 3 Phi, the quality should be less than MaxQFactor3 = 3. The track is refitted, eliminating R or Phi hits with a χ 2 greater than MaxChi2PerHit. If less than 6 hits are left, the track is discarded. Else, the track is stored as final track!

Phi tracking with unused hits

This search is performed only if the number of already found tracks is less than MaxRZForExtra = 200, if we have found already a number of space tracks higher than half the number of RZ tracks or if there are less than 20 RZ tracks. These tests are there to ensure that we will not have too high a combinatorics. Only forward tracks are searched for.

One searches for 3 R and 3 Phi hits in 3 consecutive stations making a segment. The seed sensor is taken in the last 8 stations, starting from the most downstream one. For each unused hit in s0, for each unused hit in s1 in a compatible ring (the track cannot go from the outside ring to the inside ring, while the reverse is true), a first compatibility check is made: The distance in the XY plane of the centre of strip in s0 to the strip in s1 should be less than PhiUnusedPhiTol = 5 mm. This tolerance is increased to 20 if the event is clean, i.e. has less than 10 RZ tracks. For each unused Phi hit in s3, the distance in the XY plane to the average (s0, s1) strip centre is computed and has to be less than PhiUnusedSecondTol = 10 mm, 20 mm for clean events. From these 3 Phi clusters, one computes the 45° zone in the first R sensor r0. Hits in this R zone compatible with the Phi ring are scanned, and for each a track candidate is built: There are 4 measurements, enough to compute a space line. The hit closest to this space line in r1 is selected, provided its χ 2 is below 400. The zone in this R sensor is selected according to the track extrapolation: This is how tracks crossing R zone boundaries can be found. The track LHCb-PUB-2011-001 parameters are updated with this 5 th hit, and the best hit in r2 is added, provided its χ 2 contribution is below MaxChi2ToAdd. If two or 3 of these R hits are already used, the candidate is ignored. The track is then extended, in R and Phi, until the extrapolation is outside the sensitive area, or until more than MaxMissed = 1 hits are missed. This is to recover tracks not found earlier, one can be quite strict! A full fit is performed, removing R and Phi clusters until the highest χ 2 contribution is below MaxChi2PerHit, and if there are 6 or more hits left, the track is stored and hits are tagged as used.

Final storage

Before converting the working tracks to final LHCb tracks, clone cleaning is performed on forward tracks. For each pair of tracks, the number of shared hits is computed. If it is higher than FractionForMerge = 0.70 times the shortest track, these tracks are clone candidates. The fraction of shared R hits is computed, if this is below FractionForMerge of the smallest number of R hits, both tracks are kept: Two space tracks sharing the same RZ were already compared as indicated in section 3.4.5 . Else the shortest one, or the one with worst quality factor, is removed.

Finally, the track is converted to an LHCb track, setting the proper parameters. A first state is provided at the point of closest approach to the beam line if StateAtBeam = true, else the first state is at first measurement. For forward tracks, a state is also provided at the end of the Velo, for extrapolation. Contrary to previous implementations, the track fit is not redone, as it was shown that the reweighting factor used previously generated no real improvement, while it had a cost in CPU.

Performance

The performance of the code is tested on Monte-Carlo simulation, where the efficiency and ghost rates are easy to measure. Tests were also done on real data, but there one needs a physics signal (e.g. J/Ψ→μμ) to measure the efficiency, and these signal events were selected using the previous version of the code. Tracks from a clean decay channel that weren't found by the old code can't be found by the new one: The events were not selected!

On 8000 event at ν = 3, which is close to the 2010 data taking conditions of μ≈2.5, the following performance is measured.

Item

Refitting the track: FastVeloFitLHCbIDs

In order to be able to reproduce the track fit in DaVinci, it is required to provide the same state on a track, starting from the list of LHCb Ids. This is implemented in this tool, which is used in the refit of tracks. Basically, it re-generates the FastVeloTrack and fits it, and then adds the states to the LHCb track.

When processing a track, the first step is to get the list of FastVeloHit from the LHCbID, using a dedicated method of the FastVeloHitManager. The R zone is computed by a majority vote of the R hits. The FastVeloTrack is built from the R clusters, the Phi clusters are added, R parameters updated and the track fitted. The track states are then created with the same code as in the original algorithm. It is assumed that the previous states have been removed from the track. It has been checked that this procedure give results essentially identical to the original ones. The source of difference is that the R parameters are updated from a different status of the track. These R parameters are essentially the equation of the tangent to the R strips once the space track is built. As in the original pattern code there can be extra hits that are

Figure 1 Figure 2

 12 Figure 1 Layout of the Velo sensors

Table 1

 1 Performance of the FastVeloTrackingIt is clear that there is a serious efficiency improvement, 1 to 1.5 %, with less ghost, better purity and efficiency, and in a time shorter by one third.

		Old code	FastVeloTracking	Difference
	Ghost rate (track)	9.8 %	8.1 %	-1.7%
	Efficiency long tracks	96.0 %	97.4 %	+1.4 %
	Efficiency long tracks > 5 GeV	96.7 %	98.2 %	+1.5%
	Eff. Long B daughters	96.1 %	97.2 %	+1.1%
	Eff. Long B daughters > 5 GeV	96.8 %	98.1 %	+1.3%
	Eff. long good B daughters	97.5 %	98.6 %	+1.1%
	Eff. long good B daughters > 5 GeV	97.6 %	98.9 %	+1.3%
	Eff. long Ks/Λ	87.8 %	90.4 %	+2.6%
	Eff. long Ks/ Λ > 5 GeV	90.3 %	93.3 %	+3.0%
	Clone rate all long tracks	1.7 %	1.7 %	-
	Hit Purity all long tracks	99.4 %	99.6 %	+0.2%
	Hit efficiency all long tracks	95.9 %	97.7 %	+1.8%
	Total time per event	5.65 ms	3.77 ms	67%

3.

8 List of all properties of FastVeloTracking

	Property name	Default value	Comment
	OutputTrackName	"Rec/Track/Velo"
	OnlyForward	false	Don't search for backward tracks
	OnlyBackward	false	Don't search for forward tracks
	HLT1Only	false	Only quadruplets and space
	HLT2Complement	false	Triplets and unused phi
	MaxRZForExtra	200	Don't do unused phi if that many RZ tracks
	StateAtBeam	true
	ResetUsedFlags	false	Clear used flags, for reprocessing the same event
	ZVertexMin	-170 mm
	ZVertexMax	+120 mm
	ZVertexMaxBack	+1200 mm	Define the maximum angle for backward tracks
	MaxRSlope	0.450
	rMatchTol4	1.0
	rMatchTol3	1.1
	rExtraTol	4.0
	rOverlapTol	1.0
	MaxMissed	1
	MinToTag	4
	PhiMatchZone	0.410
	PhiCentralZone	0.040
	MaxDelta2	0.05
	FractionFound	0.70
	MaxChi2PerHit	12.
	MaxChi2ToAdd	40.
	MaxQFactor	6.
	MaxQFactor3	3.
	DeltaQuality	0.5
	FractionForMerge	0.70
	PhiUnusedFirstTol	5.
	PhiUnusedSecondTol	10.
	DebugToolName	""	Use "PatVeloDebugTool" for debugging
	WantedKey	-100	If ≥ 0, debug the MCParticle with that key
	MeasureTime	false	Measure detailed timing

Table 2

 2 List of job options for FastVeloTracking

"FastVelo, a fast and efficient pattern recognition package for the Velo"

"FastVelo, a fast and efficient pattern recognition package for the Velo" LHCb-PUB-2011-001

LHCb-PUB-2011-001 removed afterwards, the exact values are not guaranteed. But the differences are well below the accuracy of any measurement; the states have all their parameters equal to better than 10 -3 .

Tracking when the Velo is open

It is important to be able to find primary vertices when the Velo is open, with the same algorithm as the HLT should not be reconfigured when the Velo is closed, this would lose precious luminosity. The FastVeloTracking code was tested in simulation and on real data with the Velo open. Compared to the previous code, the performance is better in all aspects.

Item

Old

Tracking beam gas events

In the HLT, it is also important to select beam gas events, to measure the beam profile and normalize the luminosity. These events have vertices far from the luminous region. Due to the optics of the machine, the interesting vertices should not be too far, the zone ±1 meter around the nominal interaction point is enough. The use case in the HLT is to detect that an event has a good chance to have a beam gas interaction, and then run a dedicated instance of the Velo reconstruction with relaxed cuts (mainly the length of the luminous region) to get the maximal efficiency on these rare events, less than 1 percent of the interaction rate.

FastVeloTracking is efficient enough to get most of the tracks due to beam gas inside the standard reconstruction. These tracks come from the beam axis, so are found by the RZ tracking.

The minimal track angle gives the maximal vertex distance. This is infinite for forward tracks as one accepts track parallel to the beam line. It is limited to ZVertexMaxBack = 1200 mm for backward tracks. Not all backward tracks (beam2-gas) are found with the standard setting, as we limit the range of Velo sensor looked at to the first few, up to the luminous region. But the code should find a few tracks, and then indicates that there is a potential beam2-gas interaction, thus triggering the execution of a dedicated version of FastVeloTracking looking at all backward tracks.

The exact implementation of this beam-gas processing in the HLT is still being worked out by the expert, Plamen Hopchev.

Summary

An entirely rewritten version of the Velo pattern recognition has been written and released in autumn 2010, and will be used for the 2011 run. It has better performance in all aspects, and supports Velo closing and beam gas search.

I want to thank for useful discussion the members of the tracking group, in particular Stephanie Hansmann-Menzemer, Silvia Borghi and Kurt Rinnert. I had interesting exchanges with Plamen Hopchev and Jaap Panman for the beam-gas validation.


Laurent Sagart

PAN MORPHOLOGY IN PHYLOGENETIC PERSPECTIVE

A recurrent problem when building phylogenetic trees from sound changes is that while they are relatively easy to identify2 and their directionality is usually clear, the objects these changes take for their targets -sounds-are not highly language-specific: neighboring languages tend to have much overlap in their sound inventories, so that when a sound change reaches a language boundary, the risk is high that it will cross it, provided that there are bilinguals and that the target sound exists on the other side. For this reason, shared phonological innovations often do not reflect shared inheritance. Blust's sound-changebased phylogeny for Formosan languages (1999) assumes a ten-branch star shape because the phylogenetic signal is hard to extract from the phonological isoglosses of a group of related languages in contact.

Lexical innovations take words as their targets. As associations of a string of phonemes and a meaning, words are far more language-specific than sounds. A lexical innovation affects the meaning of a particular word. When a lexical innovation reaches a sharp language boundary, the target word will usually not be present on the other side: the change cannot cross the boundary, even if there are bilinguals. True, a word having undergone a semantic change can be borrowed by a neighboring language, giving the appearance of lexical change in the borrowing language. The confounding effect of borrowing can be minimized by excluding cultural vocabulary from consideration.

In my ( 2004), based on a pattern of phylogenetic compatibility among six innovative characters in the An numerals from 5 to 10, I presented a new higher phylogeny of An. This was slightly modified in my communication to 10ICAL and in my (2008). I reproduce the tree in Figure 1. The story told by this tree is that of a founder group crossing from northern Fujian, where the straits are narrowest and a mountain top in NW Taiwan can be seen with the naked eye: expanding from their first settlements in NW Taiwan southward along the west coast, reaching southern Taiwan and continuing their expansion along the east coast: from where a second founder group left to establish Austronesian-speaking colonies ancestral to PMP and Tai-Kadai in the Philippines and on the south China coast.

Panel on Reconstruction of PAn morpho-syntax and implications for the An settlement on Taiwan 11ICAL -June 22-25, Aussois, France

In this paper I begin to examine the contribution of morphological innovations to early Austronesian (Formosan) phylogeny, using the same compatibility approach as in my 2004: I will discuss three nested characters, each of which is clearly innovative. The Austronesian languages which show these characters form are distributed in three nested isoglosses, forming a clique (sensu [START_REF] Meacham | Compatibility methods in systematics[END_REF], like the six numerals in Sagart 2004. When a node dominates another node in a phylogenetic tree, the innovations which define those nodes should show a nesting structure in geographical space.

1. Three early Austronesian morphological innovations.

In terms of their readiness to spread across language boundaries, morphological changes are probably intermediate between sound changes and lexical changes: on the one hand they are associations of a string of phonemes and a morphological meaning, making them more language-specific than sounds; on the other they are typically shorter than words, making them more like phonemes. Moreover, at least when closely related languages are in geographical contact, affixal inventories are likely to be shared to a significant extent: this increases the risk that a morphological change spreading from across the boundary will find its target on this side, and cross. At the same time morphological changes are reputed to be relatively resistant to borrowing. 3 In this, morphological changes are like changes affecting the basic vocabulary. [START_REF] Starosta | Rukai-Tsouic: subgroup or treetop ? Paper presented at the Seventh International Conference on Austronesian Linguistics[END_REF][START_REF] Starosta | A grammatical subgrouping of Formosan languages[END_REF] accepted that phylogenies should be built from innovations. He had the insight that because of their resistance to borrowing, morphological innovations are choice phylogenetic material. He built a phylogeny for Formosan languages based exclusively on morphological arguments. He accepted that Taiwan is the Austronesian homeland: he also accepted that MP subgroups with Amis (Harvey; Reid). He regarded Philippine morphosyntax as highly elaborate and assumed the history of morphosyntax between PAn and the Philippine subgroup to have been one of gradual complexification. He insisted that if Rukai, the Formosan language with the simplest morphosyntax, had lost morphological processes found in other languages of Taiwan or the Philippines, lexicalized traces of these processes should be found in Rukai.

The directionality of morphological change is not easy to establish on general principles. Starosta took as his guiding idea the view that in the absence of lexicalized vestiges, the Formosan language with the simplest morphology -Rukai-must be the first to branch off. My approach is different. I accept that morphological complexity can decrease over time in a language, leaving only faint traces in the lexicon. Arguments from 'absence of vestiges' are a variant of argumenta e silentio (arguments from silence) and must be used with caution: absence of evidence is not evidence of absence. My guiding principle is the idea, established in my 2004, that the PAn homeland was in NW Taiwan and that the first languages to branch off were Pazeh and Saisiat (whether as a single subgroup or not). This is because they show none of the early An innovations in the numeral system, but have the long additive expressions that gave rise to the shortened numerals *pitu 'seven', *walu 'eight' and *siwa 'nine'. In general I regard morphological characters shared by Pazeh or Saisiat and at least one other Austronesian language not in contact with them as PAn. When applied to neutral focus markers, this yields the system of PAn focus markers in Table 1:

non-actor focus markers patient V-en location V-an instrumental, beneficiary etc.

Si-V actor focus marker mu-V/<um>V On the same grounds, PAn perfective aspect was marked by <in> infixed in the verb. and the second set as an innovation at the Pituish node of the tree, as Ross notes. The path by which the second set arose out of the first can be described as follows: first, the genitive singular *ni-Cia was reduced to *nCia: *s and *C would be hard to distinguish following *n, and *nCia was reinterpreted as *nsia. In an analogically-motivated change, the base was *Cia was levelled to *sia. Finally *nsia was reduced to *nia. This sequence of changes was completed in proto-Pituish. In Pituish languages outside of Atayal and Malayo-Polynesian, the history of third-person pronouns is one of replacement of the new set by demonstratives. These replacements took place independently in each branch, as shown by the fact that the third-person pronouns in these languages are non-cognate.

1.2. Loss of -en in perfective patient focus forms.

There is no principle reason why we should expect an asymmetry in the way aspect and focus marking affixes combine in early An verb forms. Yet in the An world the perfective aspect marker (PERF) *<in> and patient focus marker (PF) *-en are found attached to the same verb stem only in four West coast Formosan languages: Saisiat, Pazeh, Thao and Siraya. Elsewhere such forms are not found.

Here is a Saisiat example [START_REF] Zeitoun | The temporal/aspectual and modal systems of the Formosan languages: a typological perspective[END_REF]: For Pazeh, Li and Tsuchida (2001: 28;38 n. 17) state that "In most Formosan languages and western Austronesian languages, the perfective form of the patient-focus does not bear the focus affix -en, but the infix <in>. Pazih may bear only the the affix <in> (e.g. b<in>aket 'to have been beaten'), -en (e.g. 'to be beaten or to have been beaten'), or both affixes (e.g. b<in>aked-en 'to have been beaten') to indicate the perfective aspect of the Patient-focus." They do not provide sentence examples.

hiza
Verb forms with both <in> and -in (< *-en) are grammatical in Thao (Blust 2003:238). Examples: in-dahip-in 'was helped', in-fari-n 'was blown by the wind', lh-in-irik-in 'was poked or pierced ', sh-in-umshun-in 'was Siraya, an extinct language of the SW coast of Taiwan known through 17th-century missionary materials, has been studied by [START_REF] Adelaar | Grammar notes on Siraya, an extinct Formosan language[END_REF]. Adelaar gives the past tense of undergoer-oriented verbs ('PF') as either ni-V or ni-V-ən, where past tense ni-is clearly the reflex of the PA *<in> (although it is prefixed to the verb stem, not infixed in it) and -ən is the reflex of *-ən. Examples of ni-V-ən from [START_REF] Adelaar | Grammar notes on Siraya, an extinct Formosan language[END_REF]:

ni-sulat-ən da PAST-write-PF ? 'It has been written...' ni-patimxa-ən tin ta vare vaung-appa PAST-punish-PF (by) him TM wind sea-also 'he rebuked the winds and the sea' Are we dealing with a retention from PAn or with an innovation by four western Formosan languages ? Under Blust's ten-branch phylogeny, an innovation is strongly indicated, since the four languages fall into two primary branches of An: a retention would take eight independent losses in order to account for the form's absence outside of the Western plains and Northwest Formosan branches. Predictably, [START_REF] Blust | A Note on the Thao Patient Focus Perfective[END_REF], to whom the Siraya, Pazeh and Saisiat facts were moreover not available at the time, came strongly on the side of the innovation view, arguing that in Thao PF -in (< *-en) could have been added to PF-PERF forms originally without it, in order to avoid homonymic clash with AF-PERF forms where the AF marker had a zero allomorph: if it was not for -in, he argued, nothing would distinguish AF-PERF forms like d<in>uruk 'stabbed' and f<in>ariw 'bought' from the corresponding PF-PERF verbs. He noted that among Thao verb stems in perfective patient focus, the PF marker -in (< *-en) is optional in some and obligatory in others; he suggested there is a tendence for those in which it is obligatory to have a zero allomorph of the AF marker <um> (like d<in>uruk 'stabbed' and f<in>ariw 'bought'), and to be optional in verbs where the allomorph of the AF marker is not zero (like s-m-iraq 'to kiss', s-m-in-iraq 'kissed' whose PF-PERF form is either s-in-iraq or s-in-iraq-in); he warned that the correlation is not perfect.

But the homonymic clash problem he outlined only exists under the view that -en is not original in Thao PERF-PF verbs. If loss of *-en is the innovation, Thao never lost it: it could have permitted zero allomorphs of <um> precisely because *-en was still there. Blust's discussion of this issue does not provide evidence that PF-PERF verb forms were without *-en in PAn; rather, it is a speculative account of why *-en could have been added to PF-PERF verbs in Thao assuming PF-PERF verbs were without that suffix in PAn.

In fact, consideration of the Saisiat, Pazeh and Siraya facts clearly shows that avoidance of homonymic clash is not what accounts for the presence of *-en in Pazeh b<in>aked-en 'to have been hit', Saisiat in-tani-in 'stopped' and Siraya ni-kita-(ə)n 'was seen', for these verbs require full allomorphs of the AF marker in their respective languages: Pazeh m<in>u-baket 'to have hit' (AF-PERF); Saisiat m-in-tani 'stopped' (AF-PERF), Siraya ni-k<m>ita 'saw' (AF-Panel on Reconstruction of PAn morpho-syntax and implications for the An settlement on Taiwan 11ICAL -June 22-25, Aussois, France PAST).

Because <in>V-en verbs are found in Saisiat and Pazeh as well as in two other languages -Thao and Siraya-, I regard them as part of PAn, and treat their absence in the rest of the An world as innovative. Judging from their low frequency of occurrence in Pazeh and Saisiat and their optionality in Siraya, it is likely that *-en was already optional in PAn when *<in> was present, in other words that loss of *-en in PF-PRF verbs was already in progress in PAn.

The motivation for this change appears to have been simplification. Removing *-en allowed to decrease morphological marking on verbs forms which were presumably very frequent in discourse. Since the *<in> perfective marker only occurred in verbs marked for verbal focus, removal of *-en could be achieved without losing a distinction, by taking advantage of a gap in the system.

1.3. Extension of ki-prefixation to verb stems.

In the following discussion, I rely in part on a draft paper by Stacy Teng and Elisabeth Zeitoun4 . There attention is drawn to a prefix *ki-that attaches to noun stems, deriving verbs with the broad meaning of "get N" (or "collect N", "harvest N", "cut N" etc.) in Saisiat, Kavalan, Rukai, Kanakanabu, Saaroa, Paiwan and Puyuma (as well as in Bunun where it assumes the specialized meaning of 'take off, remove'). Thao has a prefix kin-which serves the same function; it is unclear how it relates to *ki-. No Amis examples can be found in [START_REF] Fey | Amis dictionary. The Bible Society in the Republic of China[END_REF], and Teng and Zeitoun accordingly do not mention Amis among the languages having such forms; yet Zeng (1991) and Pourrias and Poinsot (ms) list some.

Here are some illustrative word pairs: Saisiat kaehoey 'tree, wood, brushwood' vs. ki-kaehoey 'gather brushwood', Kavalan tamun 'vegetable' vs. qi-tamun 'pick vegetables', Kanakanavu tamemi 'sweet potato' vs. ki-tamemi 'gather sweet potatoes', Rukai (Mantauran) paiso 'money' vs. i-paiso 'earn money'5 Paiwan sudju 'sweetheart', ki-sudju 'go courting, look for a sweetheart', Puyuma daqiŋ 'a share' vs. ki-daqiŋ 'claim one's share' Amis6 runaŋ 'mud pool' vs. ki-runaŋ 'wallow in mud' (of water buffaloes).

As expected of verbs which incorporate their object, denominal ki-verbs are one-argument intransitives: their unique argument is a nominative agent. Most cannot cooccur with any voice affixes. Teng and Zeitoun argue that derivation of denominative verbs by *ki-is a good candidate for PAn ancestry. This is likely from my point of view because this process is seen in Saisiat. The affix may have grammaticalized out of a verb meaning 'get, collect, gather'. A denominal verb *ki-pañay 'harvest rice' was probably part of PAn (Saisiat ki-pazay 'harvest rice', Paiwan ki-paday 'harvest rice'). 7Teng and Zeitoun describe a second construction involving the prefix *ki-, reflected in a more limited collection of languages. There the *ki-prefix attaches to a verb root to derive another verb which can be described as 'get V-ed': having reflexive or middle voice. This process is seen, with semantic variation, in Paiwan, Puyuma, Rukai, Bunun, Kavalan and, again, Amis.8 It is not seen in Saisiat, Pazeh, Siraya, Thao, Tsou, Kanakanabu and Saaroa. In

Pass-beat

Nom Kivi 'Kivi got beaten.' As Teng and Zeitoun noted, the nominative argument or subject in these constructions is the verb's patient. This is in contrast to verbs derived by ki-out of nouns (see above). Moreover Teng and Zeitoun contrasted ki-passives in Puyuma and Paiwan with patientfocus constructions in the same languages, observing that the ki-passives (where X gets V-ed according to his/her wishes) are more volitional than the Patient-focus construction. In a 'tentative conclusion' they further suggested that the rise of ki-passives is what triggered the realignment of verbal morphology around an active/passive distinction in Rukai:

"In Rukai, ki-V expanded so much that it came to replace the N[on-]A[ctor]V[oice] affixes while losing its volitional feature; the earlier NAV affixes were preserved in their nominalizing functions" While this conclusion was offered tentatively, it is attractive. Several Walu-Siwaish languages aside from Rukai have undergone wide-ranging morphosyntactic realignments: competition between verbs with ki-and NAF (Non-Actor Focus) constructions has the potential to explain why Puyuma *ki-passives have replaced the old neutral NAF constructions with affixes *-en, *Si-, *-an; 9 while the NAF "projective" constructions with *aw, *-ay and *-anay and with them, the entire NAF category were able to maintain themselves.

Ross (this panel) argues that Puyuma never had verbal uses of *-en, *Si-and *-an because no vestiges of verbs carrying these affixes can be found in that language. Here is a possible counter-example. Several languages of the northern Philippines: Isneg, Agta, Casiguran Dumagat etc. have a verb reflecting *[qʔ]unik 'to climb. 10 Puyuma (Cauquelin) has a stem qunkun 'to jump over', 11 eligible for agent and non-agent focus marking: munkun or muthe ouside by some object' ('中間有東西擋住') 9 *-an was maintained in its nominalizing function. 10 Isneg ʔumuneʔ, Casiguran Dumagat ʔunek 'climb up a tree', Agta ʔimunek [START_REF] Reid | Philippine minor languages: Word lists and phonologies[END_REF]. The last vowels in those forms reflect *i (Lawrence Reid, p.c., June 2009 13 Puyuma unkun can derive from *qunik-en, the PF form of *qunik. I suggest that when neutral patient-focus constructions were abandoned for ki-verbs, the PF suffix of *qunik-en was incorporated into the stem, making it eligible for all focus constructions and moving stress to the right. *quniken then underwent unstressed vowel syncope to *qunken and eventually (q)unkun, the modern form. I speculate that vowel syncope helped the stem to survive in recognizable form. Verb stems having incorporated *-en normally had three syllables: this made them prone either to be discarded as too cumbersome, or to have their first syllable pruned off to make them disyllabic. Pruning would however have a side-effect: to make them hard to etymologize. 14Here we see how morphological simplification can fail to leave conspicuous vestiges in the lexicon.

Like Puyuma, Paiwan developed a ki-passive construction but competition between NAF constructions and ki-passives did not lead to the loss of certain NAF markers as in Puyuma or of the entire NAF category as in Rukai. Likewise, in Amis and Bunun, the entire NAF category is maintained. Amis, Bunun and PMP settled the competition in favour of the NAF constructions by adding prefixes like ma-, pa-, pi-and mi-before ki-verbs, thereby giving those verbs the grammatical properties attached to these prefixes. A gradation can be seen among the prefix-adding languages: Bunun retains many unprefixed ki-verbs, Amis only a few: the MP languages keep none. This tree conflicts with the tree in Figure 1 in that the *enem isogloss in the numerals-based tree overlaps with the Loss of *-en isogloss in the morphology-based tree: thus Atayalic is above Thao and Siraya in the lexical-based tree, but below them in the morphology-based tree. 15 The problem can be fixed by allowing *-en to be lost twice: once in Atayalic and another time in a node corresponding to Enemish in the numerals-based tree. This is not unrealistic for a character that was already optional in PAn.

It is not surprising, in any case, that the two trees are entirely compatible, given the likelihood of contact effects and of independent innovations. What matters is that they tell very similar stories: a PAn homeland on the northwest coast of Taiwan, full settlement of the west coast before southern Taiwan and finally the east coast are settled, and an origin of the MP migration in the East coast languages. This view, now supported by seven 16 lexical and three morphological innovations, is consistent with archaeological dates -earlier on the west coast than on the east coast of Taiwan-. Population genetics (Sanchez-Mazas, this panel) give a similar picture.

15

It is impossible to be certain that extinct Favorlang did not allow *-en in at least some perfective patient-focus verb forms. 16 The Enemish node is also supported by the displacement of *kawaS by *CawiN as 'year' [START_REF] Sagart | The higher phylogeny of Austronesian and the position of Tai-Kadai[END_REF]. Why, then, do [START_REF] Gray | Language Phylogenies Reveal Expansion Pulses and Pauses in Pacific Settlement[END_REF], who work with the basic vocabulary, find a star-like phylogeny -although unlike in Blust's, and like in mine, PMP appears as coordinate with Paiwan-? this is undoubtedly because of their choice of both Old Chinese and Tai-Kadai as outgroups. Probably no Austronesianist thinks that both Tai-Kadai and Chinese are simultaneously related to An and exterior to it. Yet this is the basic assumption that underlies the analysis of Gray and his colleagues. By forcing Tai-Kadai, really a subgroup of An, to assume the role of an outgroup (that is, a language family outside of Austronesian but related to it), they are forcing their Bayesian statistics algorithm to treat every item shared by Tai-Kadai and any Austronesian language as a retention from PAn, thus preventing it from making use of post-PAn innovations like the 'standard' numerals between 5 and 10 or the shifts in 1st--and 2nd-person pronouns.17

9/11

I predict that if Gray and his colleagues give up Tai-Kadai as an outgroup and allow it to place itself in the tree rooted in Old Chinese, Tai-Kadai will branch off the Austronesian tree near where Malayo-Polynesian branches off, and the phylogeny they will find for Formosan languages will be very similar to those in Figure 1 and Figure 2.

This paper is dedicated to the memory of Stanley Starosta

  Panel on Reconstruction of PAn morpho-syntax and implications for the An settlement on Taiwan 11ICAL

  of PAn morpho-syntax and implications for the An settlement on Taiwan 11ICAL -June 22-25, Aussois, France

1. 1 .

 1 Third-person pronouns replacedRoss (2006:536-537) described two sets of third-person pronouns in Formosan languages: one, reflected in Pazeh and Saisiat, is formed on a base *Cia. The nominative third-person pronouns in Pazeh and Saisiat, sia and sia (singular and plural not distinguished) reflect that base. Singular and plural forms are distinguished in the genitive: singular ni-sia, plural n-asia in both languages, reflecting *ni-Cia (singular) and *ni-a-Cia (plural). The second set is reflected in proto-Atayal and in PMP: only the singular forms are cognate. The base form is *sia, the genitive *nia. In terms of the tree in Figure1, the first set can be regarded as PAn Panel on Reconstruction of PAn morpho-syntax and implications for the An settlement on Taiwan 11ICAL -June 22-25, Aussois, France

  Panel on Reconstruction of PAn morpho-syntax and implications for the An settlement on Taiwan 11ICAL -June 22-25, Aussois, France Paiwan, Puyuma and Rukai one sees a further development to passive verbs derived from verb s by *ki-. Teng and Zeitoun warn that this development may involve contact or parallel innovation. The following examples are reproduced from Teng and Zeitoun:

  Figure 2: Higher An phylogeny based on three morphological characters. (1): replacement of 3rd-person pronouns; (2) loss of *-en in perfective PF; (3) extension of ki-prefixation to verb roots.

  Panel on Reconstruction of PAn morpho-syntax and implications for the An settlement on Taiwan 11ICAL -June 22-25, Aussois, France

Table 1 :

 1 PAn neutral focus markers

  worshipped by someone' etc.Panel on Reconstruction of PAn morpho-syntax and implications for the An settlement on Taiwan 11ICAL -June 22-25, Aussois, France 'my foot was pierced by a slender bamboo some time ago'

	i-nay a	hulus f<in>ariw-in suma
	this LIG	shirt buy-PERF-PF someone
	'Someone has bought this shirt'
	aki	sh-in-umshum-in	iza	ita
	God worship-PERF-PF	already	1Pi-NOM
	'we worshipped God'
	nak	a	kuskus d-<in>uruk-in	sa	shkish	kahiwan
	1s-GENLIG	foot pierce-PERF-PF	by	slender-bamboo	some-time-ago

  ). 11 Listed by Cauquelin under unkun. Panel on Reconstruction of PAn morpho-syntax and implications for the An settlement on Taiwan 11ICAL -June 22-25, Aussois, France qunkun (AF), unkun-ai (PF), unkun-aw (LF). 12 Examples from Cauquelin's dictionary:

	ku	unkun-aw	na	gung
	1S	jump-over-LF the	ox
	I jump over the ox
	tu	unkun-ai	ku	Da	suan
	3S	jump-over-PF 1S	the	dog
	The dog jumps over me
	If the semantics are not judged too divergent,

Thanks go to John Wolff and Lawrence Reid for useful discussion and information.

All phonological mergers are innovations.

I make a distinction between transfer of a feature by spread and by borrowing: I speak of borrowing when a feature which is already fixed in a language -not in the process of spreading-is transferred from that language into another.

The passive ki-in Rukai, Paiwan and Puyuma: borrowing, shared innovation or parallel development? shown me by E. Zeitoun in August 2005.

This example shows that the process was still productive during the Spanish occupation of Taiwan.

This example is from[START_REF] Siqi | Taiwan Ameisi-yu yufa [A Grammar of Amis of Taiwan[END_REF]

See Sagart (2004) on the reconstruction of the PAn phoneme traditionally referred to as *j.

Zeng (1991:29) cites the pair ʔadiŋ 'keep off, fend off, shelter from' (遮擋) vs. ki-ʔadiŋ 'to be protected from

The variety of Puyuma investigated by Cauquelin loses *q word-initially. Initial q-in this word is preserved in one of the AF forms thanks to the mu-prefix. The other AF form munkun is analogically motivated.

Semantically the main difference between 'jump over' and 'climb' is the ballistic character of the former.

Puyuma has a number of disyllabic verbs ending in -un and without cognates in other languages, which are candidates for incorporation of *-en and loss of first syllable: repun 'to assemble, get together', compare Puy. reprep 'to swarm with, be infested with' (insects); kuLun 'to roll', compare root -kul 'curl, bend'; LuDun 'to sink', compare Isneg allad 'sink', Tiruray eled 'sink' (forms from ACD, under *qeled 'sink'; Blust mentions a rot -led 'sink').

Note that the choice of TK as an outgroup distorts the tree only above the point where TK really branches off: the non-Formosan part of the phylogeny in[START_REF] Gray | Language Phylogenies Reveal Expansion Pulses and Pauses in Pacific Settlement[END_REF] does not suffer from that problem.