Jean Collomb
email: jean.collomb@univ-smb.fr

Pascale Balland

Pascal Francescato

Yves Gardet

David Leh

Philippe Saffré

Y Gardet²

D Leh²

Validation of a low inertia

Keywords: Heat transfer, Profiled channels, Conformal channels, Rapid Heat Cycle Molding (RHCM), Rapid molding, Injection molding

teaching and research institutions in France or abroad, or from public or private research centers.

First experimental results are presented Table 1 and Figure 2. For these experiments, the heating of the mold is carried out with a thermal power of 12kW and a circulation of water. The temperature setpoint is set at 80°C. The improvement of thermal performances, thanks to rectangular heating channels and low inertia mold design, is confirmed with a heating rate more than 4 times greater and a temperature difference on the cavity mold 2 times lower. Moreover, a finite element model of the innovative mold, presented Figure 3, is defined with consideration of technological aspects: heating power, regulation and pressure drop. A numericalexperimental confrontation is also presented in Figure 3 and highlights the good prediction of the numerical model. This comparison of a conventional massive mold and an innovative low inertia mold with rectangular heating channels clearly confirms the improvement of thermal performances and the interest of this concept. With the numerical model, a thermomechanical multiparameter geometric optimization of this mold concept will be conducted in a future study.

Figure 1 :

 1 Figure 1: Massive reference RTM mold (left) and innovative low inertia mold with rectangular heating channels (right).

Figure 2 :

 2 Figure 2: Thermal gain of innovative mold for a production cycle compared to the reference mold.

Figure 3 :

 3 Figure 3: Numerical and experimental temperature comparison on the cavity surface for the innovative mold.

Table 1 :

 1 Experimental comparison of thermal performances

		Heating rate (°C.min -1 )	∆𝑇 𝑚𝑎𝑥 𝑐𝑎𝑣𝑖𝑡𝑦 𝑚𝑜𝑙𝑑 (°C)
	Reference massive mold	1.9	11.5
	Low inertia mold with rectangular heating channels	8.1	5.0


Introduction

Reactive polymeric systems are widely used in foam blowing or in resin transfer molding (RTM) processes. Numerical modelling of reactive processes is very difficult because the flow kinematics is coupled with the chemical reacting kinetics, which depends mainly on the reaction time (time elapsed since the reactants' mixing). The kinematics -chemical kinetics coupling -is a key point in foam blowing, being less important for example in usual RTM processes, where the chemical formulation can be adjusted in order to finish the mold filling before a significant rise in the resin viscosity.

The polymerization reaction increases the material viscosity and the expansion induces the material growth, generating a porous structure. When the polymerization takes place much faster than the foam expansion, the material becomes too rigid before the conclusion of the expansion. In this case the cellular walls can be broken with direct consequences for the mechanical properties of the conformed pieces. On the other hand, when the expansion finishes before the complete polymerization, the low material consistency is not enough to preserve the final geometry, and the structure collapses when the expansion is finished.

The numerical modelling of forming processes involving the flow of chemicallyreacting polymeric foams requires taking into account the different problem scales. Thus, in industrial applications a macroscopic approach is suitable, whereas the macroscopic flow parameters depend on the cellular structure: porosity, size shape and orientation of the cells, cellular walls properties, etc. Moreover, the shape and orientation of the cells are induced by the flow during the foam expansion (induced anisotropy). We have analyzed the influence of process conditions on the cellular shape. Elongated cells were obtained by applying a traction or a shear in the material during its expansion [START_REF] Chirivella | Experimental analysis of the induced anisotropy in foam forming processes[END_REF]. This microscopic information can be introduced in a macroscopic model by means of a homogenization technique. Thus, at each point in the material domain we can associate a characteristic volume containing some cells, whose size, shape and orientation depend on the point considered. Then, from the matrix and gas viscosities, and the cellular structure, we can compute an equivalent viscosity tensor. However, some difficulties appear in the homogenization procedure: (1) sometimes surface tension cannot be neglected, which requires a very accurate description of the inclusion geometry as well as consideration of its effect on the homogenized model; (2) the behaviour of the matrix fluid becomes quickly non-Newtonian, and its non linear character makes difficult the application of usual homogenization techniques; (3) the reaction kinetics remain uncertain for usual industrial processes; ( 4) the numerical tracking of a cell (microscopic level) which grows and deforms into the material flow is today an open problem, which in spite of some promising works [START_REF] Sussman | A level set approach for computing solutions to incompressible two-phase flows[END_REF](Sussman et a/., , 1999;;[START_REF] Sussman | An efficient interface preserving level set redistancing algorithm and its applications to interfacial incompressible fluid flows[END_REF][START_REF] Sussman | A coupled level set and volume-of-fluid method for computing 3D and axysymmetric incompressible two-phase flows[END_REF], important difficulties persist (Paredes, 2001); and (5) the physical mechanism of gas diffusion, the cells interaction and fusion, and the modelling of the final stage of the process, when the material is concentrated on the cell walls, are, among many others, important and unsolved difficulties found in the global modelling of foam growing.

In this way a macroscopic modelling where the reaction kinetics and the viscosity evolution are obtained experimentally is suitable from the point of view of its industrial application. Thus, such a numerical modelling could provide a visualization of mold filling (foam volume evolution and flow front position during the expansion reaction) and will be useful to optimize the mold and process design.

From that simulation we can determine the number, size and optimal position of the evacuation orifices to avoid internal overpressures (reducing the cleaning necessities as well as the material losses). Moreover, an accurate flow front tracking allows us to predict the interactions between different flow fronts when the curing reaction is too advanced (welding lines are fragile) as well as to avoid filling defects due to the big air bubbles retained between the flow front and the mold walls.

Simplified flow model

In order to obtain a macroscopic simplified flow model, able to carry out fast simulations of reacting foam flows, we consider the usual balance equations:

-The equilibrium equation neglecting the inertia and mass terms gives Divcr= 0 where cr is the stress tensor.

-The mass conservation

()p + Div( pv)= 0 dt - [I] [2]
where ~is the velocity field and p the foam density. The previous equation can be rewritten as

dp +pDivv+v•Gradp = 0 dt -- [3]
When the foam expansion is homogeneous p =p(t) and Gradp = Q. In this case

equation [3] becomes Divy_ = -()p _!_ = ~(t) dt p [4]
This function can be defined from the chemical reaction kinetics or from a simple experimental test. The first procedure, as previously discussed, remains very inaccurate because it is very difficult to know all the chemical reactions and their kinetics, in order to compute accurately the evolution of the density. Thus, we prefer, as proposed in Chinesta et a!. (Chinesta et a!., 2000), to proceed from a direct determination of the density evolution.

-The constitutive equation establishes the relation between the stress tensor a and the strain rate tensor D [5] where 2D = Grady+(Grady_'f. For compressible flows Tr( !2,)= Div~ ::f:. 0.

Equation [5] indicates also that the foam rheology depends on the reaction time. The simplest constitutive relation consists of [6] where A, denotes the volumetric viscosity, p the pressure field, ! the unit tensor and TJ the shear viscosity. If the volumetric viscosity contribution is neglected [START_REF] Chirivella | Experimental analysis of the induced anisotropy in foam forming processes[END_REF] where !2,' denotes the deviatoric strain rate tensor. The viscosity evolution of the matrix fluid can be measured using usual techniques [START_REF] Richter | Viscosity changes during the isothermal and adiabatic network polymerization[END_REF]. The most difficult task is the determination of the foam viscosity evolution (when the foam is considered as a pseudo-homogeneous fluid). The rheological characterization is difficult on account of the impossibility to proceed with standard rheometers which are not adapted to multiphase flows where the cell size is of the same order as the diameter or the gap of capillary or cone-plate rheometers respectively. Thus, some authors propose the use a of viscosity law in the form of a mixing rule that includes as control variables the concentration in reactants, the shear rate, the temperature and the gas volume fraction [START_REF] Lefebvre | Mathematical simulation of the flow of chemically-reacting polymeric foams[END_REF]. Other authors prefer to describe the foam rheology by means of a simple power law model [START_REF] Deshpande | Foam phenomena in sudden expansion and contractions[END_REF] which must be identified using a nonstandard device. The flow consistency and the behaviour indices of the power law equation could depend on time to describe the material polymerization.

-We can notice that thermal effects are not taken into account in the previous model in order to minimize the parameters which must be identified experimentally, although it is well known that the kinetics depends significantly on the temperature field (usually the determination of the reaction kinetics is made from the temperature evolution measured during the material curing).

Numerical modelling

When the foam expansion is heterogeneous, the density depends on the spatial coordinates, and the mass conservation equation is expressed by equation [3]. In this case the variational formulation of equations [1], [3] and [6] results in [8] and

J •( . A Gradp ) -p Dw!:-p+--•_!:' dil=O n(r) P [9]
The numerical discretisation of the flow model for an homogeneous foam expansion was discussed in a former work [START_REF] Chinesta | Evaluation d'un modele numerique simplifie de !'expansion de mousses[END_REF] From the kinematics point of view the discretisation of the flow equations governing the heterogeneous foam expansion, equations [8] and [9], do not introduce particular difficulties, if we can assign to each point in the foam domain and for each time the corresponding values of p, /..., ' 11 and p. Usually, the volumetric viscous term can be neglected, p derives directly from the density p , whose time evolution is obtained following the procedure described in the previous section. Finally, in some simple rheological models, the foam shear viscosity depends on the equivalent strain rate as well as on the time, due to the polymerization and blowing kinetics.

In this form, for each point in the fluid domain at any time, we need to know, for the fluid located at that point, the time elapsed from the beginning of the chemical reactions (reaction time). Knowing that time we can assign the density and viscosity values.

In the case of homogeneous expansions we assume that the reactants are dumped in the mould simultaneously, and in consequence the reaction time is the same for any point, which corresponds to the usual Eulerian time when we consider that the reaction starts at t = 0 . However, the manufacture of big pieces requires a great amount of reactants that must be distributed in large areas. In this case, we can not neglect the time elapsed in the reactants' dumping, so that, sometimes, the blowing starts before the complete reactant dumping, giving rise to a heterogeneous foam expansion. In that case, we need to label each point of the reactant fluid volume with its dumping time Td , and transport this field during the foam expansion. Thus, we must solve the following linear advection equation

1 dT dT -=-+vGradT=O dt dt - T( !_, t = 0) = Td ( !_) [10]
For example, if we suppose that reactants are dumped over a plate of one meter, from the left (x = 0) to the right (x = 1), in one second, the time label Td could be for example Td = X -I . If we take the origin of time ( t = 0) when the total amount of reactants is dumped, the reaction time t r will be given by t r = t -T .

When an Eulerian technique, using a fixed mesh, is adopted to solve equation [10] some difficulties appear, related to the field updating in elements whose filling stage is just starting. To illustrate this fact we consider an empty finite element at time t , which starts its filling process. In this case, as we are using a fixed mesh strategy, the field T updating in an element ge, i.e. the computation of Te at time t + l'!t , makes use of its value at time t . Because this element starts its filling process at time t, the field T(t) is not defined. If the simulation is carried out without paying attention to the different field initializations, great deviations appear, mainly in the neighborhood of the flow front, whose incidence cannot be neglected. For example, in the welding of several flow fronts we require active chemical kinetics to assure the material cohesion. In order to predict this kind of material defect, we need carry out simulations taking into account the hyperbolic character of advection equation, keeping good accuracy in the flow front treatment.

The main original contribution of this work is the establishment of a new general strategy able to carry out the field transport accurately even through the flow front. This technique minimizes the numerical diffusion and implicitly transports all fields with the flow front movement. The proposed technique can be also applied to general models involving advection or advection-diffusion equations.

Transport problems

As discussed in the last section, a transport problem appears associated with the advection of the field T. which is governed by the following linear advection equation [ 1 0]. If we consider a first order discontinuous finite element discretisation then [11] where re (t + L1t) is the value ofT in the element Qe at time t + L1t, nthe fluid volume corning to the element ne from the upstream element ne-(where the value ofT is given by Te-), n+ is the fluid volume leaving the element Qe and l nel is the volume of ne.

A first problem appears suddenly if we consider the time t for which an element Qe starts its filling process from its upstream element Qe-. To illustrate this limitation we consider the simplest situation of an incompressible fluid, i.e. p = 0. In this case, with n+ = 0 , we obtain [12] which establishes that solution at t +!it is dominated by the solution existing in the element at time t, at which this element is empty, and consequently the property T is not defined. The question is: what is the value to assign to T before the beginning of its filling process?

An accurate and consistent discrete form of the conservation equation can be obtained by considering a discontinuous finite element discretisation of the equation which results from the addition of equation [10] multiplied by the fluid fraction I, whose evolution is governed by the advection equation (Chinesta eta!. 2000) di a1

-=-a +vGradl=O dt t - [ 13]
with equation [13] multiplied by the field T, as proposed in Chinesta et al. [START_REF] Chinesta | Some di fficulties in the flow front treatment in fixed mesh simulations of composites forming processes[END_REF]

d(IT) = CJ(IT) + v Grad(IT) = 0 dt dt - [14]
which results in [15] where n~ is the volume of fluid coming to ne from the upstream element n,ei ;

Conclusions

From these results we can draw conclusions about the accuracy and low numerical diffusion of the numerical strategy proposed in the discretisation of equation [ 1 0]. The results are significantly improved, in relation to standard discontinuous finite element discretisations, mainly in the neighborhood of the flow front In this way the proposed strategy can be used to compute accurately the reaction times in welding areas in 3D foam expansion process simulation, to predict defects in the material welding.

Figure I . 0 "Figure 2 .

 I02 Figure I. Foam presence and velocity field during the expansion and welding


Joris Guerry

Alexandre Boulch

Bertrand Le Saux

Julien Moras

Aurelien Plyer

David Filliat

Aurélien Plyer

SnapNet-R: Consistent 3D Multi-View Semantic Labeling for Robotics

come

L'archive ouverte

Introduction

Visual scene understanding is a key capability to let intelligent robots evolve and interact in their environment. After decades of progress, we have reached a point where the acquisition of complex 3D scenes with fine details is now possible through multiple ways: precise laser scanners, commodity Red-Green-Blue-Depth (RGB-D) sensors or reconstruction techniques based on stereo and Simultaneous-Localization and Mapping (SLAM). Theses kinds of reconstruction are sufficient for simple navigation and collision avoidance, but a new step is required for the robot to act purposefully: semantic analysis of the 3D data.

In this paper we present a novel approach for semantic labeling of the scene perceived by a robot. It is built on an efficient multi-view approach called SnapNet [START_REF] Boulch | Unstructured point cloud semantic labeling using deep segmentation networks[END_REF] which established state-of-the-art results for urban, remote sensing data such as semantic3d [START_REF] Hackel | Semantic3d.net: A new large-scale point cloud classification benchmark[END_REF]. We propose new multi-view sampling strategies and better image and geometry integration. But most importantly, we show how 3D structure reconstruction and 2D semantics can mutually benefit from each other.

3D ¨2D : Instantaneously, when the robot gets a single RGB-D capture of its environment, we generate new virtual views which are consistent with the 3D structure of the perceived scene. It improve 2D semantic labeling at training (by data augmentation) and at prediction (by voting procedure).

2D ¨3D : At reconstruction time after capture of multiple images, adding semantics to the images allows to filter the points used for 3D reconstruction and to obtain more precise point clouds, which are in turn easier to label with semantics.

Overall, the contributions of this paper are:

• SnapNet-R, an improved multi-view Convolutional Neural Network (CNN) with simpler image and geometry integration as compared to SnapNet [START_REF] Boulch | Unstructured point cloud semantic labeling using deep segmentation networks[END_REF];

• 3D-consistent data augmentation for 2D semantic labeling with state-of-the-art performances;

• massive 3D-consistent sampling of 2D virtual views for 3D semantic labeling of reconstructed point clouds, validated on challenging data.

The paper is organized as follows. Section 2 presents the related work on image and point cloud semantic labeling. The description of our approach can be found in section 3: the core method SnapNet-R is explained in section 3.1 while implementations for RGB-D data and point clouds are detailed in section 3.2 and section 3.3 respectively. Finally, we evaluate our segmentation method in section 4.

Related work

Semantic labeling consists in point-wise classification of 2D or 3D data. In 2D, it became popular in the 2000s with challenges such as the CamVid dataset [START_REF] Brostow | Semantic object classes in video: A high-definition ground truth database[END_REF] or the Pascal Visual Object Challenge [START_REF] Everingham | The pascal visual object classes (voc) challenge[END_REF], a task consisting in giving a pixel-wise class label. In 3D, the related task is to identify the class membership of each 3D point. It gained interest with the availability of laser scanners which made possible the acquisition of large point clouds for application such as urban cartography [START_REF] Haala | 3D urban GIS from laser altimeter and 2D map data[END_REF] or robotics [START_REF] Anguelov | Discriminative learning of markov random fields for segmentation of 3d scan data[END_REF]. In the following, we describe the most recent and efficient approaches for tackling this problem, focusing on neural networks and robotics.

In 2D, fully-convolutional networks (FCN) [START_REF] Long | Fully convolutional networks for semantic segmentation[END_REF] were a milestone in the field of segmentation with dense prediction: by keeping the 2D spatial structure all along the network, they offered a simple and highly efficient approach. They were followed by several encoder-decoder architectures which exploit the symmetry of input and output spaces: U-Net [START_REF] Ronneberger | U-Net: Convolutional networks for biomedical image segmentation[END_REF], SegNet [START_REF] Badrinarayanan | SegNet: A Deep Convolutional Encoder-Decoder Architecture for Image Segmentation[END_REF], etc. Alternatively, Lin et al. [START_REF] Lin | Exploring context with deep structured models for semantic segmentation[END_REF] obtained state-of-the-art results using a multi-scale patch-based architectures with Conditional-Random-Fields (CRFs). While these networks perform in the Red-Green-Blue (RGB) domain, robotics commonly use Depth (i.e. distance from the sensor, denoted by D) as an additional 3D information, as reflected by reference datasets such as Sun-RGBD [START_REF] Song | SUN RGB-D: A RGB-D scene understanding benchmark suite[END_REF]. To exploit RGB-D, Gupta et al. [START_REF] Gupta | Learning rich features from rgb-d images for object detection and segmentation[END_REF] proposed an object detection method based on Region-CNN [START_REF] Girshick | Rich feature hierarchies for accurate object detection and semantic segmentation[END_REF] with depth encoding and semantic labeling output. Recently, the FuseNet architecture [START_REF] Hazirbas | Fusenet: Incorporating depth into semantic segmentation via fusionbased CNN architecture[END_REF] combined an encoder-decoder with early integration of depth at encoding stage. Later developments then lead to the Multi-View Consistent network (MVCNet) [START_REF] Ma | Multi-view deep learning for consistent semantic mapping with rgb-d cameras[END_REF] which takes benefit of unlabeled frames captured by a mobile robot by warping semantic predictions on these images in a common reference view (for which labels are known) based on a trajectory estimate. With respect to these approaches, our method uses the depth to perform 3D-consistent data augmentation, which allows us to transform a FuseNet (chosen for its state-of-the-art performances) in a multi-view CNN. Unlike MVCNet, our method does not require to be trained on video sequences to extract real, nearby frames but instead it creates virtual views from single RGB-D images, and still imposes 3D consistence.

In 3D, designing the most discriminating features for training a supervised classifier has been a standard approach for long. For example, in [START_REF] Charaniya | Supervised parametric classification of aerial Lidar data[END_REF] expert features such as normalized height or luminance were selected and ag-gregated. Generic descriptors able to represent the points and their neighborhood were also proposed: for example the fast point feature histograms [START_REF] Rusu | Fast geometric point labeling using conditional random fields[END_REF] or the signature histograms [START_REF] Tombari | Unique signatures of histograms for local surface description[END_REF]. But now, by using a deep learning framework, representations and classifier are learned all at once. In this field, three approaches compete.

First, voxel-based methods use a voxelization of the 3D space to create 3D tensors in order to feed a 3D convolutional neural network (CNN) [START_REF] Lai | Unsupervised feature learning for 3D scene labeling[END_REF][START_REF] Wu | 3D shapenets: A deep representation for volumetric shapes[END_REF][START_REF] Maturana | Voxnet: A 3D convolutional neural network for real-time object recognition[END_REF], mainly for object classification. Following this idea of 3D encoding, which may be computationally expensive, Hackel et al. [START_REF] Hackel | Semantic3d.net: A new large-scale point cloud classification benchmark[END_REF] proposed to use local, multiscale 3D tensors around the actual 3D points. In the voxel semantic labeling task (which is slightly different from point labeling) of the ScanNet benchmark [START_REF] Dai | Scannet: Richly-annotated 3d reconstructions of indoor scenes[END_REF], the proposed baseline network predicts labels for a whole column of voxels at once according to the voxels' neighborhood. Cherabier et al. [START_REF] Cherabier | Multi-label semantic 3d reconstruction using voxel blocks[END_REF] reduce computational cost of this kind of approach by dividing the space and working on smaller 3D blocks.

Second, the multi-view strategy consists in applying neural networks to 2D tensors which are collections of images of the scene. For retrieving and classifying shape models of objects, Multi-View CNN (MVCNN) [START_REF] Su | Multiview convolutional neural networks for 3D shape recognition[END_REF] takes several pictures all around a 3D meshed object and then performs image classification using a deep network. The PANORAMA representation [START_REF] Sfikas | Exploiting the PANORAMA Representation for Convolutional Neural Network Classification and Retrieval[END_REF] introduces another trick: projections on bounding cylinders oriented following the 3 principal directions of the volume. SnapNet [START_REF] Boulch | Unstructured point cloud semantic labeling using deep segmentation networks[END_REF] randomly takes numerous views all around the scene, creates virtual RGB and geometry-encoded images, and process it through U-net networks. Our approach has common features with these works: we generate snapshots of the 3D scene in order to use a 2D CNN with images as input. But instead of assigning a single label per 3D shape as in MVCNN, we compute dense 3D point labeling. And unlike SnapNet, we directly process RGB-D data in a single network.

Third, point-based methods work directly on unordered point sets, using architectures with fully-connected and pooling layers instead of convolutional layers. Thus, Point-Net [START_REF] Su | Pointnet: Deep learning on point sets for 3d classification and segmentation[END_REF] can output classes for the whole 3D shape or perform semantic segmentation of a 3D scene. However, it lacks the ability to capture local context at different scales. To go around this drawback and improve generalization to complex scenes, PointNet++ [START_REF] Qi | Pointnet++: Deep hierarchical feature learning on point sets in a metric space[END_REF] introduces a feed-forward network which performs alternatively hierarchical grouping of points and PointNet layers optimisazion.

Approach: multi-view segmentation networks for 3D semantics

Our approach consists in 3D-consistent view generation for improving the quality of semantic labeling. It may be applied in two cases of semantic labeling : 2D (image and depth map) or 3D point cloud.

In the first case, the objective is semantic segmentation of a unique 2D image for which some depth information is available (such as RGB-D data captured by Microsoft Kinect or Realsense R200 devices). For a single image, we generate other views of the scene such that the views are geometrically consistent with the 3D nature of the observed scene. Thus, all views correspond to what would be seen of the scene from a different point of view. It may be considered as data augmentation, but it differs from standard data augmentation (such as crop, translations, flip...) because it does not benefit only from the information contained in the 2D image plane but can also extract knowledge from the 3D spatial structure through the set of different appearances.

In the second case, the objective is semantic labeling of each point of a point-cloud. We focus in this study on point cloud reconstructed from a sequence of RGB-D images or couples of stereo images, for example using SLAM [START_REF] Engel | Direct sparse odometry[END_REF] which have been extensively studied in the last years. Compared to the RGB-D mono-image case, we exploit the higher completeness of the point cloud, allowing us to massively generate viewpoints very different from original cameras positions. We extend the SnapNet approach [START_REF] Boulch | Unstructured point cloud semantic labeling using deep segmentation networks[END_REF] with a new network architecture and a two-step procedure to generate a labeled point cloud: first 2D labeling of RGB-D images generated from stereo images for filtering a large part of outliers and then, 3D labeling using SnapNet on the points.

SnapNet-R

In this section we present SnapNet-R, a novel version of the previous SnapNet framework [START_REF] Boulch | Unstructured point cloud semantic labeling using deep segmentation networks[END_REF]. The SnapNet approach consists in independently doing the 2D semantic segmentation of each generated image from the 3D scene and then efficiently re-projecting the per-pixel attributed label back to its corresponding point in the 3D point cloud. This method is based on three main parts : the 2D image generator, the semantic segmentation technique and the efficient 3D back-projector.

Concerning the 2D image generator, the sampling strategy to pick viewpoints and create virtual views has to be adapted to the available data and to the application. In a complete 3D reconstructed scene (given as a point cloud) the generator has to take into account the observation scale (e.g. getting closer to small objects like cars or far enough from big objects like houses in order to have enough context information) and the orientation angle (e.g. not to look toward the empty sky), but can globally evolve freely in the environment. However, such a strategy no longer holds in single RGB-D image scenes like those from the SUNRGBD dataset [START_REF] Song | SUN RGB-D: A RGB-D scene understanding benchmark suite[END_REF]. The virtual viewpoints must be chosen close to the original camera pose. Indeed, a RGB-D image contains very sparse information: only points directly seen from a single point of view are captured, and every azimuth or at-titude shift from this original pose can lead to significant holes in the newly generated image as seen in Figure 2.

For semantic segmentation, the original approach uses two SegNet [START_REF] Badrinarayanan | SegNet: A Deep Convolutional Encoder-Decoder Architecture for Image Segmentation[END_REF], one for each modality (RGB and geometric features), and fuses them by a residual correction module [START_REF] Audebert | Semantic Segmentation of Earth Observation Data Using Multimodal and Multi-scale Deep Networks[END_REF]. The geometric features are handcrafted: normals, local noise and depth to the virtual camera. Thus, this protocol implies three different steps : training the RGB expert, training the geometric expert and finally training the fusion module. Here, we replace all theses components by a FuseNet SF5, introduced by Hazirbas et al. [START_REF] Hazirbas | Fusenet: Incorporating depth into semantic segmentation via fusionbased CNN architecture[END_REF], which takes directly RGB and depth map as input and so only needs one training stage.

Finally, the process of back-projection of semantic segmentation information back into 3D point cloud stays unchanged.

RGB-D geometric data augmentation.

The SUNRGBD dataset is composed only of single view images (an overview is available in Figure 5). This implies a different sampling strategy than for a full dense point cloud scene. Besides, the SnapNet approach was first designed to work on few big scenes and extract as many images from each scene as wanted. But SUNRGBD is composed of 5285 training images and 5050 testing images where each one of them is actually a scene. Thus, we can't sample hundreds of samples in each scene without taking the risk of exploding the final number of generated images. Therefore, we chose a sampling strategy that only generates 5 images for each scene and thus results in 26425 images for training and 25250 images for testing.

As previously said, because of the very sparse aspect of the 3D point cloud generated from RGB-D data, the sampling strategy cannot be random. Bad examples of random selection of viewpoints in SUNRGBD dataset are available in Figure 2. Also, since we decided not to produce too many images per scene we arbitrarily defined 5 camera poses:

• the first viewpoint is the original camera viewpoint • from the first viewpoint we get 2 more viewpoints by changing the azimuth angle of +10°and -10°• from the first viewpoint we get 2 more viewpoints by changing the attitude angle of +10°and +20°.

This sampling strategy is illustrated in Figure 1. No uplooking viewpoints are picked because objects mainly stand on the floor and are best viewed from above. Besides, the scenes are inside-rooms so objects are at close range and thus, looking upper could make us miss them. Furthermore, the depth sensors suffer from bright surfaces like mirrors, windows, screens or glass objects and the measurement can be very noisy resulting in black artifacts in depth map and thus holes in the RGB images resulting from 3D re-projection. That is why we applied inpainting preprocessings : based on Navier-Stokes method [START_REF] Bertalmio | Navier-stokes, fluid dynamics, and image and video inpainting[END_REF] for the depth images and by using mean color of non-black neighbors on RGB images (before training and testing).

3D scene reconstruction and labeling

We focus here on standard scenarios for robotics, where a robot equipped with various sensors moves around, as in the 3D Reconstruction Meets Semantics (3DRMS) dataset (cf. section 4.1). The image sets for reconstruction of 3DRMS consist of RGB and Gray stereo pairs. A pose estimation is given for each acquisition. In this section, we propose three different reconstruction and labeling pipelines, each of them corresponding to a different robotic use case: from a SLAM fashion reconstruction to global multi-view estimation. These pipelines are presented on Figure 3.

Reconstruction and semantic from navigation

The first pipeline -Classif 2Drepresents a type of reconstruction that can be used in robotic navigation. Labels and points are accumulated for each acquisition image pair. This pipeline follows 4 steps:

Depth estimation: For each camera pair, we rectify the RGB/Gray images and the ground-truth label images. These rectified data will be used for the whole following tasks. Then, we use the Efficient Large-scale Stereo Matching (ELAS) algorithm [START_REF] Geiger | Efficient large-scale stereo matching[END_REF] in order to compute the disparity map and add a depth channel to the RGB images, resulting in raw RGB-D data.

2D semantic segmentation: Using the train sets of these rectified RGB-D data we train a Fusenet SF5 for semantic segmentation (with same condition as experiment 1 in Table 1). Despite the low number of training images, we are able to generate rough label prediction maps on all sets.

Point cloud reconstruction: For both train and test sets, by taking into account the given global position of cameras and the computed rectifications, we accumulate the points from all RGB-D images into a global coordinates system (a voxel grid with voxel size of 0.01 coordinate unit). Due to the small baseline between acquisition cameras, such an accumulation produce point clouds with a lot of outliers (see Fig. 4 up left). Therefore, we implemented a filtering strategy: first, for each image we filter the points that were too close (< 1m) or too far (> 3.5m) from the camera. The closest points are mainly due to the on-board vehicle used during recording and the farthest were too noisy because of the small baseline. Then, we filter each image-relative point set according to two filters: label based and geometry based. Based on labels provided by the 2D semantic segmentation, we discard unlabeled and background points and points in a 3 pixel margin between two different labeled objects. Based on geometry: we compute the normals using [START_REF] Boulch | Fast and robust normal estimation for point clouds with sharp features[END_REF] and reject points which are too aligned with camera observation axis (dot product between normal and direction vector from camera center to considered point > 0.7) A snapshot of the result on test set is presented in Figure 4 up right.

Direct labeling: like in SnapNet-R, the label at each point is computed by voting from the FuseNet semantic map predictions.

Reconstruction from navigation and 3D labeling

This method -Classif 2D-3Dis a variation from the previous pipeline. The reconstruction step is still done in a SLAM fashion but here the 2D semantic maps are only used for filtering. At inference time, to predict labels the former direct labeling is replaced by SnapNet-R. At training time, labeled training point clouds are computed using classif 2D where the predicted semantic maps are replaced by the provided ground truth.

SnapNet-R: a label is given to each point using the SnapNet-R framework. We adapt SnapNet-R to exploit the point cloud topology, more complete than RGB-D data. We compute new RGB-D pairs by randomly picking view points. The camera position is selected such that: a) it looks in the direction of the scene (particularly, the center of the image is directed towards an existing point) and b) it is positioned over ground level. For each camera direction, we take two snapshots distant of 2 and 5 meters from the targeted point.

Off-line reconstruction and 3D labeling

The last approach -Classif 3Dis purely post processing of the data. It is a two step pipeline : a global point cloud reconstruction and a 3D semantic labeling.

Global multi-view reconstruction: we used the Agisoft Photoscan software. By using all images at once for reconstruction, it allows to search for point matching between cameras with a greater baseline than reconstruction from stereo pairs. It yields point clouds with less estimation artifacts such as continuous transitions between objects in depth map estimation (cf. Figure 4, bottom).

SnapNet-R: we labeled the point cloud with the previously trained model (from section 3.3.2).

Experiments and results

Datasets

We validate our approach on three datasets used in robotics which rely on different types of sensors to get 3D information.

Sun-RGBD is a dataset of images captured by lowcost RGB-D cameras [START_REF] Song | SUN RGB-D: A RGB-D scene understanding benchmark suite[END_REF]. It combines previously existing smaller RGB-D datasets (NYUDv2 [START_REF] Silberman | Indoor segmentation and support inference from RGBD images[END_REF], Berkeley B3DO [START_REF] Janoch | A category-level 3d object dataset: Putting the kinect to work[END_REF] & SUN3D [START_REF] Xiao | SUN3D: A database of big spaces reconstructed using sfm and object labels[END_REF]) to reach a size of more than 10k images from 4 various RGB-D sensors. It is completed with several types of 2D and 3D annotations over the whole dataset which allows training and evaluating algorithms for various task such as scene classification, semantic segmentation or 3D object detection. All 37 classes are given in Fig. 7 but to name just a few there are: table, bed, chairs, curtain, fridge, mirror, sink, floor, wall, etc.. We also evaluate on NYUDv2 [START_REF] Silberman | Indoor segmentation and support inference from RGBD images[END_REF] alone, considering the 40 classes and 13 macro-classes segmentation tasks as in [START_REF] Eigen | Predicting depth, surface normals and semantic labels with a common multi-scale convolutional architecture[END_REF][START_REF] Long | Fully convolutional networks for semantic segmentation[END_REF].

The 3DRMS Reconstruction Challenge provides us with a dataset of image sequences captured by a robotic platform.

This lawn-mower-looking ground robot is equipped with two stereo rigs for RGB and grayscale images and a Leica Lidar sensor for point cloud acquisition. The dataset contains 5 sequences of calibrated images with the corresponding camera pose: 4 for training (totalizing 108 views) and 1 for testing (with 125 views). For the training part only, there are also ground truth 2D semantic annotations and a semantically annotated 3D point cloud depicting the area of the training sequence (cf. Fig 6). Classes are Unknown, Grass, Ground, Pavement, Hedge, Topiary, Rose, Obstacle, Tree, Background.

Evaluation criteria

Let C be the set of classes (labels) and X be the input points to classify. Let X c ⊂ X be the points classified with label c ∈ C. Finally, we note X * c the objective classification (i.e. the ground truth) for label c.

Overall accuracy: (OA) is the proportion of well labeled points:

OA = 1 |X| c∈C |X c ∩ X * c |.

It does not take into account the unbalance between classes but gives a good view of the global behavior of the classifier.

Mean Accuracy (M A) is computed with the per class accuracies :

M A = 1 |C| c∈C |Xc∩X * c | |X * c | .
Intersection over union (IoU ) penalise also false negative predictions :

IoU = 1 |C| c∈C |Xc∩X * c | |Xc∪X * c | .

Semantic segmentation of RGB-D data with 3D-consistent multi-view SUNRGBD:

The semantic segmentation task of SUN-RGBD is a hard task. The best method to our knowledge, the Context-CRF of Lin et al. [START_REF] Lin | Exploring context with deep structured models for semantic segmentation[END_REF], only achieves 42.3% of IoU. To this purpose they use a very rich and complex set of neural networks, united in a conditional random field framework. They first use a neural network to get a feature map of the input image, another pair of networks to predict patchto-patch semantic relation and at last, after bilinear upsampling of the low resolution prediction, they apply a boundary refinement post-processing with a Dense CRF [START_REF] Krähenbühl | Efficient inference in fully connected crfs with gaussian edge potentials[END_REF]. By comparison, the results we show in Table 1 do not make use of any refinement post-processing neither conditional random field inference. The complete processing pipeline of our strategy consists in wrapping several images from a single scene to 224x224 patches, infer a low resolution semantic map through a FuseNet SF5 for each of them, upsample each prediction map by nearest neighbor interpolation, reproject all information in the original point of view, and finally vote for the label to assign to every pixel.

Our fully implemented method corresponds to the experiment 9 in Table 1. We achieve a new state of the art result with 58.13% of mean accuracy on the SUNRGBD semantic segmentation task, compared to the 53.4% of the Context-CRF [START_REF] Lin | Exploring context with deep structured models for semantic segmentation[END_REF]. We also obtain 78.04% of overall accuracy against 78.4% for the Context-CRF and 39.61% IoU (vs 42.3%), which are the second best known results on this dataset.

Compared to the Fusenet SF5 framework alone, our multiview strategy leads to a pure gain of performance for all metrics, especially on mean accuracy going from experiment 1: 54.81% (52.61% in the original article [START_REF] Hazirbas | Fusenet: Incorporating depth into semantic segmentation via fusionbased CNN architecture[END_REF]) to 58.13%.

If the Context-CRF method was to be used as the basic Table 1. SUNRGBD [START_REF] Song | SUN RGB-D: A RGB-D scene understanding benchmark suite[END_REF] quantitative results. All numbered experiments are using a FuseNet SF5 trained on all SUNRGBD train data (on the contrary of [START_REF] Hazirbas | Fusenet: Incorporating depth into semantic segmentation via fusionbased CNN architecture[END_REF] who removed RealSense images for depth quality reason). We replicated the result of FuseNet SF5 in experiment 1 and applied our full multi-view strategy in experiment 9. For each criterion, best values are emphasized in bold, second best values in bold italics. Results are discussed in section 4.3.

segmentation method of our SnapNet-R approach, the results would be probably even better for each accuracy metrics.

Because all the information is going through a 3D step, we had to preprocess the depth maps not to loose RGB information. Thus we present the results of this processing in Table 1. Experiment 5 uses pre-processing at training and testing time and leads to a little improvement over the FuseNet SF5 alone (experiment 1). But experiment 6 shows that multiview testing reduce IoU performance (-9 points). This implies that multiview generated images are slightly different from original view images and must be seen at training time. Some qualitative results are shown in Fig. 7. SnapNet-R approach is less sensitive to occlusion situation (see row 5). Besides, looking at the chairs in row 6 shows a much better edge segmentation for SnapNet-R than for FuseNet alone. The poor results of FuseNet (row 1) can be explained by the depth map quality. In this specific situation SnapNet-R takes fully advantage of its inpainting processings.

We also tried to directly train at full resolution (experiment 10 and 11) by picking random crops of 224x224 pix-els. Even if this approach is not able to achieve the same performances as the low resolution full frame experiments, we observe the same improvement tendency by using our multi-view data augmentation.

NYUDv2 labeled dataset is extracted from RGB-D video sequences. Our method only generates artificial new data based on a single frame on the contrary of MVCNet method [START_REF] Ma | Multi-view deep learning for consistent semantic mapping with rgb-d cameras[END_REF] who uses more data with the adjacent video frames for training and testing. Results on NYUDv2 (see Table 2) show the same trend as for SUNRGBD. The SnapNet-R approach reacts well to more than 7 times fewer data at training time. About the 40 classes semantic segmentation task we outmatch the state of the art method M A by +8.77 points and reach the second best score for OA and IoU . Concerning the 13 classes task we achieve new state of the art results, reaching 81.95% OA (+2.82), 77, 51% M A (+6.92) and 61.78% IoU (+2.71).

3DRMS challenge

As the segmentation task is part of the ongoing challenge, the test ground truth is not available yet. Besides, the small quantity of training data would lead to over-fitting if we had practiced cross-validation on them only. We settle then for qualitative results on the testing set with models trained on the training sets.

Reconstruction. The figure 4 compares the point sets generated with stereo data without filters (up left), with outliers removal (up right) and using a global approach (bottom). We have to filter rough data before using SnapNet-R: small objects would not be seen otherwise. However stereo baseline is small and leads to depth estimation errors that are not corrected at accumulation. It results in multiple occurrences of small or thin objects such as tree trunks. In the case of 3DRMS challenge, global approaches provides smoother point cloud with less artifacts than incremental stereo. Thanks to normals estimation it is visible on Figure 8 through shadow which reflects local orientations.

Labeling. Visual results provided in Figure 8 represent the predicted classes with artificial colors. Compared to Classif 2D (left), SnapNet-R leads to more consistent labeling (Classif 2D-3D middle). This is due to the synthetic cameras which provide views of the point cloud inaccessible to the ground robot. As an example, views looking at the ground with a small incidence angle allows to capture the global geometry of the pavement and ground areas. The experiments also show the robustness of the approach for 3D point clouds. We trained the SnapNet-R on the 2D classification training point clouds. The regularization induced by the CNN and the prediction smoothing with multiple snapshots reduce the imperfect labels due to reconstruction artifacts.

Perspectives. The main perspective concerns the global reconstruction. In this paper, we used the models learned from the Classif 2D-3D pipeline. The next step is to train directly on reconstructed data and labels using the 3D classification framework without intermediary steps. As presented in section 3.2 dealing with RGB-D data, the closer the data you train on are to the test data, the better the results are. This gives us a glimpse of improvement for the Classif 3D pipeline.

Conclusion

We presented in this paper an extension to SnapNet for robotics applications. We changed the core neural network architecture and proposed a snapshot strategy adapted to RGB-D single view data. Using the later as data augmentation improve significantly the segmentation performances on SUNRGBD and NYUDv2, leading to new state of the art results.

We also applied SnapNet-R on the 3D Reconstruction Meets Semantics challenge data. The three pipelines described in the this paper correspond to different use cases from robotic navigation to offline reconstruction. The qualitative results are promising and illustrate the potential of multiview strategy for 3D data labeling.

Figure 1 .

 1 Figure 1. Illustration of the sampling strategy over single view data from SUNRGBD (real distances and angles are not respected for illustration purpose).

Figure 2 .

 2 Figure 2. Illustration of bad sampling strategies with single view data. Upper left: good sample close to real camera point of view, upper right: seen through dresser, down left: seen from behind, down right: seen from bedside lamp.

  Pipelines and processing blocks are described in section 3.3. Qualitative results on 3DRMS are discussed in section 4.4.

Figure 3 .

 3 Figure 3. 3D reconstruction and labeling pipelines for the 3DRMS challenge.

Figure 4 .

 4 Figure 4. 3DRMS challenge point cloud (test set), incremental accumulation without filter (up left), with filter (up right) and global reconstruction (bottom).

Figure 5 .

 5 Figure 5. Sample of the SUNRGBD dataset: the 3D semantized point cloud and on right column RGB image, depth map and grounth-truth (see Fig.7 for label legend).

Figure 6 .

 6 Figure 6. Overview of the 3DRMS Challenge dataset: training sequence "boxhood row" with the 3D semantized point cloud and on right column 2D semantic annotations over image, front-left RGB image, front-right grayscale image.

Figure 7 .

 7 Figure 7. Qualitative segmentation results on SUNRGBD [41]. The first three columns contain depth, RGB and ground-truth images. Last two columns present the results obtained by Fusenet SF5 [23] then SnapNet-R.

Figure 8 .

 8 Figure 8. Results of labeled reconstruction for the 3DRMS challenge, test set. Left, labels obtained with direct labeling (Classif 2D), middle is SnapNet-R applied on stereo reconstruction (Classif 2D-3D) and right is SnapNet-R applied on multiview reconstruction (Classif 3D).

Table 2 .

 2 NYUDv2[START_REF] Silberman | Indoor segmentation and support inference from RGBD images[END_REF] quantitative results. For each criterion, best values are emphasized in bold, second best values in bold italics. Results are discussed in section 4.3.


Y W A N G

C Wang

Z Z H A N G

Segmentation of clustered cells in negative phase contrast images with integrated light intensity and cell shape information

Keywords: Cell segmentation, clustered cells, grey-weighted distance transform, light intensity, negative phase contrast images, shape

Automated cell segmentation plays a key role in characterisations of cell behaviours for both biology research and clinical practices. Currently, the segmentation of clustered cells still remains as a challenge and is the main reason for false segmentation. In this study, the emphasis was put on the segmentation of clustered cells in negative phase contrast images. A new method was proposed to combine both light intensity and cell shape information through the construction of grey-weighted distance transform (GWDT) within preliminarily segmented areas. With the constructed GWDT, the clustered cells can be detected and then separated with a modified region skeletonbased method. Moreover, a contour expansion operation was applied to get optimised detection of cell boundaries. In this paper, the working principle and detailed procedure of the proposed method are described, followed by the evaluation of the method on clustered cell segmentation. Results show that the proposed method achieves an improved performance in clustered cell segmentation compared with other methods, with 85.8% and 97.16% accuracy rate for clustered cells and all cells, respectively.

Introduction

Observation and analysis of live cells play a vital role in characterisation of cell behaviours for both biology research and clinical practices, including wound healing, immunoreaction, and cancer cell metastasis [START_REF] Piccinini | Automated image mosaics by non-automated light microscopes: the MicroMos software tool[END_REF][START_REF] Das | Computational microscopic imaging for malaria parasite detection: a systematic review[END_REF]. In these applications, individual cells are segmented in optical images to implement cell morphological analysis, cell classification and even tracking [START_REF] Xing | Robust nucleus/cell detection and segmentation in digital pathology and microscopy images: a comprehensive review[END_REF]. Cell image segmentation is a process of separating individual cells from background in images containing cells. This could be done Correspondence to: Yuliang Wang, School of Mechanical Engineering and Automation, Robotics Institute, Beihang University, 37 Xueyuan Rd., Haidian District, Beijing 100191, China. Tel: +86-18612525756; fax: +86-10-82314554; e-mail: wangyuliang@buaa.edu.cn manually, which is tedious and time consuming. Therefore, automated cell image segmentation techniques become important for high throughout cell morphological and dynamics studies [START_REF] Irshad | Methods for nuclei detection, segmentation, and classification in digital histopathology: a review-current status and future potential[END_REF][START_REF] Wang | Quantitative characterization of cell behaviors through cell cycle progression via automated cell tracking[END_REF].

Different methods have been proposed to implement automated cell image segmentation, including thresholding method [START_REF] Chalfoun | Segmenting time-lapse phase contrast images of adjacent NIH 3T3 cells[END_REF], morphological erosion method [START_REF] Chowdhury | Cell segmentation by multi-resolution analysis and maximum likelihood estimation (MAMLE)[END_REF], watershed method [START_REF] Lin | A hybrid 3D watershed algorithm incorporating gradient cues and object models for automatic segmentation of nuclei in confocal image stacks[END_REF][START_REF] Wahlby | Combining intensity, edge and shape information for 2D and 3D segmentation of cell nuclei in tissue sections[END_REF][START_REF] Yi | Automated segmentation of multiple red blood cells with digital holographic microscopy[END_REF][START_REF] Kim | Nucleus extraction of uterine cervical pap-smears using marker information and watershed algorithm[END_REF], level set method [START_REF] Ersoy | Cell spreading analysis with directed edge profile-guided level set active contours[END_REF] and active contour models [START_REF] Chenyang | Snakes, shapes, and gradient vector flow[END_REF][START_REF] Cosatto | Grading Nuclear Pleomorphism on Histological Micrographs[END_REF][START_REF] Ali | An integrated region-, boundary-, shapebased active contour for multiple object overlap resolution in histological imagery[END_REF][START_REF] Seroussi | Segmentation and tracking of live cells in phase-contrast images using directional gradient vector flow for snakes[END_REF][START_REF] Vink | Efficient nucleus detector in histopathology images[END_REF], or the combination of several different algorithms. Most methods can achieve automated cell segmentation to some extent, especially for individual cells. However, the segmentation of clustered cells still remains as a challenge.

The methods used to segment clustered cells can be categorised into shaped-based [START_REF] Kumar | A rule-based approach for robust clump splitting[END_REF][START_REF] Al-Kofahi | Improved automatic detection and segmentation of cell nuclei in histopathology images[END_REF][START_REF] Kong | Partitioning histopathological images: an integrated framework for supervised color-texture segmentation and cell splitting[END_REF][START_REF] Wienert | Detection and segmentation of cell nuclei in virtual microscopy images: a minimum-model approach[END_REF][START_REF] Veta | Automatic nuclei segmentation in H&E stained breast cancer histopathology images[END_REF][START_REF] Chen | Accurate segmentation of touching cells in multi-channel microscopy images with geodesic distance based clustering[END_REF][START_REF] Song | Accurate segmentation of cervical cytoplasm and nuclei based on multiscale convolutional network and graph partitioning[END_REF][START_REF] Gharipour | Segmentation of cell nuclei in fluorescence microscopy images: an integrated framework using level set segmentation and touching-cell splitting[END_REF] and light intensity-based approaches [START_REF] Ali | Advanced phase-based segmentation of multiple cells from brightfield microscopy images[END_REF][START_REF] Wang | Segmentation of the clustered cells with optimized boundary detection in negative phase contrast images[END_REF]. Shape-based approaches use distance transform [START_REF] Al-Kofahi | Improved automatic detection and segmentation of cell nuclei in histopathology images[END_REF][START_REF] Chen | Accurate segmentation of touching cells in multi-channel microscopy images with geodesic distance based clustering[END_REF][START_REF] Song | Accurate segmentation of cervical cytoplasm and nuclei based on multiscale convolutional network and graph partitioning[END_REF], concave points along contours [START_REF] Kumar | A rule-based approach for robust clump splitting[END_REF][START_REF] Wienert | Detection and segmentation of cell nuclei in virtual microscopy images: a minimum-model approach[END_REF][START_REF] Gharipour | Segmentation of cell nuclei in fluorescence microscopy images: an integrated framework using level set segmentation and touching-cell splitting[END_REF], and radial symmetry transform [START_REF] Kong | Partitioning histopathological images: an integrated framework for supervised color-texture segmentation and cell splitting[END_REF][START_REF] Veta | Automatic nuclei segmentation in H&E stained breast cancer histopathology images[END_REF][START_REF] Gharipour | Segmentation of cell nuclei in fluorescence microscopy images: an integrated framework using level set segmentation and touching-cell splitting[END_REF] to conduct cell segmentation, assuming cells have regular (elliptical, for example) and convex shapes. Light intensity-based methods segment cells using the special distribution of light intensity over cell bodies, assuming the light intensity has direct correlation with cell height. The commonly used methods are level set method [START_REF] Ali | Advanced phase-based segmentation of multiple cells from brightfield microscopy images[END_REF] or peak detection method [START_REF] Wang | Segmentation of the clustered cells with optimized boundary detection in negative phase contrast images[END_REF].

Both shape-based and light intensity-based methods have pros and cons in clustered cell segmentation. On one hand, light intensity, like phase contrast images, represents height information of cells. As a result, the light intensity information is useful in determining cell locations and exact boundaries between two contacting cells within an area of clustered ones. On the contrary, the shape-based methods only use the boundary information, which is normally obtained from a binarisation operation. This inevitably loses some useful information of the clustered cells. On the other hand, the light intensitybased methods require high image contrast. However, the actual cell images, including phase contrast images, suffer from various noises. The methods may lead to oversegmentation or undersegmentation. From above analysis, one can see that light intensity and cell shape information should be integrated to perform a more efficient segmentation, rather than being utilised solely.

Some efforts have been put to combine the shape and light intensity information [START_REF] Lin | A hybrid 3D watershed algorithm incorporating gradient cues and object models for automatic segmentation of nuclei in confocal image stacks[END_REF][START_REF] Wahlby | Combining intensity, edge and shape information for 2D and 3D segmentation of cell nuclei in tissue sections[END_REF][START_REF] Ali | An integrated region-, boundary-, shapebased active contour for multiple object overlap resolution in histological imagery[END_REF][START_REF] Dewan | A method for automatic segmentation of nuclei in phase-contrast images based on intensity, convexity and texture[END_REF]. These methods normally rely on some assumptions, either from shape or from light intensity point of view. For example, in the gradientweighted distance transform method, it is assumed that internuclei gradients are higher than intranuclei gradients [START_REF] Lin | A hybrid 3D watershed algorithm incorporating gradient cues and object models for automatic segmentation of nuclei in confocal image stacks[END_REF]. In the synergistic boundary and region-based active contour model, nuclei are assumed to be roughly elliptical in shape [START_REF] Ali | An integrated region-, boundary-, shapebased active contour for multiple object overlap resolution in histological imagery[END_REF]. The assumptions limit their applications and reduce their robustness.

Aiming at developing a clustered cell segmentation method coping with the complexity of live cells either on shape or light intensity in negative phase contrast images, we focus on fusing the light intensity and cell shape information. It gives an improved performance for clustered cell segmentation. Specifically, the fusion of the two information solves two issues existed in negative phase contrast image segmentation. One is the pseudo light intensity peaks caused by intracellular organelles, which could cause oversegmentation. The other is the reversed image contrast for cells with larger height (like cells under division), which could cause undersegmentation.

The rest of the paper is organised as follows. In the Experimental section, experimental details is introduced, including cell culturing, cell imaging and programming software. In the section of Principle of the clustered cell segmentation, the detailed process of clustered cell detection and separation are presented. In the Results and discussion section, the clustered cell segmentation using the proposed method is demonstrated, followed by the comparison and evaluation of clustered cell segmentation among different methods.

Experimental

In this study, an established nontumorigenic breast epithelial cells MCF 10A (Michigan Cancer Foundation 10A) were cultured for the experiments. The cells were maintained in 47.5% Dulbecco's modified Eagle's medium and 47.5% F-12 medium supplemented with 5% horse serum, EGF (Epidermal Growth Factor) (20 ng mL -1 ), hydrocortisone (0.5 µg mL -1 ), Cholora toxin (0.1 µg mL -1 ), NaHCO 3 (0.2 nM), and 1% penicillin/streptomycin.

Negative phase contrast images were acquired using a phase contrast microscope (IX71, Olympus, Japan) equipped with a 10× negative phase contrast lens. The image size is 1344 × 1024, corresponding to the actual field of view of 900 µm × 686 µm.

All algorithms developed in this study, including preprocessing, suspected clustered cell detection and skeleton-based clustered cell separation, were implemented in Matlab (Version 2016b, The Mathworks, Inc., USA).

The principle of the clustered cell segmentation

In this study, the shape and light intensity information are integrated to segment the clustered cells. The schematic of the segmentation is shown in Figure 1. For a simulated image shown in Figure 1(A), the preliminary cell shape information is obtained by thresholding operation to the cell image. Then a binarised image which indicates the preliminary shape of cells is obtained, as shown in Figure 1(B). The boundary of the mask can then be used to represent cell shape. Meanwhile, light intensity of cell bodies within the preliminarily segmented mask area can also be extracted (Fig. 1C). After that, the light intensity information and distance transform information are fused (Fig. 1D).

For the sake of the clustered cell separation, a skeletonbased separation algorithm is employed. The Skeleton (green line in Fig. 1E) of the preliminarily segmented area represents the distribution of cells in the region. The mask area can then be separated with assistant of the graph theory. The final step is the refinement of cell boundary detection using a contour expansion method [START_REF] Wang | Segmentation of the clustered cells with optimized boundary detection in negative phase contrast images[END_REF], as shown in Figure 1(F).

Results and discussion

In this section, the detailed procedure of the proposed method is first demonstrated, followed by the evaluation of the proposed method.

Construction of grey-weighted distance transform

In this study, the grey-weighted distance transform (GWDT) [START_REF] Soille | Generalized geodesy via geodesic time[END_REF] was operated on the raw images to integrate cell shape and light intensity information. For the preparation of GWDT, cell shape information was first extracted. Instead of using the thresholding method with a global fixed threshold value, here the adaptive binarisation method was applied [START_REF] Roth | Adaptive thresholding using the integral image[END_REF]. A raw cell image and corresponding result after adaptive binarisation are shown in Figures 2(A) and (B), respectively. After that, the obtained binary images were processed by morphological close operation to remove the burrs near cell borders. The areas of less than 200 pixels were filtered out. The areas selected by rectangles in Figure 2 The detailed process of GWDT construction is demonstrated as follows. In a grey scale image f, take a path P going from pixel p 0 in a foreground area to pixel p m in a background area as a m tuple (p 0 , p 1 , . . . , p m ). The time t f (P) is the accumulation of the mean intensity for any two neighbouring pixels along the path P, given by:

t f (P ) = m i =1 f ( p i -1 ) + f ( p i ) 2 = f ( p 0 ) 2 + f ( p m ) 2 + m-1 i =1 f ( p i ). (1)
The geodesic time t f (p 0 , p m ) between p 0 and p m is defined as the smallest value of time for all the possible paths from pixel p 0 to pixel p m , given as

t f ( p 0 , p m ) = min{t f (P )|P links p 0 to p m }.
(2)

With the definition of geodesic time, the value g(p 0 ) at pixel p 0 in a foreground area to the background Y is defined as the smallest value of geodesic time linking p 0 to any pixel y in Y, given as

g( p 0 ) = t f ( p 0 , Y) = min y∈Y t f ( p 0 , y) .
(3)

For three clustered cells shown in Figure 3(A), the cell shape information was first extracted through the thresholding method. A mask area containing the three clustered cells was obtained, as shown in Figure 3(B). Within the mask area, the light intensity can also be extracted (Fig. 3C). The GWDT can then be obtained using Eq. ( 3), as shown in Figure 3(D).

Clustered cells detection

Before the operation of clustered cell separation, one needs to determine clustered cells. In raw negative phase contrast images, the light intensity achieves its higher value around nucleus area and gradually decreases towards cell boundaries. Therefore, cell locations can be determined by detecting the maximum value (peaks) [START_REF] Wang | Segmentation of the clustered cells with optimized boundary detection in negative phase contrast images[END_REF]. However, due to the existence of cell organelles, which have higher light refractive index and hence cause higher light intensity in cell images, the pseudo-peaks may be detected if only light intensity is considered, as shown in Figure 3(C). This will lead to oversegmentation.

In this study, instead of using light intensity, the GWDT was applied to detect clustered cells through peak detection. The regional maximum points were detected and taken as peaks. The peaks were then used to represent locations for individual cells. For the clustered cells shown in Figure 3(A), three peaks were detected with GWDT, as shown in Figure 3(D). One can see that the integration of shape and light intensity information helps to reduce oversegmentation.

The detection of clustered cells is demonstrated in Figure 4. For a raw image shown in Figure 4(A), the thresholding method was first applied to get the preliminary segmentation. The GWDT was then constructed in the preliminarily detected mask areas. The mask areas of clustered cells were detected through peak detection, as shown in Figure 4(B).

Masks containing only one peak were taken as individual cell regions, whereas the ones containing multiple peaks were taken as clustered cell regions. In the figure, the yellow and green mask areas correspond to the detected individual and clustered cells, respectively.

Region skeleton-based clustered cell separation

After clustered cell detection, the next step is to perform separation of them. Here a region skeleton-based method was proposed to separate clustered cells, including skeleton extraction, skeleton thinning, graph model representation, splitting points searching, dividing path searching and contour expansion.

A region skeleton is a curved line going through a mask area and normally has the equal distance to the boundary of the mask area. It is used to represent the topological features of clustered cells. In this study, augmented fast marching method (AFMM) [START_REF] Telea | An augmented fast marching method for computing skeletons and centerlines[END_REF] was used to extract skeletons in the mask areas. The method relies on the evolution of the parameterised boundaries of the mask areas. In the method, pixels along a boundary of a mask area are first coded to be parameterised. The parameterised boundary pixels are then evolved inward the mask area under a velocity field constructed by solving an Eikonal equation. When two pixels in the evolving boundary encounter each other and the difference of their codes is larger than a threshold value, the encountered pixels are transformed into a skeleton point. By doing so, the parameterised boundary will be converted into a skeleton for the mask area.

Take the clustered cells shown in Figure 3(A) and the preliminary segmentation shown in Figure 3(B) as examples. For the mask area, the skeleton was first obtained using the AFMM method, as shown in Figure 5(A). The skeletons directly obtained by AFMM are normally not single-pixel width, as shown in the inset enlarged figure . Here a modified thinning algorithm [START_REF] Zhang | A fast parallel algorithm for thinning digital patterns[END_REF][START_REF] Chin | A one-pass thinning algorithm and its parallel implementation[END_REF][START_REF] Holt | An improved parallel thinning algorithm[END_REF] was used to transform the preliminarily obtained skeletons into 8-connected components. In the thinning processing, each foreground pixel is compared with eight different 3 × 3 templates (Fig. 5B). In the templates, P is the centre of the templates, 0 and 1 are the pixels in the background and foreground, respectively, and X means the value can be either 0 or 1. When a pixel in the preliminarily detected skeletons pixel is superimposed with the centre of a template and the eight neighbourhood pixels in the image are also matched with that of the template, the pixel will be deleted. The eight templates scan the foreground pixels in the image until no more foreground pixel can be deleted. The remaining skeleton will be taken as the final one for the mask area. Figure 5(C) shows the obtained skeleton after the thinning processing.

With the obtained skeletons, the graph theory is applied to better represent the relationship among different regions in the mask area. The graph is consisted of vertices and edges. The contiguous vertices are linked by edges. In the construction of graphs, the pixels in the extracted skeleton are divided into three types: end-points, junction-points and connecting-points [START_REF] Bai | Path similarity skeleton graph matching[END_REF]). An end-point (V 1 , V 2 , and V 3 in this case) is a pixel at the end of a skeleton branch with only one pixel from the skeleton in its eight neighbourhoods. A junction-point (V 4 ) is the intersection of different skeleton branches, with at least three pixels in its eight neighbourhoods. All the rest pixels in the skeleton are defined as connecting-points. The end-points and junction-points are called vertices. The skeleton branches connecting any two vertices are defined as edges (e 1 , e 2 and e 3 ).

The constructed graphs can be represented by three elements, numbered vertex set V = {V 1 , V 2 , V 3 }, edge set E = {e 1 ,e 2 ,e 3 } and relations ϕ among vertices and edges, which is given as ϕ:E→E(V) = {(u,v)|u,v ࢠ V}. In this study, a modified graph G is proposed. In the graph, the values of the GWDT along the detected skeleton are taken as the fourth element. The modified model G can then be expressed as

G = (V, E , ϕ, W).
(4)

In the equation, W = {w(e 1 ), w(e 2 ), w(e 3 )}, where w(e i ) represents the distribution of GWDT along a skeleton branch linking two vertices:

w(e 1 ) = {d (V 1 ), d ( p 1 ), d ( p 2 ), . . . , d ( p n ), d (V 4 )}, (5)
where {p 1 ,p 2 , . . . ,p n } are the connecting-points along the skeleton path, and d(p i ) is the value of GWDT.

In order to suppress the noise of the region boundary, the width vector W is smoothed by a local regression and a second order polynomial model. The local minimums of the smoothed width vector w(e i ) correspond to the turning points of the GWDT values along the skeleton e i and are likely to be the contacting points of cells. The distribution of GWDT values along the three edges (e 1 , e 2 and e 3 ) is shown in Figure 5(E). Two local minimums were found along the curves of w(e 1 ) and w(e 2 ), corresponding to two splitting points in the cell region. Because cells have certain size, here the skeleton branches with length smaller than 10 pixels are ignored, like the red curve w(e 3 ) in Figure 5(E).

With the detected splitting points, the dividing paths were searched to split the clustered cells by searching slope-lines [START_REF] Serrat | On ridges and valleys[END_REF]. In hydrology, the water flow runs downhill along the steepest slope, which is called slope-lines. In morphology, slope-lines refer to the paths along which the light intensity gradient achieves its highest values along the path. Assuming that the GWDTs are mountains, slope-lines split the mountains and realizes the optimal separation of contacting cells. In this paper, slope-lines are determined by searching the pixels with minimal intensity level in their eight-neighbourhood. The path searching starts from the splitting points. Two opposite pixels in the eight-neighbourhood which have smallest intensity value will be taken as the next point on the path, leading the directions of the path. Along the two directions, pixels with minimal intensity in eight-neighbourhood of current pixels will be added onto the path until it reaches the background, as shown in Figure 5(G). The separation results for the clustered cells in Figure 3(A) are shown in Figures 5(F) and (H) for 3D and 2D images, respectively.

Compared with the methods which solely rely on the shape or light intensity, the combination of GWDT and cell skeleton reduces oversegmentation and improves cell separation accuracy. Here we take a cell in Figure 6(A) as an example. The preliminarily detected mask area is shown in Figure 6(B). The mesh plot of the distance transform is shown in Figure 6(C). One can see that there is a local minimum of the width distribution along the detected skeleton (Fig. 6D), which induces oversegmentation (Fig. 6E). However, in the proposed method, which combines GWDT and cell skeleton, the oversegmentation can be avoided. The mesh plot of the GWDT of the mask area for the cell is shown in Figure 6(F). Within the mask area, only one peak is detected. Accordingly, the result shows that there is no local minimum for the grey-weighted distance transform value along the skeleton (Fig. 6G). Therefore, the cell can be successfully detected (Fig. 6H).

After the preliminarily segmented mask areas were separated with the proposed method, the contours of the separated subareas are first extracted, as shown in Figure 7(A). To get the optimised cell boundary detection, the contour expansion method was applied to the contours, as what we did previously [START_REF] Wang | Segmentation of the clustered cells with optimized boundary detection in negative phase contrast images[END_REF]. In the contour expansion method, the contours are first parameterised. An energy function involving light intensity and geometries of the contours are constructed. The final contour is obtained by iteratively minimising the total energy of the evolving contour. Driven by the field of light intensity, the contour converges towards and will stop at the actual cell boundary. The optimised cell boundary detection is shown in Figure 7(B).

Comparison with other methods

To validate the proposed method, a comparison of clustered cell segmentation was implemented with several other popular methods, including the watershed method [START_REF] Meyer | Topographic distance and watershed lines[END_REF], an intensity-based method [START_REF] Wang | Segmentation of the clustered cells with optimized boundary detection in negative phase contrast images[END_REF], a markercontrolled watershed method [START_REF] Song | Accurate segmentation of cervical cytoplasm and nuclei based on multiscale convolutional network and graph partitioning[END_REF], and the proposed method with 393 clustered cells. The result is shown in Table 1.

Among the four methods, the proposed method achieves 85.8% accuracy rate, followed by the intensity-based method, which gives 80.2% accuracy rate. The other two methods provide a lower accuracy rate of around 70%. Authors want to note here that all cells involved in the evaluation in this study (Table 1 and whole cell evaluation in Fig. 9) were manually confirmed one by one. Because all the sample images were captured using time lapse phase contrast microscope with a certain time interval, all cells can be manually distinguished with the assistance of their adjacent frames. This guarantees the validity of the ground truth selected for evaluation. [START_REF] Meyer | Topographic distance and watershed lines[END_REF], marker-controlled watershed method [START_REF] Song | Accurate segmentation of cervical cytoplasm and nuclei based on multiscale convolutional network and graph partitioning[END_REF], intensity-based method [START_REF] Wang | Segmentation of the clustered cells with optimized boundary detection in negative phase contrast images[END_REF] and the proposed method, respectively.

Compared with the other three methods, especially the watershed method and the shape-based method, the proposed method has much lower oversegmentation rate. Figure 8 shows the comparison of some typical examples of clustered cell separation. The first column shows the raw images of typical clustered cells, including mitotic cells, elongated cells and clustered cells.

Among the four methods, the watershed method [START_REF] Meyer | Topographic distance and watershed lines[END_REF] and the marker-controlled watershed method [START_REF] Song | Accurate segmentation of cervical cytoplasm and nuclei based on multiscale convolutional network and graph partitioning[END_REF] can easily induce oversegmentation. Both methods are shape-based methods. The light intensity information within the mask areas is ignored. This will cause oversegmentation, especially when it comes to the mitotic cells which normally have much higher image contrast. The light intensity-based method [START_REF] Wang | Segmentation of the clustered cells with optimized boundary detection in negative phase contrast images[END_REF] took only the light intensity information into account. When there are multiple pseudo-peaks, especially to elongated cells, the method will cause oversegmentation. The proposed method integrates both light intensity information and the shape information and achieves the improved segmentation results, as shown in the fifth column in Figure 8.

Whole image evaluation

The proposed method was implemented on the whole image of 40 negative phase contrast images. Undersegmentation, over-segmentation, corresponding positive predictive value (PPV) and recall were shown in Figure 9. The overall accuracy rate, PPV and recall of cell segmentation using the proposed method are 97.16%, 98.82% and 98.64%, respectively.

The image segmentation of four negative phase contrast images is demonstrated in Figure 10. For four raw phase contrast images shown in Figure 10(A), corresponding segmentation results are given in Figure 10(B). Under segmentation and over segmentation are marked by yellow and green arrows, respectively.

The accuracy rate of the algorithm is limited by the initial segmentation of binarisation and the skeleton-based separation method. The performance of adaptive binarisation in the initial segmentation needs to be improved, especially when cell density is high. A possible way to further improve the accuracy rate of the clustered cells could be the application of the deep learning or machine learning methods.

Conclusion

In summary, we have developed a new method for clustered cell segmentation in negative phase contrast images of MCF 10A human breast cells. This is done by fusing cell shape and light intensity information through the construction of greyweighted distance transform (GWDT). The method consists of two steps. In the first step, the clustered cells are detected by peak detection in GWDT. This avoids the pseudo-peaks induced by intercellular organelles and reversed image contrast for cells with large height. The second step is the separation of clustered cells. To do this, a modified region skeleton method was proposed. The value of GWDT along the detect skeletons was used to separate clustered cells. Compared to the widely applied shaped-based or light intensity-based methods, the proposed method provides more accurate detection and hence better separation of clustered cells. The achieved accuracy rate is 85.8% and 97.16% for clustered cells and all cells, respectively, in negative phase contrast images.

  (A) are two examples of the clustered cells.

Fig. 1 .

 1 Fig. 1. Schematic diagram of the proposed segmentation scheme. (A) A simulated image showing three contacting cells in an area. (B) The binary image of (A), indicating the preliminary shape of the clustered cell area. (C) Mesh plot of light intensity in the area, where a spurious peak is detected. (D) Mesh plot of the grey-weighted distance transform by integrating shape and light intensity information, where the spurious peak is eliminated. (E) The skeleton of the clustered cell region and corresponding cell separation result. (F) Refined contours are obtained with an active contour method.

Fig. 2 .

 2 Fig. 2. Adaptive binarisation of cell images. (A) Raw negative phase contrast cell image. (B) Binary image with the adaptive threshold method. Cell shape information is obtained by the preliminary segmentation.

Fig. 3 .

 3 Fig. 3. Demonstration of the construction of grey-weighted distance transform. (A) A raw image of clustered cells. (B) Binary image of the clustered cells with a preliminarily detected mask area used to represent cell shape information. (C) Light intensity within the preliminarily detected mask area for the three clustered cells. Within the areas, two pseudo-peaks were detected, which normally cause oversegmentation. (D) Grey-weighted distance transform obtained by integrating cell shape and light intensity information. With the method, peaks can be correctly detected for each cell.

Fig. 4 .

 4 Fig. 4. Detection of the clustered cells. (A) A raw negative phase contrast image. (B) Result of clustered cell detection in the grey-weighted distance map. Green mask areas refer to suspected clustered cells whereas yellow areas refer to individual cells.

  Figure 5(D) shows the diagrammatic form of the graph corresponding to the skeleton shown in Figure 5(C).

Fig. 5 .

 5 Fig. 5. Demonstration of the region skeleton-based clustered cell separation. (A) A skeleton of the mask area in Figure 3b is obtained (orange curve) with AFMM method. The inset in the bottom left corner is the enlarged area of the green dashed box. (B) Templates for the skeleton thinning operation. In the templates, P is the centre of the templates, 0 and 1 are the pixels in the background and foreground, respectively, and X means the pixel can be either 0 or 1. (C) The skeleton after thinning operation. V 1 , V 2 and V 3 are end-points of the skeleton, and V 4 is the node-point. The inset in the bottom left corner is the enlarged area of the green dashed box. (D) Graph diagram of (C), where vertices are linked by edges of e 1 , e 2 and e 3 . (E) Width distribution curves of lines e 1 , e 2 and e 3 . Two local minimums are detected on curves w(e 1 ) and w(e 2 ), whereas w(e 3 ) is ignored due to its short length. (F) Dividing paths are searched on the smoothed grey-weighted distance map. (G) Enlarged area selected by the dashed box in (F). The star mark represents the start pixel corresponding to the local minimum of w(e 1 ) in (E), and the points of the path are searched one by one, as shown by arrows. The right half of the curve represents the path behind the surface. (H) 2D plot of smoothed grey-weighted distance map and dividing paths. Clustered cells are separated by the skeleton-based method with the fused information of intensity and shape.

Fig. 6 .

 6 Fig. 6. Comparison of cell segmentation with grey-weighted distance transform and Euclidean distance transform. (A) Raw negative phase contrast image of an elongated cell. (B) Binary image of the elongated cell with only shape information. (C), (F) Mesh plot of Euclidean distance transform (C) and grey-weighted distance transform (F). (D), (G) Distribution of values along the skeleton in Euclidean distance transform (D) and grey-weighted distance transform (G), which corresponds to the red curves in (C) and (F). A local minimum was detected in (D). (E), (H) Segmentation results with distance transform and grey-weighted distance transform, respectively. The distance transform causes oversegmentation, whereas the grey-weighted distance transform gives a correct segmentation.

Fig. 7 .

 7 Fig. 7. Optimisation of cell boundary detection. (A) Initial contours extracted from preliminarily segmented mask areas. Yellow contours refer to the initial contours of individual cells, while the green ones refer to the initial contours of clustered cells after the skeleton-based separation method. (B) Result of the optimised contour detection after contour expansion operation.

Fig. 8 .

 8 Fig. 8. Comparison of cell segmentation with different methods. (A) Raw negative phase contrast images of candidate clustered cells. (B)-(E) Segmentation results by watershed method[START_REF] Meyer | Topographic distance and watershed lines[END_REF], marker-controlled watershed method[START_REF] Song | Accurate segmentation of cervical cytoplasm and nuclei based on multiscale convolutional network and graph partitioning[END_REF], intensity-based method[START_REF] Wang | Segmentation of the clustered cells with optimized boundary detection in negative phase contrast images[END_REF] and the proposed method, respectively.

Fig. 9 .

 9 Fig. 9. Whole cell evaluation of the proposed method. The stacked bars show cell number of correct segmentation, under segmentation and over segmentation, scaled by the left vertical axis. The line plots show the PPV and recall, scaled by the right vertical axis.

Fig. 10 .

 10 Fig. 10. Cells segmentation results on whole cell images. (A) Raw negative phase contrast images of MCF 10A cells. (B) Segmentation results with the proposed method. Under segmentation and over segmentation are pointed by yellow and green arrows, respectively.

Table 1 .

 1 Comparison of clustered cell separation with four different methods.

		Number of	Correct seg-
	Methods	cells	mentation	Oversegmentation	Undersegmentation	Accuracy rate	PPV	Recall
	Watershed method (Meyer,		279	93(23.7%)	20(5.1%)	71.0%	75.00%	93.31%
	1994)
	Intensity-based method		315	37(9.4%)	40(10.2%)	80.2%	89.49%	88.73%
	(Wang et al., 2015)	393
	Marker-controlled		269	97(24.7%)	29(7.4%)	68.5%	73.50%	90.27%
	watershed method (Song
	et al., 2015)
	The proposed method		337	24(6.1%)	31(7.9%)	85.8%	93.35%	91.58%

Acknowledgements

This study is financially supported by National Natural Science Foundation of China (Grant No. 51775028) and Beijing Natural Science Foundation (Grant No. L150002).