Alireza Asvadi

Andrei Mitriakov

Christophe Lohr

Panagiotis Papadakis

Digital Twin Driven Smart Home: A Feasibility Study

Keywords: Ambient assisted living, cyber physical system, living-lab, semantics

We aim to facilitate the daily-life activities of frail or elderly people in collaboration with mobile assistive robots through the means of a digital twin-powered smart home. Being able to quickly and efficiently produce a digital twin of the human user's environment, can help to further develop personalized assistive solutions. As our first investigation toward this goal, we describe our proof-of-concept "digital twindriven smart home" implementation. It consists of a virtual representation, robot navigation and environment semantics using open-source software. The initial obtained results on the building process of the digital twin are encouraging and suggest the possibility of integration of digital twin for smart spaces.

Introduction

In recent years, the availability of low-cost sensors and open-source middle-ware software have opened up interesting research areas in the robotics field. In particular, the next generation of simulation models called Digital Twin (DT), which represents a continuous virtual replica of a physical system, has gained increasing attention. It can be used to create a simulation of a smart home including assistive/service robots and human users. It has applications in the optimization of robots and smart home settings. For example, finding the optimal number and configuration of sensors especially when new robots or users are introduced. As another examples, monitoring in real-time, further analysis and learning of edge cases and rare situations in DT for the safety of human users, e.g., by pushing those events from simulation to real-world and vice-versa (thanks to DT that enables a bi-directional link between simulation and real-world). In practice, such optimization should be carried out over a variety of houses and users.

In this context, this paper aims to put forward the use of DT technology in the development of a smart home system. The main contributions are:

. This work looks at the practical aspects of the building process of DT. . Provides analysis and insights into the difference in the performance between the actual and simulated environments.

The paper is organized as follows. Section 2 gives a definition and recent related works. Section 3 details the proposed DT development process. Experimental results are described in Section 4, and Section 5 brings some concluding remarks and proposes directions for further research.

Background

With the progress of the internet of things and artificial intelligence, several strategies were developed to advance smart technologies, including Industry 4.0 [START_REF] Lasi | Industry 4.0[END_REF] and Society 5.0 [START_REF] Deguchi | From Smart City to Society 5.0[END_REF]. At the core, these technologies are based on a cyberphysical system (CPS) that integrates physical and digital entities. In this context, DT is an approach to implement cyber-physical integration. The DT is a digital representation of a physical system that precisely models and can predict how the actual system will perform. In this respect, it serves as a tool to bridge the gap between physical and digital worlds. Definition. Some of the commonly accepted concepts in the definition of the DT are the following : it is a digital representation of a physical entity (including geometric, functional, and usage models), being ultra-realistic and fully consistent with the physical entity, continuous and automatic update to be able to adapt to changes, semantic content to support intelligent decision making, and implementation in a decentralized structure to collaborate with other DTs [START_REF] Tao | Digital twin driven smart manufacturing[END_REF]. DT for Robotics. While there is a rich literature on DT concept [START_REF] Fuller | Digital twin: Enabling technologies, challenges and open research[END_REF][START_REF] Jones | Characterising the Digital Twin: A systematic literature review[END_REF], there have been few reports on the 3D modeling aspect in robotics and smart home applications. Kousi et al. [START_REF] Kousi | Digital twin for adaptation of robots' behavior in flexible robotic assembly lines[END_REF] combined multi-sensor data and CAD models to create a virtual representation (including semantic and dynamic update components) of a production environment. Havard et al. [START_REF] Havard | Digital twin and virtual reality: a co-simulation environment for design and assessment of industrial workstations[END_REF] proposed a co-simulation approach between DT and virtual reality in a human-robot collaborative workplace. Another notable work is Habitat 2.0 [START_REF] Szot | Habitat 2.0: Training home assistants to rearrange their habitat[END_REF] which is a platform for training simulated robots in interactive environments. Phanden et al. [START_REF] Kumar Phanden | A review on simulation in digital twin for aerospace, manufacturing and robotics[END_REF] recently reviewed existing simulation software for DT development. For robotics, they found Gazebo and V-REP [START_REF] Rohmer | V-REP: A versatile and scalable robot simulation framework[END_REF] for modeling sensors and image-based systems, ROBOTRAN to model multi-body systems, and ROBCAD for multi-device robotic and automated manufacturing processes. DT for Smart Home. As to digital models of buildings, To et al. [START_REF] To | Drone-based AI and 3D Reconstruction for Digital Twin Augmentation[END_REF] performed drone-based reconstruction for DT augmentation of buildings. Calderita et al. [START_REF] Luis V Calderita | Designing a cyber-physical system for ambient assisted living: a use-case analysis for social robot navigation in caregiving centers[END_REF] proposed a CPS for ambient assisted living including humans and assistive robots in a care-giving center. A related term to DT is Building Information Modeling (BIM). Czerniawski and Leite [START_REF] Czerniawski | Automated digital modeling of existing buildings: A review of visual object recognition methods[END_REF] reviewed approaches for BIM. They identified the following steps for model creation: collecting sensor data and 3D reconstruction, semantic and geometric modeling for the recognition of semantic classes in sensor data and describing the instances' shapes, and finally the creation of BIM that integrates semantic and geometric components. In comparison with BIM (which is one-way modeling), in a DT there is a bi-directional data flow between the virtual and physical entities. There are a number of other works and projects that discuss related concepts on ambient assisted living and the synergy of ambient intelligence with robots. Examples are GIRAFF+ that consists of a network of smart home sensors and a telepresence robot to monitor elderly people, and STRANDS that explores adaptation to changes over time for mobile robots in 3D dynamic human environments. Previous work of the authors concentrated on the potential services that could be offered from such systems [START_REF] Papadakis | System design for coordinated multi-robot assistance deployment in smart spaces[END_REF].

DT-Powered Smart Home

The DT of a smart home is defined as: using a priori information (e.g., architectural plans, input from inhabitants) and given input sensory data (e.g., home automation and robot sensors), build a continuous model including robots and user models. As a first step toward this goal, a prototype DT model of the smart home has been developed. Figure 1 demonstrates an abstraction of our understanding of what a DT workflow should consist of. Ideally, such a workflow has limited dependence on human supervision so as to reduce errors and favor portability and generalization. The steps are defined as follows. The smart home is scanned using a hand-held RGB-D camera and a SLAM algorithm to build a reconstruction of it. A modeling software is used to construct the CAD models of the home and objects inside it from the scan. The CAD models and a robot model are imported into a robotic simulator. A proof-of-concept is developed to enable semantic understanding of objects.

The DT can be conceptualized in data, algorithmic and model. The data layer represents the communicated data including point cloud, navigation (velocity and pose), etc. The algorithmic layer includes software to develop the modules. The model layer comprises the geometric, functional and semantics including the CAD models, a TurtleBot1 , and a pre-trained object detection model. More information can be found on the project website2 .

3D Reconstruction

To have a metric measurement of the actual environment a hand-held scanning is performed. The Orbbec Astra Pro RGB-D camera is used which is shown to be one of the best low-cost sensors for the reconstruction of indoor spaces [START_REF] Da | Comparison of RGB-D sensors for 3D reconstruction[END_REF]. The RTAB-Map (Real-Time Appearance-Based Mapping) [START_REF] Labbé | RTAB-Map as an open-source lidar and visual simultaneous localization and mapping library for largescale and long-term online operation[END_REF] algorithm is used. It consists of loop closure and proximity detection (using a bag of words), graph optimization (to decrease odometry drift) and global map assembling (using local occupancy grid). It requires an odometry estimation, i.e., motion estimation between consecutive scans. This is performed by computing visual odometry using RTAB-Map and the camera. The process outputs a dense point cloud. The CAD models are created from point cloud data. A modeling software is used to craft the CAD model of the home structure (e.g., walls) from the point cloud, also online public repositories are used to obtain the models of the real home items. It should be noted that an automatic transformation of 3D scans into CAD is an open question (see e.g., Scan2CAD project) that will not be addressed in our work.

Physics-based Simulation

The Gazebo is used to simulate the dynamics of the system. A digital 3D model is developed to represent the smart home environment and the objects in the simulator. The CAD models are used to define the visual and collision properties of the items in this model. Other physical properties can be defined, for example pose, mass, being static (e.g., walls) and dynamic (e.g., objects). A TurtleBot model is spawned in the simulated system. The use of Gazebo with ROS enables communication between the simulated and the real system.

Semantics

The detection of the categories of objects of interest (e.g., tables and chairs) and their position can have several benefits in a DT framework. The You Only Look Once (YOLO)-v3 algorithm [START_REF] Redmon | Yolov3: An incremental improvement[END_REF] is used for object detection. It reads camera image data and returns the 2D bounding boxes. The pre-trained model is able to detect 80 classes including chair, sofa, table, bottle, microwave, person, etc. The darknet ros 3d package3 is used to add the 3D bounding boxes of objects to the YOLO. It combines the detected 2D bounding boxes with point cloud data to calculate the 3D bounding boxes.

Results and Discussion

Sections 4.1 and 4.2 describe the process results for the creation of the digital 3D model and a comparison of the DT and the actual system. Section 4.3 provides discussion on the semantics aspect. The approach is evaluated in a small experimental flat including a kitchen, a bathroom and a living room with furniture (see Fig. 2-TopLeft).

Creation of the Digital 3D Model

3D Reconstruction A detailed colored 3D point cloud of the flat is achieved which has about 700k points (see Fig. 2-Bottom). The following are some of the observed limitations.

. The windows need to be covered in sunlight. The Astra Pro camera works based on infrared technologies. The light coming from the infrared projector can be outshined by the sunlight and no point cloud can be computed. . The mirror, window glass, TV surface, and black colored objects degrade the quality of the depth image and the point cloud.

. The textureless walls, ground and ceiling must be avoided or rectified (e.g., by adding picture frames to plain walls). Because the visual odometry in RTAB-Map is based on extracted features from the RGB-D images. Therefore, in environments without enough features, the odometry cannot be computed.

Digital 3D Model An idealistic 3D model of the flat was made by a team of students. The model was built by first measuring different dimensions of the flat and then recreating it using a modeling software. A quantitative evaluation is performed with a comparison of our CAD model against this idealistic model (see Fig. 3-Top). To perform the evaluation the models are converted to point cloud data using ray tracing operations. The root mean squared error (RMSE) and Hausdorff distance metrics (which are common general-purpose metrics for comparing 3D shapes) are used to estimate the differences (see Table 1). The accuracy of about 9 cm average RMSE is obtained which is encouraging and demonstrates an interesting perspective for an automated reconstruction-based model creation. As it can be seen in Fig. 3-Top, the perpendicularity of the walls in our developed model (shown in red) are preserved well with only a slight slope on the right side. Figure 3-Bottom shows the furnished smart home with a simulated TurtleBot in the Gazebo.

Robot Navigation

Experiments were performed to compare the traversed trajectory and the generated 2D maps in the actual and simulated conditions.

Trajectory Evaluation A velocity message with the linear and angular values of +0.25 and -0.1 is published (with the frequency of 10 Hz for 4 seconds) to control the translational and rotational speed of the robot. The objective is to perform the navigation in the DT and reality, and to compare the trajectories.

To obtain the trajectory of the actual Turtlebot, a visual marker was used. We position the marker on the upper part of the Turtlebot in the field of view of a high-resolution camera (see Fig. 4). The robot poses can be retrieved by detecting the marker in the camera image. A toolkit tracks marker poses which enables tracing of the Turtlebot trajectory in the camera coordinate frame. The experiments were performed in 5 iterations. The average and standard deviation of the trajectories in both real and simulated cases were reported in Fig. 5. In addition a simplified analytical solution for the expected trajectory of a robot performing uniform circular motion is computed. It is equal to traverse of 0.4 radians on the perimeter of a circle with radius of 2.5 m. It should be mentioned that the complete Turtlebot model is a differential drive consisting of two wheels and two caster wheels which was not considered here. It was noticed that the turtlebot in the real condition traversed a shorter distance in the trajectory (about 25 cm difference with the simulated one). This was probably because of the uncertainties on the friction parameters between the robot and the ground surface. Due to using a simplified model, there were also differences in the recorded trajectories in comparison with the analytical solution.

Understanding the difference between the virtual model and the real world and adjusting such parameters can improve the robot's behavior in reality. It helps to accurately estimate and plan the robot trajectory and in the application level to provide both reliable and consistent service for elderly users in the smart home.

Mapping Experiment A 2D mapping process is performed and compared in the DT and reality. The turtlebot rotates 360 degrees (with the angular speed of 10 degrees per second) to cover the whole space. The gmapping4 package from ROS is used [START_REF] Grisetti | Improved techniques for grid mapping with rao-blackwellized particle filters[END_REF]. It is a laser-based SLAM. To be able to work with a RGB-D camera, it converts a depth image to a fake laser scan. It uses particle filters (each particle carries an individual map) and several adaptive techniques to learn the occupancy grid map of the environment. Figure 6 shows examples of the generated maps. It was observed that the actual turtlebot could not cover the full profile of the home using the rotation message of 360 degrees. The rotation of 540 degrees was published to the actual turtlebot to fulfill the coverage of the whole area.

Semantics

Object detection can be used as a coupling mechanism of the real and the digital worlds. As a feasibility study, 2D and 3D versions of YOLO object detection were employed to detect the objects in the smart home. Quantitative Analysis of 2D Object Detection Figure 7 shows 2D object detection results in different areas of the smart home. An evaluation is performed to quantitatively measure the object detection performance in the smart home. A number of 30 detection frames were selected randomly and the precision, recall, and F1 scores computed for the two most frequent items seen in the smart home: the chairs and TV monitors, with the total number of 50 chairs and 41 monitors observed in all images. As it can be seen in the Table 2, the detection rates were different for each item. For example, the TV monitors were detected with higher precision whereas the detector has a better recall value on chairs. Qualitative Analysis of 3D Object Detection Figure 8 shows an example of the 3D detection results. The top shows the 2D detections. The yellow boxes in the bottom show projected detections into the point cloud in the 3D space.

Conclusions

In this paper, we presented an approach for the generation of a digital 3D model of a smart home using open-source software and a low-cost sensor. The methodology and initial results (including the obtained models) are promising and appear to have a potential for the design of DT-driven smart homes. However, further research is needed to explore the integration of the digital profiles of the connected home automation sensors (e.g., thermostat and pressure sensors) into the DT. The next challenge is to extend the developed DT to include other entities: sophisticated robots (than a simple TurtleBot) and human users (e.g., to sense where the user is and what is doing). Another interesting directions are benchmarking, further analysing and quantifying the different steps of DT creation.

Fig. 1 .

 1 Fig.1. The conceptual digital twin of the smart home. The dashed-red arrows mark the process that are not implemented in a fully automated way. The white box shows the component that is developed primarily for a feasibility study.

Fig. 2 .

 2 Fig. 2. Top-left: Top panoramic view of the smart home. top-right: Illustration of an area in the reconstructed map (from the view shown by a red arrow in the left image). Bottom: the full reconstructed map (from top view).

Fig. 3 .

 3 Fig. 3. Top: Illustration of the representations. Superimposed simplified mesh (obtained from the point cloud), the model created from point cloud (red), and the idealistic model (green) for a visual inspection. Bottom: The simulated smart home with a Turtlebot model.

Fig. 4 .

 4 Fig. 4. Left: The trajectory evaluation experiment. Right: Left and right show the Turtlebot at the start and end of the experiment. Top and bottom show the simulated and actual setup.

Fig. 5 .

 5 Fig. 5. The average and standard deviation of the obtained trajectories in the real, simulated and analytical situations (the starting point is on the left-side of the curves).

Fig. 6 .Fig. 7 .

 67 Fig.6. Illustration of the 2D occupancy grid maps. From left to right: the generated grid map in the simulated environment, the created map in the actual home, the obtained map after the additional rotation in the real case.

•Fig. 8 .

 8 Fig. 8. 3D Object detection using YOLO-v3 and Darknet ROS 3D. Left: shows a person standing and moving in the smart home. Right: shows results in another setting with different furniture.

Table 1 .

 1 A comparison of the quality of the developed smart home's CAD model against the idealistic model.

	Metrics	RMSE Hausdorff
	Error (in m) 0.088	1.354

Table 2 .

 2 A comparison of 2D object detection in the smart home.

	Object/Metrics Precision Recall F-Score
	Chair	0.53	0.88	0.66
	Monitor	0.97	0.68	0.80

https://www.turtlebot.com/turtlebot2/

https://sites.google.com/view/heron-project

https://github.com/IntelligentRoboticsLabs/gb visual detection 3d

http://wiki.ros.org/gmapping

Acknowledgments. The authors would like to thank Maroua El Houicha, Zhi Zhang, Jeanne Thévenin, and Damien Bouchabou, students of IMT Atlantique -Brest, for providing the idealistic 3D smart home model. The authors would also like to thank Fabien Dagnat and Jean-Christophe Bach for helpful discussions on the digital twin architecture. The work has been performed in the context of project HERON -Habitat intelligent Et RObotique d'assistance basés sur le jumeau Numérique -financed by the region of Brittany (Région Bretagne) and the department of Finistere (CD29).


Emilie Yu
email: emilie.yu@inria.fr

J Andreas Baerentzen

Rahul Arora

Karan Singh

Adrien Bousseau
email: adrien.bousseau@inria.fr.

Piecewise-Smooth Surface Fitting onto Unstructured 3D Sketches

Keywords: Computing methodologies → Shape modeling sketch-based modeling, surfacing, VR sketching

Fig. 1. 3D sketches created in Virtual Reality (VR) or with sketch-based modeling systems often depict piecewise-smooth surfaces (a), but lack proper inter-stroke connectivity to detect the individual surface patches, as illustrated in the insets where pen strokes do not intersect precisely, and where detail strokes lie on the imaginary surface without being connected to other strokes. State-of-the-art surfacing algorithms only produce smooth surfaces from such sparse and unstructured 3D data (b). Our algorithm segments such an initial smooth surface into regions aligned with the pen strokes to produce a piecewise-smooth surface that better captures the intended shape. ©James Robbins, used with permission.

We propose a method to transform unstructured 3D sketches into piecewise smooth surfaces that preserve sketched geometric features. Immersive 3D drawing and sketch-based 3D modeling applications increasingly produce imperfect and unstructured collections of 3D strokes as design output. These 3D sketches are readily perceived as piecewise smooth surfaces by viewers, but are poorly handled by existing 3D surface techniques tailored to wellconnected curve networks or sparse point sets. Our algorithm is aligned with human tendency to imagine the strokes as a small set of simple smooth surfaces joined along stroke boundaries. Starting with an initial proxy surface, we iteratively segment the surface into smooth patches joined sharply along some strokes, and optimize these patches to fit surrounding strokes. Our evaluation is fourfold: we demonstrate the impact of various algorithmic parameters, we evaluate our method on synthetic sketches with known ground truth surfaces, we compare to prior art, and we show compelling results on more than 50 designs from a diverse set of 3D sketch sources.

INTRODUCTION

Drawing, be it on cave walls, paper, or touch screens, has been a quintessential tool of visual communication for millennia. Beyond the established importance of drawing on surfaces, sketching midair in AR/VR environments, like waving a magic wand, provides artists the unprecedented freedom to draw directly in immersive 3D [Adobe 2020;Google 2016;[START_REF] Gravitysketch | Gravity Sketch[END_REF]Smoothstep 2021]. Such unstructured 3D curve sketches, while beautiful to behold and effectively perceived by viewers, are non-trivial to transform into 3D surface models for further design communication and processing. We present, to our knowledge, the first approach to transform such unstructured sketches into piece-wise smooth 3D surface models that preserve salient sketched features (Fig. 1).

Our method is motivated by several key characteristics of 3D sketches, highlighted as insets in Fig. 1(a).

-Artists often depict piecewise-smooth surfaces by strategically placing strokes at sharp surface discontinuities, consistent with traditional sketching practices of representing an object through its feature curves -ridges and valleys [START_REF] Cole | Where do people draw lines?[END_REF][START_REF] Gryaditskaya | OpenSketch: a richly-annotated dataset of product design sketches[END_REF]]. -While the main patches that form the envisioned surface are often delimited by strokes, other strokes lie within individual patches to depict sub-parts and decorative details. -3D sketches are often imprecise, exhibiting over-sketching or gaps and missing intersections between strokes.

While the relative importance of these characteristics in 3D sketches may vary by artist and the drawing tool used, they underlie both AR/VR sketching and traditional sketch-based 3D modeling systems [Bae et al. 2008;[START_REF] Kim | SketchingWithHands: 3D sketching handheld products with first-person hand posture[END_REF][START_REF] Xu | True2Form: 3D curve networks from 2D sketches via selective regularization[END_REF]], making our approach well-suited in general to surfacing any sketches, that are provided as an unstructured collection of 3D curves.

Existing surfacing methods only partly account for these characteristics. Point cloud surfacing algorithms target densely sampled surfaces and fail on sparse stroke clouds [START_REF] Hoppe | Surface reconstruction from unorganized points[END_REF][START_REF] Michael Kazhdan | Poisson surface reconstruction[END_REF]] or cannot capture sharp features of the sketch (Fig. 1b).

Methods dedicated to 3D drawings focus on well-connected curve networks [START_REF] Abbasinejad | Surface patches from unorganized space curves[END_REF][START_REF] Bessmeltsev | Design-driven quadrangulation of closed 3d curves[END_REF][START_REF] Orbay | Sketch-based surface design using malleable curve networks[END_REF][START_REF] Sadri | Flow-complex-based shape reconstruction from 3d curves[END_REF][START_REF] Zhuang | A general and efficient method for finding cycles in 3D curve networks[END_REF], where surface patches are bounded by closed cycles of curve segments. Due to their strong requirements on curve network topology, such methods cannot process the imprecise, unstructured stroke clouds we target.

To reconstruct piecewise-smooth surfaces from inaccurate 3D sketches, we must jointly determine which parts of the sketch correspond to different surface patches, and recover the geometry of these patches away from the strokes.

A first challenge is to determine where surface patches should lie in the empty space between the strokes. This task is particularly ambiguous for sparse sketches, since the inter-stroke Euclidean proximity does not necessarily denote geodesic proximity on the envisioned surface. We reduce this ambiguity by complementing the input sketch with a low-fidelity proxy surface, which we obtain by applying the smooth point cloud surfacing method of [START_REF] Huang | Variational implicit point set surfaces[END_REF], or in cases where this automatic method fails, by asking users to create a simple proxy with low-poly modeling tools (see Section 5). The proxy surface provides us with a manifold domain on which we project nearby strokes, such that our problem of locating the surface patches amounts to segmenting the proxy into regions roughly bounded by strokes.

Our second challenge is to shape each proxy region into a surface patch that represents well the geometry of the strokes that project in that region, which is especially difficult when the strokes themselves are sparse and approximate. We address this challenge by representing each patch as an implicit surface defined as the zero level-set of a low degree polynomial -which we refer to as an implicit polynomial surface. Such implicit surfaces offer multiple benefits in our context. They have infinite support, allowing the surface regions they capture to grow or shrink arbitrarily as the algorithm progresses, and they are fast to evaluate and fit to stroke points.

Equipped with the concepts of a proxy surface and implicit polynomial surface patches, we cast finding the piecewise-smooth surface that best represents the sketch as a multi-model fitting problem [START_REF] Isack | Energy-based geometric multi-model fitting[END_REF]. In a nutshell, our algorithm alternates between refining a segmentation of the proxy into regions well represented by a set of implicit surfaces, and improving the implicit surfaces by fitting them over the strokes within each region. After convergence, we generate the final 3D triangle mesh by projecting the proxy surface onto the zero-level sets of the implicit polynomial surfaces, effectively recovering sharp features at the intersections between patches.

While our method automatically produces piecewise-smooth surfaces aligned with the input strokes, it also supports user indications to further improve the quality of the surface; for instance, to force the creation or removal of surface discontinuities. We also rely on user indications to discard strokes that are not supported by our approach, notably internal strokes that lie inside the outermost depicted surface, small disconnected parts, and so-called skeletal strokes that depict thin cylindrical features (see Fig. 15).

In summary, we present the first surfacing method for 3D sketches that produces piecewise-smooth meshes while being oblivious to stroke connectivity. We evaluate our method by surfacing more than 50 sketches created with a variety of VR and sketch-based modeling systems, by comparing to surfacing methods tailored to well-connected curve networks, by measuring deviation from ground truth surfaces from which synthetic sketches were generated, and by studying the impact of the main components of our algorithm. Finally, we provide an interactive visualization of all our results as supplemental materials for close inspection.

RELATED WORK

Surfacing curve networks. The problem of surfacing 3D sketches emerged with the advent of practical user interfaces, algorithms and devices to create such sketches. Early 2D interfaces deduce the depth of the strokes from the intersections they form with other strokes [START_REF] Schmidt | Analytic drawing of 3D scaffolds[END_REF][START_REF] Xu | True2Form: 3D curve networks from 2D sketches via selective regularization[END_REF] or with sketching planes [Bae et al. 2008], effectively producing well-connected curve networks by construction. As a consequence, a number of surfacing algorithms strongly rely on the connectivity of the curve network to identify closed cycles delimiting surface patches [START_REF] Abbasinejad | Surface patches from unorganized space curves[END_REF][START_REF] Orbay | Sketch-based surface design using malleable curve networks[END_REF][START_REF] Sadri | Flow-complex-based shape reconstruction from 3d curves[END_REF][START_REF] Zhuang | A general and efficient method for finding cycles in 3D curve networks[END_REF]]. Each such patch can then be surfaced by propagating geometric information from the boundary curves [START_REF] Bessmeltsev | Design-driven quadrangulation of closed 3d curves[END_REF][START_REF] Hao Pan | Flow aligned surfacing of curve networks[END_REF][START_REF] Stanko | Smooth interpolation of curve networks with surface normals[END_REF]. Curve connectivity information provides strong geometric hints, since surface normals can be estimated at each intersection -as done by [START_REF] Hao Pan | Flow aligned surfacing of curve networks[END_REF] to detect sharp features and determine which curves are flow lines. However, the reliance on accurate curve network topology prevents these methods to work on raw, unstructured stroke clouds as typically produced via freehand VR sketching [Google 2016;[START_REF] Gravitysketch | Gravity Sketch[END_REF]Smoothstep 2021].

Recent studies have shown that precise sketching in VR is more challenging than in 2D due to the lack of a supporting surface for the hand and the need for finer motor control to position strokes in 3D [START_REF] Arora | ILoveSketch: as-naturalas-possible sketching system for creating 3d curve models[END_REF][START_REF] Donaji | The Effect of Spatial Ability on Immersive 3D Drawing[END_REF]. Automatic snapping can help We take as input a 3D sketch (a) and a proxy mesh (b), for example obtained with VIPSS [START_REF] Huang | Variational implicit point set surfaces[END_REF]. We iteratively improve a segmentation of the proxy mesh (c, left) and parameters of a set of surface models (c, right) to obtain a decomposition of the surface into smooth patches that approximate the input strokes well. We represent each surface model as the zero level set of an implicit polynomial, which offers infinite spatial support and a good balance between expressivity and complexity. Finally, we project the vertices of the proxy mesh onto the surface models to recover the final piecewise-smooth surface with sharp features (d). When appropriate, we trim the proxy mesh to match user-annotated boundary strokes (a, blue), and we leverage sketch symmetry by surfacing only half of the sketch and mirroring the result.

correct for inaccuracy [START_REF] Machuca | Multiplanes: Assisted freehand VR Sketching[END_REF]Yu et al. 2021a,b], but such assistance can be detrimental to user creativity, as mentioned by users of the recent CASSIE system [Yu et al. 2021a] who criticized being forced to think of their designs in terms of curve networks. Processing 3D sketches to form well-connected curve networks is both challenging and often undesirable. First, many strokes in 3D sketches are not meant to connect to others -see the details on author1_bulbasaur (Fig. 16), and on the car (Fig. 1 insets). Simply ignoring such disconnected strokes can alter the original design intent (Fig. 12, red strokes). Secondly, 3D sketches may exhibit oversketching and imprecisions (vr_controller Fig. 16), which make forming clean curve networks from 3D strokes a challenging problem, akin to 2D vectorization. We instead propose an approach that is oblivious to stroke connectivity, achieving high robustness to inaccuracy and to detail strokes that lie on the envisioned surface but are not connected to other strokes.

Our work also differs from interactive systems designed to model 3D surfaces by sketching in 2D [START_REF] Dvorožňák | Seamless Reconstruction of Part-Based High-Relief Models from Hand-Drawn Images[END_REF][START_REF] Li | BendSketch: Modeling Freeform Surfaces Through 2D Sketching[END_REF][START_REF] Nealen | FiberMesh: Designing Freeform Surfaces with 3D Curves[END_REF]]. Users of these systems draw in dedicated interfaces and provide annotations of surface discontinuities. In contrast, we take as input 3D sketches created with a variety of tools, and we propose a multi-model fitting algorithm to automatically locate sharp features where and only where they are needed.

Surfacing point clouds. The 3D sketches we consider can easily be converted to point clouds by sampling points along each stroke, which would allow leveraging the wealth of methods developed for surfacing such unstructured 3D data [START_REF] Berger | A survey of surface reconstruction from point clouds[END_REF]]. Unfortunately, most existing methods have been designed to process dense point clouds acquired using 3D scanning technology [START_REF] Hoppe | Surface reconstruction from unorganized points[END_REF][START_REF] Michael Kazhdan | Poisson surface reconstruction[END_REF], and are doomed to fail on 3D sketches that only provide a very sparse, non-uniform sampling of the envisioned surface. While some methods are robust to missing data, they only fix holes that are relatively small compared to the scale of the overall surface [START_REF] Hornung | Robust reconstruction of watertight 3D models from non-uniformly sampled point clouds without normal information[END_REF], or guide surface completion by fitting geometric primitives on otherwise dense parts of the point cloud [START_REF] Schnabel | Completion and reconstruction with primitive shapes[END_REF][START_REF] Tagliasacchi | Curve skeleton extraction from incomplete point cloud[END_REF]]. In contrast, 3D sketches are dominated by large holes, and only contain 3D information along thin, 1-dimensional subspaces.

The recent VIPSS algorithm of [START_REF] Huang | Variational implicit point set surfaces[END_REF] achieves impressive resilience to sparsity and non-uniformity of the point cloud, and has even been demonstrated on unstructured stroke clouds similar to our target sketches. But this robustness is obtained thanks to a global smoothness energy that misses sharp surface discontinuities. Nevertheless, we use VIPSS to initialize our method, and focus on the problem of segmenting its result into individual patches forming a more faithful reconstruction of the input sketch.

Our approach relates to algorithms that recover piecewise-smooth surfaces from point clouds by identifying locally-smooth patches and by detecting sharp discontinuities as the intersections of these patches [START_REF] Fleishman | Robust moving least-squares fitting with sharp features[END_REF][START_REF] Huang | Edge-Aware Point Set Resampling[END_REF][START_REF] Jenke | Patch-Graph Reconstruction for Piecewise Smooth Surfaces[END_REF]]. However, these methods rely on a dense sampling of the surface away from sharp features, while 3D sketches are mostly empty in such regions and are only densely sampled along feature curves.

We also take inspiration from methods that reconstruct or approximate 3D shapes as a collection of parametric geometric primitives, such as planes [START_REF] Cohen-Steiner | Variational shape approximation[END_REF][START_REF] Monszpart | RAPter: rebuilding man-made scenes with regular arrangements of planes[END_REF][START_REF] Trung | Geometrically Consistent Plane Extraction for Dense Indoor 3D Maps Segmentation[END_REF], cones, cylinders and spheres [START_REF] Li | Globfit: Consistently fitting primitives by discovering global relations[END_REF][START_REF] Wu | Structure Recovery via Hybrid Variational Surface Approximation[END_REF], generalized cylinders [START_REF] Zhou | Generalized Cylinder Decomposition[END_REF], and quadric surfaces [START_REF] Yan | Variational mesh segmentation via quadric surface fitting[END_REF]]. Similarly to these methods, we cast the discovery of representative surface patches as a multi-model fitting problem. Our originality is to leverage specifics of 3D sketches to guide this optimization, notably by encouraging large, uniform patches bounded by the input strokes, and by representing free-form surfaces as 4th-order implicit polynomial surfaces.

Note that the resulting segmentation is purely driven by geometry and does not carry any semantic meaning other than when semantic parts happen to be well-represented by different smooth patcheseg, the hood of the car and the windshield (Fig. 1c).

Surfacing stroke clouds. Our work focuses on sparse 3D sketches, which contrast with the dense VR paintings targeted by [START_REF] Rosales | SurfaceBrush: from virtual reality drawings to manifold surfaces[END_REF] that are created by covering the surface with large, overlapping ribbons that provide position and orientation samples all over the surface [Google 2016;[START_REF] Rosales | AdaptiBrush: Adaptive General and Predictable VR Ribbon Brush[END_REF]].

In the computer vision community, a few methods have been proposed to reconstruct curve clouds by matching edges in a multi-view stereo algorithm [START_REF] Fabbri | 3D Curve Sketch: Flexible Curve-Based Stereo Reconstruction and Calibration[END_REF]. [START_REF] Usumezbas | The Surfacing of Multiview 3D Drawings via Lofting and Occlusion Reasoning[END_REF] proposed a surfacing algorithm dedicated to such unstructured 3D data, where candidate surface patches are lofted between pairs of curves. But this method largely relies on the availability of multiple photographs of the shape to select valid candidate patches based on occlusion reasoning.

Closest to our problem statement, Batuhan Arisoy et al.

[2012] surface sparse and imprecise 3D sketches by smoothly deforming an initial low-fidelity surface of correct topology, using a discrete guidance vector field that points towards the closest stroke point. This approach produces globally smooth surfaces and requires user intervention to specify strokes that should be inserted into the mesh as sharp edge polylines (see Fig. 17 and 18 in their paper for a visual comparison of their results on similar sketches to ours). In contrast, our multi-model fitting formulation produces piecewise-smooth surfaces automatically.

OVERVIEW

Fig. 2 illustrates the main steps of our method, which takes as input a 3D sketch, along with an approximate proxy surface obtained with an automatic surfacing algorithm [START_REF] Huang | Variational implicit point set surfaces[END_REF]] or an existing lowpoly modeling tool (see the accompanying video for a demonstration of this workflow). After projecting the sketch onto the proxy, our goal is to segment the proxy into regions, each associated with a smooth surface model, such that the resulting piecewise-smooth surface satisfies the following desiderata.

-Reproducing stroke geometry. The surface models should run close to the input strokes, both for strokes that depict sharp surface discontinuities and for strokes that depict details within smooth areas. -Aligning patch boundaries with strokes. The boundary between neighboring regions should lie along strokes that depict sharp surface discontinuities, yet not all strokes in the sketch depict a discontinuity. -Keeping the reconstruction simple. The surface should be composed of a small number of smooth patches, rather than many intricate patches that would overfit to inaccuracy in the input strokes.

We formulate these competing requirements as energy terms in an optimization. A first energy term measures the distance between each stroke and the surface model it is assigned to. A second term measures the smoothness of the segmentation away from the strokes, encouraging the transitions between models to occur along strokes. This smoothness term also penalizes small, isolated regions, which contributes to satisfy our third desiderata. Finally, a third term measures the complexity of the reconstruction by counting the degrees of freedom of the models used.

While each of our terms has an intuitive interpretation, their combination yields a challenging optimization problem that combines discrete variables (which surface model should be assigned to which region of the proxy) and continuous variables (parameters of these models). We tackle this challenge by expressing our problem within PEARL, a general algorithm to solve multi-model fitting problems [START_REF] Isack | Energy-based geometric multi-model fitting[END_REF]. In our context, the algorithm alternates between optimizing the discrete variables describing the segmentation while keeping the model parameters fixed, and optimizing for the continuous model parameters while keeping the segmentation fixed (Fig. 2c).

After convergence, we obtain a segmented proxy, as well as a set of surface models represented by implicit polynomial surfaces. The last step of our method aims at converting this representation into a triangle mesh suitable for downstream applications (Fig. 2d). Since the proxy provides us with a mesh of correct topology relatively close to the surface models, we recover the final surface by projecting this mesh onto the zero level-sets of the polynomials.

METHOD

We now formulate the main terms of our optimization problem, before describing how we solve this problem by alternating segmentation and fitting. In what follows, we use the terms points and segments to refer to the stroke polylines, and the terms vertices and edges to refer to the triangles of the proxy mesh.

Problem formulation

While our segmentation algorithm relies on the proxy mesh as a computational domain, several of our energy terms are defined as functions of the stroke points. To reason about stroke geometry on the proxy domain, we project each stroke to its nearest location on the proxy and associate it with nearby mesh elements. Specifically, we detect intersections between the projected stroke segments and triangle edges, and associate the corresponding stroke points with the two vertices of the intersected edges. Thanks to this association, our segmentation algorithm assigns one surface model per mesh vertex, yet can fit several surface models to a single stroke point -which is critical for positioning boundaries of surface patches along strokes.

stroke proxy projection

Denoting P the set of all stroke points, we denote P 𝑣 ⊂ P the subset associated with a given mesh vertex 𝑣 ∈ V, which can contain zero or many stroke points (see inset). This association is built by a closest-point projection of the strokes onto the proxy surface (see Appendix A for more details).

Our goal is to assign each mesh vertex to a surface model 𝑓 𝛼 , associated with the label 𝛼 ∈ L. We denote the corresponding vertex labeling as 𝑙 : V → L. Furthermore, we denote as 𝜽 𝛼 the vector of real-valued parameters of the surface model 𝑓 𝛼 . Our surface model representation, detailed in Section 4.2, is the zero level-set Z (𝑓 𝛼 ) of an implicit polynomial. Finally, we denote as Θ = 𝛼 ∈ L 𝜽 𝛼 the concatenation of the parameter vectors 𝜽 𝛼 for all 𝛼 ∈ L.

Given these definitions, we cast our problem as finding the labeling and the associated model parameters that minimize

𝐸 L (𝑙, Θ) = 𝐸 fidelity (𝑙, Θ) + 𝑤 smoothness • 𝐸 smoothness (𝑙) + 𝑤 simplicity • 𝐸 simplicity (𝑙), (1)
where 𝐸 fidelity , 𝐸 smoothness , and 𝐸 simplicity capture the three desiderata listed in Section 3. Note that the number of labels required for a given sketch is also unknown, since we do not know how many smooth patches are necessary to faithfully represent the 3D surface a priori. We next describe each of the three energy terms.

Fidelity to input strokes. We seek a piecewise-smooth surface that reproduces well the input stroke geometry. We measure this property by summing over all mesh vertices 𝑣 the distance between their surface model Z 𝑓 𝑙 (𝑣) and each of their associated stroke points p ∈ P 𝑣 :

𝐸 fidelity (𝑙, Θ) = 1 𝜖 2 𝑁 𝑝 ∑ 𝑣 ∈V ∑ p∈ P 𝑣 dist Z 𝑓 𝑙 (𝑣) , p 2 , (2)
where the normalization term 𝑁 𝑝 is the total number of stroke points involved in the summation. We detail in Section 4.2 how to compute dist (Z(𝑓 ), p) 2 efficiently. The scale factor 𝜖 controls the sensitivity of 𝐸 fidelity to small deviations of the strokes from the surface models. We experimentally set this parameter to 1% of the input sketch's bounding-box diagonal for all results presented in the paper, but we show the impact of alternate values in Fig. 9.

Smoothness of the labeling. Our smoothness term seeks to concentrate changes of labels along strokes, such that sharp surface discontinuities appear at these locations when transitioning from one surface model to another. Inspired by classic work on edgeaware image segmentation [START_REF] Boykov | Interactive graph cuts for optimal boundary and region segmentation of objects in N-D images[END_REF], we encourage neighboring vertices to share the same label, unless they are separated by a stroke. We achieve this goal by assigning a weight 𝑤 𝑒 to each edge 𝑒 ∈ E of the mesh, using a low weight 𝑤 𝑒 = 1 if the edge is crossed by a projected stroke segment, and a high weight 𝑤 𝑒 = 100 otherwise. Here we assume that all edges have approximately the same length, see Appendix C for a scaling term in the case of anisotropic meshes. We then define the smoothness energy:

𝐸 smoothness (𝑙) = 1 𝑊 ∑ {𝑢,𝑣 } ∈ E 𝑤 𝑢𝑣 (1 -𝛿 (𝑙 (𝑢), 𝑙 (𝑣))) , (3)
with 𝛿 being the Kronecker delta function, 𝛿 (𝑖, 𝑗) = 1 if 𝑖 = 𝑗, 𝛿 (𝑖, 𝑗) = 0 otherwise. The normalization term 𝑊 is the sum of all edge weights. We set 𝑤 smoothness = 10.

Surface simplicity. To be resilient to the typical inaccuracy of 3D sketches [START_REF] Arora | ILoveSketch: as-naturalas-possible sketching system for creating 3d curve models[END_REF], we encourage the surface reconstruction to be as simple as possible. We measure the complexity of a solution as the total number of surface parameters involved. Denoting dim(𝜽 𝛼 ) the number of parameters of a surface model 𝑓 𝛼 , we define the simplicity energy as:

𝐸 simplicity (𝑙) = 1 𝐷 ∑ 𝛼 ∈ {𝑙 (𝑣) |𝑣 ∈V } dim(𝜽 𝛼 ), (4)
where we set 𝐷 = 35, the number of parameters of a degree 4 model, as a normalizer. Note that by summing over assigned models only, this energy also pushes for using a small number of models to explain the 3D sketch, to the point where a single model can be used to represent multiple non-adjacent regions of the surface. For all results, we set 𝑤 simplicity = 0.01 (and demonstrate varying values in Fig. 8).

Energy minimization

Algorithm 1 adapts the general PEARL multi-model fitting algorithm [START_REF] Isack | Energy-based geometric multi-model fitting[END_REF] to estimate the labeling 𝑙 and surface model parameters Θ that locally minimize Equation 1. After an initialization phase, the algorithm alternates between improving the

L ←-L ∪ L new Θ ←-Θ ∪ Θ new until 𝑙 = 𝑙 prev ;
labeling and refining the model parameters, and terminates when the labeling 𝑙 no longer changes. We next detail each step.

Segmentation. At each iteration, we first optimize the current labeling 𝑙 while keeping all surface model parameters Θ fixed. We perform 𝛼-expansion [Boykov et al. 2001] -which finds the optimal change of labeling where some nodes are assigned a label 𝛼 -for every label 𝛼 ∈ L to efficiently find a local minimum of 𝐸 L . 1 While 𝐸 smoothness and 𝐸 simplicity are straightforward to compute for a given labeling 𝑙, 𝐸 fidelity requires computing the distance dist (Z(𝑓 ), p) 2 between each surface model and the corresponding stroke points. Since the exact distance between a point and the zero level-set of an implicit polynomial surface cannot be computed exactly with a direct method, we employ the first-order approximation proposed by Taubin [1993],

dist (Z(𝑓 ), p) 2 ≈ 𝑓 (p) 2 ∥∇𝑓 (p)∥ 2 , ( 5
)
which can be computed in linear time with respect to the number of stroke points p. Since the resulting labeling might leave some labels unassigned, we update L to only keep the selected models.

Model representation and fitting. After each segmentation step, we improve the surface model of each label 𝛼 ∈ L by optimizing its parameters 𝜽 𝛼 to best fit the stroke points p ∈ P 𝛼 associated with the vertices 𝑣 ∈ V 𝛼 . We define each surface model as the zero where 𝜽 𝛼 = 𝜃 1 . . . 𝜃 𝑡 𝑇 is the 𝑡 × 1 vector of coefficients and X(𝑥, 𝑦, 𝑧) a 𝑡 × 1 vector of monomials, X(𝑥, 𝑦, 𝑧) = 1 𝑥 𝑦 𝑧 𝑥 2 . . . 𝑥 𝑑 𝑦 𝑑 𝑧 𝑑 𝑇 .

To avoid overfitting to approximate strokes, we limit the expressivity of the surface models by keeping their degree 𝑑 low. In practice we set 𝑑 to be at most 4, which corresponds to 𝑡 = 35 parameters. We describe at the end of this section how we introduce lower-degree models at the end of each iteration. Fig. 3 shows how surface models of different degrees capture shapes of varying complexity.

Given the set of stroke points p = (𝑥, 𝑦, 𝑧) ∈ P 𝛼 , we seek the vector of coefficients 𝜽 𝛼 that minimize the geometric error between the stroke points and the zero level set of the function 𝑓 . However, the geometric distance from a point to the zero level-set of a polynomial has no closed form expression, and minimizing it requires an iterative approach [START_REF] Joon Ahn | Orthogonal distance fitting of implicit curves and surfaces[END_REF]. Similarly, minimizing Equation 5 is costly since its derivatives are non-linear with respect to the model parameters 𝜃 . A common alternative consists in minimizing the algebraic distance to the implicit polynomial surface, which yields a linear least-squares regression:

arg min 𝜽 ∑ p∈ P 𝛼 X(p) 𝑇 𝜽 2 .
As shown by [START_REF] Tasdizen | Improving the stability of algebraic curves for applications[END_REF], this minimization problem can be made more stable by also encouraging the gradient of the implicit function to align with prescribed normals at a set of points

P ′ 𝛼 : arg min 𝜽 ∑ p∈ P 𝛼 X(p) 𝑇 𝜽 2 + 𝜇 ∑ p∈ P ′ 𝛼 1 -n p • ∇𝑋 (p) 𝑇 𝜽 2 ,
where ∇𝑋 (p) is a 𝑡 × 3 matrix denoting the gradient of the monomial vector 𝑋 (p), and n p denotes the normal vector at the point p. Since we do not have normal information at the stroke points, we encourage alignment with the proxy mesh normals at the mesh vertices. Furthermore, we exclude normals from vertices that lie close to the strokes, since strokes often denote surface discontinuities for which the smooth proxy mesh is only a crude approximation. We thus define the set of points for the alignment term as P ′ 𝛼 = {𝑣 ∈ V 𝛼 | P 𝑣 = ∅}. Finally, we also include an 𝐿 2 regularization term to penalize large polynomial coefficients [START_REF] Tasdizen | Improving the stability of algebraic curves for applications[END_REF]], yielding:

θ𝛼 = arg min 𝜽 ∑ p∈ P 𝛼 (X(p) 𝑇 𝜽 ) 2 +𝜇 ∑ p∈ P ′ 𝛼 (1-n p •∇𝑋 (p) 𝑇 𝜽 ) 2 +𝜆 ∥𝜽 ∥ 2 , (6)
with 𝜇 = 0.1 and 𝜆 = 1. Fig. 4 illustrates the impact of each of these regularization terms. We provide additional details on how to solve this minimization problem in Appendix B.

Note that θ𝛼 does not necessarily decrease 𝐸 fidelity (Equation 2) given a fixed labeling, since it minimizes the algebraic distance augmented with regularization terms -and not the geometric distance.

To guarantee convergence to a local minimum of Algorithm 1, we only update the parameter vector 𝜽 𝛼 if this yields a decrease of 𝐸 fidelity for the vertices labeled as 𝛼. If not, we keep the previous parameter vector and introduce the new parameter vector as a different model.

Initialization. Starting the algorithm with a good initial guess leads to higher quality solutions (see evaluation in Fig. 7). We obtain this initial guess by leveraging the observation that some of the strokes depict boundaries of smooth surface patches. Assuming that this is the case for all strokes, we compute an over-segmentation of the proxy surface by running spectral clustering [START_REF] Shi | Normalized cuts and image segmentation[END_REF] on the mesh, using the same edge weight 𝑤 𝑒 as in 𝐸 smoothness . We then fit an implicit polynomial surface to the strokes associated to the vertices of each cluster to obtain our initial set of models.

Spectral clustering requires the number of clusters to be specified in advance. Yet, the appropriate number varies among sketches, depending of their level of details. We address this challenge by running spectral clustering with different numbers of clusters (in practice, 20, 30, 40, and 50 clusters), which produces surface models of different scales.

For this initial surface fitting, we boost the regularization weights 𝜇 and 𝜆 by a factor 10 to limit complexity of the models and prevent the appearance of sub-optimal large models across sharp features.

New model proposal. We end each iteration by proposing new surface models, which helps decrease the energy 𝐸 L in subsequent iterations [START_REF] Isack | Energy-based geometric multi-model fitting[END_REF]. First, we introduce models of lower degree by fitting a polynomial of degree 𝑑 -1 to the set of stroke points P 𝛼 for each active model 𝑓 𝛼 of degree 𝑑 > 1. Second, we propose new models by merging pairs of neighboring regions and fitting a polynomial to the union of their stroke points. We prioritize merging regions that are separated by edges with a high weight 𝑤 𝑒 , as this strategy is more likely to yield a decrease in 𝐸 smoothness . In practice, we select three pairs of regions to be merged per iteration.

Extracting the surface mesh

The outcome of our multi-model fitting algorithm is a set of implicit surfaces that extend infinitely beyond the strokes they capture (Fig. 2c). Our goal is now to extract a triangular mesh from this arrangement of surfaces. One option to reach this goal would be to trim each surface along its intersection with other surfaces, for example by meshing the isosurfaces, resolving intersections [START_REF] Cherchi | Fast and Robust Mesh Arrangements using Floating-point Arithmetic[END_REF], and then selecting the set of patches that best cover the strokes while forming a closed manifold [START_REF] Bauchet | Kinetic shape reconstruction[END_REF][START_REF] Du | Boundary-Sampled Halfspaces: A New Representation for Constructive Solid Modeling[END_REF]]. Unfortunately, the halfspaces defined by our implicit polynomial surfaces are not guaranteed to bound the desired shape, as two surfaces that describe adjacent regions might not intersect. Even when neighboring surfaces do intersect, they can be nearly tangential to each other, which makes the detection of intersections prone to numerical inaccuracies. We bypass all these difficulties by leveraging the proxy mesh, as it provides us with a good estimate of the mesh we are looking for. We formulate our problem as the local minimization of an energy that projects the mesh towards the implicit surfaces assigned to each region, regularized by a smoothness term:

𝐸 mesh = 𝐸 project + 𝐸 regularize (7)
Insertion of segmentation boundaries. Before attempting to minimize Equation 7, we first need to insert edges that split the mesh triangles crossed by segmentation boundaries, so that these boundaries can emerge as sharp surface discontinuities in the optimized mesh. The inserted vertices (blue circles in inset) inherit the two or three labels of the regions they separate, yielding a multi-labeling function that we denote l.

Projection to the implicit surfaces. We project each mesh vertex onto its associated surface patches by minimizing its distance to the zero-level sets of the corresponding polynomials:

𝐸 project (V) = ∑ 𝑣 ∈V ∑ 𝛼 ∈ l (𝑣) 1 | l (𝑣)| dist (Z(𝑓 𝛼 ), 𝑣) 2 , ( 8
)
where we use the first-order approximation from Equation 5 to compute the distance.

Regularization. Minimizing Equation 8 sometimes pulls neighboring vertices in opposite directions, yielding a distorted surface.

We achieve smoother results by regularizing the optimization with a 2D Laplacian term on the mesh triangles, and a 1D Laplacian to favor smooth sharp feature curves:

𝐸 regularize (V) = 𝛾 2D ∑ 𝑣 ∈V w 𝑣 ||Δ 2D 𝑣 || 2 + 𝛾 1D ∑ 𝑣 ∈V seams ||Δ 1D 𝑣 || 2 , (9)
where V seams is the set of vertices inserted along segmentation boundaries, and Δ 2D and Δ 1D denote the discrete graph Laplacian operator for triangles and edges respectively [START_REF] Nealen | Laplacian mesh optimization[END_REF]].

We set 𝛾 2D = 100 and 𝛾 1D = 500. The weight w 𝑣 controls the strength of the surface regularization, which we want to vanish along segmentation boundaries that correspond to sharp surface discontinuities. We set w 𝑣 = 1 for vertices away from segmentation boundaries, and w 𝑣 = 0 for boundary vertices that are shared by more than two segmentation regions. For vertices that lie in-between two regions, we adjust their weight according to the angle formed by the implicit surface normals on each side of the boundary, with w 𝑣 = 0 when the angle is greater than 𝜋/3, w 𝑣 = 1 when the angle is smaller than 𝜋/8, and w 𝑣 varies linearly in-between. Effectively, this preserves sharpness of the boundary between regions that intersect sharply, while favoring smoothness elsewhere.

In Appendix D, we describe the special case of open meshes, which necessitate an additional term to prevent the mesh boundary to contract under the effect of the regularization.

We minimize 𝐸 mesh with L-BFGS, by precomputing the gradient matrices for the quadratic terms.

IMPLEMENTATION DETAILS

Proxy creation. Our method requires a proxy surface that has the same topology as the desired result, and that lies close to the envisioned surface, so that projecting the strokes onto the proxy yields the same embedding as it would on the envisioned surface.

The automatic VIPSS surfacing algorithm [START_REF] Huang | Variational implicit point set surfaces[END_REF]] produced suitable proxies for many of our results. We additionally re-meshed the output of VIPSS to obtain a high-resolution, uniform mesh [START_REF] Wenzel | Instant Field-Aligned Meshes[END_REF]]. However, VIPSS sometimes fails to surface complex sketches that exhibit many sharp features, a very sparse sampling, or thin features. In these cases, we experimented with multiple casual modeling tools to create our proxies in just a few minutes, including low-poly modeling with Blender [START_REF] Blender | Blender[END_REF]] (Fig. 5c), rough VR sculpting with Adobe Medium [Adobe 2020] (Fig. 5d), assembly of simple geometric primitives (Fig. 17 author1_bulbasaur). We provide an example modeling session in our accompanying video.

Manually-created proxy meshes are often approximate, and would not yield satisfying results if projected directly on the implicit surface models (Section 4.3). We obtain better results by first attracting the proxy towards the strokes, which we achieve by projecting the stroke points onto the proxy and by using these points as control vertices for a Least-Squares mesh defined by the initial proxy connectivity [START_REF] Sorkine | Least-squares meshes[END_REF]. Fig. 5 shows how this attraction brings approximate proxy meshes much closer to the intended surface. We automatically improve the proxy mesh by remeshing, smoothing, and attracting it to the strokes (middle column). When the sketch is symmetric, we only process half of the proxy mesh. Note that despite a drastically different original appearance, all three proxy meshes yield a similar result (e), faithful to the input sketch (a).

Symmetric sketches. Many VR sketching and sketch-based modeling tools provide a mirror plane to ease the creation of symmetric sketches. When this is the case, we only surface half of the sketch and mirror the result, which speeds up computation and yields surfaces that are symmetric by construction.

Open surfaces. We support open surfaces by asking users to annotate strokes that denote the boundary of the envisioned surface (colored in blue in all input sketches), and by trimming the surface mesh along these strokes. We provide additional details about this procedure in Appendix D.

User control. Most results shown in this paper were obtained automatically. Nevertheless, an additional strength of our formulation is that it admits user control naturally. Fig. 6 illustrates two scribblebased controls that we offer users to refine the segmentation:

-Inserting a new region. Users can scribble over a part of the mesh to trigger the insertion of a new region in that part (Figure 6, yellow scribbles). This interaction is particularly useful to recover details that might have been over-smoothed. -Merging regions. Users can scribble over multiple regions to merge them into a single one (Figure 6, green scribbles). This interaction can be used to correct for over-segmentation.

We implement region insertion by introducing a new model, which we fit to the stroke points associated with the scribbled-on mesh vertices. We implement region merging in a similar manner, except that we fit the new model to the stroke points associated with every vertex that shares a label with any scribbled-on vertex and lies on the same connected component. For both edits, we also penalize the use of other models by adding a term to Equation 1: where V scribbled denotes the set of scribbled vertices, 𝑙 new denotes the label of the new model, and the penalty 𝜌 is set to (10𝜖) 2 , with 𝜖 from Section 4.1. We then update the solution by re-running algorithm 1 from step (1) until convergence. This combination of new model proposals and a penalties on other models effectively pushes the optimization towards a different local minimum that satisfies the user indications. We describe an edge-collapse algorithm in Appendix C that reduces the complexity of the mesh on which we compute the labeling, making our method react faster to user edits.

𝐸 penalty (𝑙) = 𝜌 ∑ 𝑣 ∈V scribbled 1 -𝛿 (𝑙 (𝑣), 𝑙 new ) (10)

EVALUATION AND RESULTS

Algorithm evaluation

Initialization strategy. We initialize our method by over-segmenting the proxy mesh using spectral clustering (Section 4.2). Fig. 7 compares this strategy with two baselines, over multiple runs of the method. The first baseline strategy creates the initial surface models by sampling random sets of points among all stroke points in the sketch. The second baseline applies spectral clustering but keeps the default regularization weights rather than scaling their values during the initial surface fitting.

This evaluation reveals that our strategy is more stable, as it makes the optimization converge to approximately the same energy for 1) of the converged state over 10 runs with different random seeds. Initializing the surface models with a random selection of stroke points yields high variance across runs (left column of each plot). Spectral clustering identifies local patches that form good candidate models, but the default regularization weights used to fit these models can result in complex surfaces that cover large portions of the sketch (middle column of each plot, red circles in segmentation results). Boosting the regularization weights favors simpler models at the start of the optimization, yielding better and more stable results after convergence (right column of each plot).

different random seeds. In addition, our strategy is more effective, as evidenced by the lower or comparable median energy values at convergence. Finally, visualizing the segmentation produced by each strategy highlights the benefit of boosting the regularization when fitting the initial models, as it prevents the creation of complex models covering large regions of the sketch, which would be difficult to remove in subsequent steps of the optimization. While our initialization strategies guarantees low variance and similar results regardless of random seed (see bottom two inset results), for all subsequent results we reduce the influence of initialization by running our method 10 times and selecting the solution with the lowest energy.

Segmentation parameters. Fig. 8 shows the effect of our simplicity term (Equation 1) by varying its weight. The result is oversegmented when we omit the simplicity term altogether (Fig. 8a). Increasing this parameter (Fig. 8b,c) leads to smaller model counts 𝑚, and to models of lower degrees (see the spaceship cockpit), at the cost of an increase in mean fitting residual ε. The simplicity energy also encourages using the same model across disconnected regions (Fig. 8c, yellow ellipsoid), which cannot be achieved with the smoothness energy alone.

The other tunable parameter of our method is the factor 𝜖 that defines the scale of the fitting error we tolerate (Equation 2). Increasing this factor leads to results that follow the strokes more loosely, which can be useful if the input sketch is imprecise and needs to be smoothed (Fig. 9).

Influence of proxy mesh. Fig. 10 illustrates the influence of the proxy resolution and shape on our results. Since the proxy serves as a discrete domain for our segmentation algorithm, its resolution impacts the size of the regions we can detect. At low resolution, neighboring strokes will be lumped together and sharp features will be misplaced Fig. 8. We demonstrate the impact of the simplicity energy by varying the weight 𝑤 simplicity in Equation 1 on two sketches. We measure the number 𝑚 of models used in the segmentation, and the mean deviation ε between stroke points and the models of the vertices they are associated with, as a % of the bounding box diagonal. (Fig. 10a, third row). In addition, a low-resolution mesh might lack degrees of freedom to follow the curvature of the implicit surfaces, and is thus unable to reproduce the shape of smooth high curvature parts of the sketch (Fig. 10a, fourth row).

The proxy also serves to embed the strokes into a manifold domain representative of the envisioned surface. When the proxy is too far from the desired surface, strokes that should be disjoint might end up projecting to the same location on the proxy and be approximated by the same surface patch (Fig. 10b, top). A slightly more precise proxy mesh resolves the issue (Fig. 10b,bottom).

0% > 5%

Fig. 11. We evaluate our method on synthetic sketches extracted from ground truth surfaces using FlowRep [START_REF] Gori | FlowRep: Descriptive Curve Networks for Free-Form Design Shapes[END_REF]] (bottle, boat) or from the boundary representation of CAD models (iron, bishop). The deviation is measured as a % of the sketch bounding-box diagonal, and is low almost everywhere except in places where the sketch lacks strokes to disambiguate the intended surface (interior of the iron handle, bottom of the bishop). We provide the median and maximum deviation for each sketch.

Table 1. For each sketch in the paper, we provide the runtimes for the initialization (1) -which corresponds to the strokes-proxy association (Appendix A) and step 0 of algorithm 1; for the iterative segmentation (2), and for the final mesh optimization (3). We also give the number of iterations of algorithm 1 necessary to converge to a stable solution. We use a light gray material to shade proxy meshes computed automatically with VIPSS [START_REF] Huang | Variational implicit point set surfaces[END_REF], and a dark gray material for the ones created manually. Note how our method recovers the fine details and sharp edges depicted in the sketch even from very smooth, approximate proxy meshes. Our supplementary website displays each result with an interactive 3D viewer. We detail in Table 1 the performance of our algorithm on all results shown in the main paper. Timings vary from a few seconds on simple sketches up to a few minutes on complex sketches composed of many strokes. The bottlenecks reside in the initialization of the algorithm and in the final mesh optimization, while the segmentation algorithm takes less than 10 seconds in most cases. Users can thus refine the segmentation multiple times with relatively short wait times and only compute the final mesh once satisfied.

Comparison against ground truth surfaces. We evaluate our method quantitatively by comparing our results to ground truth surfaces for which curve networks are available. We obtain such data from FlowRep [START_REF] Gori | FlowRep: Descriptive Curve Networks for Free-Form Design Shapes[END_REF]], a method to generate descriptive curve networks from 3D meshes, as well as from CAD models for which we extract the boundary representation (B-rep) [START_REF] Koch | Abc: A big cad model dataset for geometric deep learning[END_REF]. Fig. 11 visualizes the results of this experiment, where the color map indicates that our results are very close to ground truth except in ambiguous regions devoid of strokes. Sketchified network Ours Ours [START_REF] Huang | Variational implicit point set surfaces[END_REF] Ours [START_REF] Huang | Variational implicit point set surfaces[END_REF] + high frequency noise Ours [START_REF] Huang | Variational implicit point set surfaces[END_REF] [ [START_REF] Huang | Variational implicit point set surfaces[END_REF] Fig. 12. Our method can surface well-connected curve networks, similarly to [START_REF] Hao Pan | Flow aligned surfacing of curve networks[END_REF]. More importantly, it also achieves similar results from imprecise sketches that contain duplicate strokes, gaps (middle row) or high-frequency noise (bottom row). Our method also accounts for disconnected parts, like the red stroke that depicts a concavity on the car, which [START_REF] Hao Pan | Flow aligned surfacing of curve networks[END_REF] ignored.

Comparison against curve network surfacing. Fig. 12 compares our method to the one by Pan et al. [2015], which is the most recent method for surfacing well-connected curve networks. The two methods produce very similar results on the original connected network. However, a unique strength of our method is to also produce similar results from disconnected, noisy sketches, which we synthesized by perturbing the original network (displacing and duplicating strokes, introducing gaps). Our method also accounts for strokes that are not connected to the main network, which [START_REF] Hao Pan | Flow aligned surfacing of curve networks[END_REF] had to ignore (Fig. 12, red strokes at the front of the car).

In Fig. 13, we compare to Pan et al.

[2015] on a more complex input. By balancing fidelity with smoothness and simplicity, our method tends to miss small details, such as the buttons on top of the machine. Our smooth surface models also do not capture well the generalized cylinder that forms the nozzle. In contrast, by assuming that the input strokes form a clean curve network, Pan et al. can trust every curve to create interpolating surfaces and can leverage connectivity information to detect which curves are sharp features.

Limitations

Relying on geometric criteria only. We emphasize that our method relies purely on geometric criteria to place sharp features in the final surface. As a consequence, our method can miss semanticallyimportant features if they do not contribute significantly to the shape, such as the round headlights of the car in Fig. 1.

Internal strokes Non manifold

Separate in parts Ignore strokes

Skeletal strokes

Fig. 15. Our current implementation does not automatically support nonmanifold surfaces, as present on the boat and blender (red strokes). Furthermore, our assumption that strokes lie on the surface they depict is not always true, as there may be some strokes inside the outermost depicted surface (green) or strokes that depict tubular structures (cyan). A workaround is to manually separate the sketch into parts and surface them individually, or mark out some strokes to be ignored by our method (second row).

Sharp features not depicted by the strokes. Our algorithm assumes that all sharp features of the intended surface occur along some of the input strokes. Fig. 14 illustrates a limitation of this assumption, where the artist implicitly indicated that the headlight presents a sharp feature by sharp corners in neighboring strokes. In such cases, users need to draw a few additional strokes to obtain the intended surface (Fig. 14b). While our assumption that sharp features are explicitly represented seems to hold in our dataset, further analysis of how people draw in 3D is needed to quantitatively evaluate this hypothesis. In the future, more global cues such as stroke tangent continuity or sketch symmetry could be leveraged to make surfacing better aligned with viewer perception.

Other types of strokes. Since we assume that all strokes lie on the intended surface, our method cannot handle strokes that lie inside the shape, or that depict the skeleton of a tubular structure (Fig. 15). While we simply discarded such strokes to produce our results (see supplemental materials for all original input sketches), future work could attempt to identify these strokes and surface them with dedicated representations, such as generalized cylinders for skeletal stokes [START_REF] Zhou | Generalized Cylinder Decomposition[END_REF].

Non manifold configurations. In theory, our method could recover non-manifold surfaces if provided with a non-manifold proxy mesh. In practice, we only surfaced manifold shapes because we implemented our algorithm using a manifold data-structure to represent the mesh. As a work-around, it is typically possible to manually separate the sketch into multiple manifold pieces that can be surfaced separately (Fig. 15, bottom row and coffee cup in Fig. 13). An exciting direction for future work would be to perform topology analysis of the unstructured sketch to detect and process non-manifold parts automatically.

Iterative sketch creation and surfacing. We did not consider the problem of surfacing a sketch as it is being created and refined live by an artist. While our method can be made more efficient to reduce computation time with more engineering effort, it is not trivial how to take into account new strokes or edits on an existing stroke without globally affecting the surface. Moreover, partial sketches introduce more ambiguity to the surfacing problem since some regions may be altogether undefined at a given stage of the process.

CONCLUSION

After decades of research, a number of robust algorithms now exist to surface dense point clouds [START_REF] Berger | A survey of surface reconstruction from point clouds[END_REF]]. In contrast, very few methods have been proposed to surface the sparse stroke clouds that designers produce when sketching in VR or with sketchbased modeling systems. Inspired by the unique characteristics of this emerging form of 3D data, we have proposed an approach to locate smooth patches in unstructured 3D sketches, and to optimize their geometry to produce a piecewise-smooth surface aligned with salient strokes. The resulting 3D meshes can benefit numerous tasks. For instance, we used them to occlude hidden strokes in all figures of this paper, which helps perceive the correct shape from the sketch via relative depth cues [START_REF] Arora | ILoveSketch: as-naturalas-possible sketching system for creating 3d curve models[END_REF]. Additionally, our sketch-aligned surfaces are compatible with all downstream 3D processing and modeling tasks. This opens exciting avenues to integrate 3D sketching with other creation modalities such as sculpting [Adobe 2020;[START_REF] Blender | Blender[END_REF]Pixologic 2016]. As we have demonstrated, sculpting can be used to create a rough proxy surface, that is then automatically refined by our method to align with a sketch. It is easy to further refine the surface by sculpting, and in the future we would like to investigate tighter coupling between surface and curve-based editing methods, from both interaction design and geometric representation perspectives.

B REGULARIZED FITTING OF IMPLICIT POLYNOMIALS

From Equation 6, we obtain the parameters θ𝛼 that best fit the stroke points by solving the linear system (𝑀 𝑇 𝑀 + 𝜇𝑀 𝑇 align 𝑀 align + 𝜆𝐼 ) θ𝛼 = 𝜇𝑀 𝑇 align .

(11)

𝑀 is a (|P 𝛼 | × 𝑡) matrix composed of rows of the monomial vector 𝑋 (p) for every point in P 𝛼 . And given the set of points P ′ 𝛼 for the alignment term, 𝑀 align is a (1 × 𝑡) vector defined as:

𝑀 align = ∑ p∈ P ′ 𝛼 n p 𝑇 ∇X(p). (12)
To keep the 𝐿 2 regularization position and scale independent, and resilient to non uniform scaling in the point set P 𝛼 , we apply a transformation to the space of monomials composed of a translation -to center the data -and a non uniform scaling -to normalize it.

In practice, we compute the column-wise mean values 𝜇 and standard deviations 𝑆 of 𝑀, and consider the transformed polynomial function:

𝑓 (𝜽 𝛼 ; 𝑥, 𝑦, 𝑧) = X(𝑥, 𝑦, 𝑧) -𝜇 𝑆 𝑇 θ𝛼 ,

where the division by the vector 𝑆 is a column-wise division. As a consequence, the vectors 𝑋 (p) that compose matrix 𝑀 are transformed by subtracting 𝜇 and dividing by 𝑆, and the vector 𝑀 align is divided by 𝑆.

C GRAPH SIMPLIFICATION

Algorithm 1 has a complexity that increases linearly with the number of nodes in the graph to be labeled [Boykov et al. 2001]. Yet, we note that only the vertices that are associated to stroke points contribute to 𝐸 fidelity , other vertices being solely determined by the smoothness and simplicity terms. We leverage this observation to simplify the graph away from the strokes, which we achieve by performing greedy edge collapses [START_REF] Garland | Surface simplification using quadric error metrics[END_REF][START_REF] Hoppe | Surface reconstruction from unorganized points[END_REF]]. In our context, we are not interested in preserving the 3D shape of the original mesh, but rather its connectivity. Therefore, we prioritize collapses that yield vertices of low valence, and stop collapsing edges whenever any further collapse would yield a vertex of high valence (> 12 in practice). We exclude edges that have a vertex associated with stroke points to preserve high resolution near the strokes. This simplification produces a mesh with spatially-varying resolution, which we account for by adjusting the edge weight to w𝑒 = 𝑤 𝑒 𝑐 𝑒 , where 𝑤 𝑒 is the original edge weight as defined in Section 4.1, and 𝑐 𝑒 scales this weight according to the size of the triangles adjacent to the edge, as measured by the Euclidean distance between the vertices opposite the edge on the two adjacent faces.

After running algorithm 1 on the simplified graph, we propagate the labeling to the vertices of the original graph based on proximity. This simple strategy speeds up Algorithm 1 by a factor of 2.9 on average, while yielding qualitatively similar results. where V boundary is the set of vertices associated to the boundary strokes of the sketch. However, we do not know a priori where an open surface should be trimmed to align with the sketch boundary. We thus ask users to annotate the boundary strokes of sketches depicting open surfaces (shown in blue in all sketches). After projecting these strokes over the proxy mesh, we segment the mesh into an interior and an exterior region separated by the boundary strokes (Fig. 19c), and trim the exterior region to obtain a mesh whose boundary V boundary aligns with the boundary strokes. We perform this segmentation with a graph cut, where we encourage cuts along edges crossed by boundary strokes by setting their cost to zero (Fig. 19b). The graph cut source node is linked to all boundary vertices of the input mesh (red in Fig. 19b) and the sink node is linked to vertices of the mesh that have associated stroke points but are geodesically far from the boundary strokes (green in Fig. 19b).

D BOUNDARIES OF OPEN MESHES.

For ease of comparison, we also apply this boundary cutting step to results of [START_REF] Huang | Variational implicit point set surfaces[END_REF] against which we compare.

Fig. 2 .

 2 Fig.2. Overview of our method. We take as input a 3D sketch (a) and a proxy mesh (b), for example obtained with VIPSS[START_REF] Huang | Variational implicit point set surfaces[END_REF]. We iteratively improve a segmentation of the proxy mesh (c, left) and parameters of a set of surface models (c, right) to obtain a decomposition of the surface into smooth patches that approximate the input strokes well. We represent each surface model as the zero level set of an implicit polynomial, which offers infinite spatial support and a good balance between expressivity and complexity. Finally, we project the vertices of the proxy mesh onto the surface models to recover the final piecewise-smooth surface with sharp features (d). When appropriate, we trim the proxy mesh to match user-annotated boundary strokes (a, blue), and we leverage sketch symmetry by surfacing only half of the sketch and mirroring the result.

Fig. 3 .

 3 Fig.3. We represent surface patches as zero level-sets of implicit polynomials of degree up to 4. Low-degree polynomials can represent planes, developable and simple doubly-curved surfaces, while higher-degree polynomials can represent complex freeform surfaces.

Fig. 4 .

 4 Fig. 4. Fitting an implicit polynomial surface to a subset of strokes (thick red strokes). We visualize a slice of the unsigned distance field in the second row (distances consider a sketch bounding box diagonal of unit length). (a) Minimizing the algebraic distance without regularization yields level sets with spurious turns close to the data points. (b) By encouraging the gradient to align with the proxy normals [Tasdizen et al. 2000], the implicit function is well behaved in a neighborhood around the data points. (c) Adding 𝐿 2 regularization avoids overfitting to noise in the data points and yields a more stable zero level-set.

Fig. 5 .

 5 Fig. 5. For the same sketch, we demonstrate how a proxy mesh can be created with VIPSS (b), low-poly modeling with Blender [Blender 2020] (c), and rough VR sculpting in Adobe Medium [Adobe 2020] (d).We automatically improve the proxy mesh by remeshing, smoothing, and attracting it to the strokes (middle column). When the sketch is symmetric, we only process half of the proxy mesh. Note that despite a drastically different original appearance, all three proxy meshes yield a similar result (e), faithful to the input sketch (a).

Fig. 6 .

 6 Fig. 6. Our automatic method might produce over-segmented regions or miss geometric details (b). The user can indicate desired changes by scribbling on a preview of the result (c). Yellow scribbles trigger region insertion, green scribbles trigger region merging. Encouraging the segmentation to respect these constraints yields the desired result (d). Left column sketch ©Jacopo Colò, used with permission.

Fig. 7 .

 7 Fig. 7. Comparison of different initialization strategies. Each plot shows the energy (Equation1) of the converged state over 10 runs with different random seeds. Initializing the surface models with a random selection of stroke points yields high variance across runs (left column of each plot). Spectral clustering identifies local patches that form good candidate models, but the default regularization weights used to fit these models can result in complex surfaces that cover large portions of the sketch (middle column of each plot, red circles in segmentation results). Boosting the regularization weights favors simpler models at the start of the optimization, yielding better and more stable results after convergence (right column of each plot).

Fig. 9 .

 9 Fig. 9. We demonstrate the impact of the fidelity energy by varying the error scale factor 𝜖 in Equation 2 on the spaceship sketch (top row). In bottom rows, we show concrete examples of sketches with imprecise strokes where a higher 𝜖 gives better results compared to the default 𝜖 = 1% parameter. The error scale factor is measured in % of the sketch bounding box diagonal.

Fig. 10 .

 10 Fig. 10. Our method depends on the resolution of the input proxy surface (a) as well as on its shape (b). The proxy should have sufficient resolution to capture fine details, and be close enough to the envisioned surface to avoid thin parts of the sketch to collapse. Left column sketch ©Arturo Tolentino, used with permission.

6. 2

 2 Fig. 16 and Fig. 17provide a gallery of results obtained by surfacing 3D sketches from varied sources, ranging from clean curve networks created with 2D sketch-based modeling interfaces [Bae et al. 2008; Xu et al. 2014], to imprecise and over-sketched stroke clouds created with 2D and VR ideation tools [Colò 2021; Kim and Bae 2016; Yu et al. 2021a]. All these results were obtained without user correction.We use a light gray material to shade proxy meshes computed automatically with VIPSS[START_REF] Huang | Variational implicit point set surfaces[END_REF], and a dark gray material for the ones created manually. Note how our method recovers the fine details and sharp edges depicted in the sketch even from very smooth, approximate proxy meshes. Our supplementary website displays each result with an interactive 3D viewer.We detail in Table1the performance of our algorithm on all results shown in the main paper. Timings vary from a few seconds on simple sketches up to a few minutes on complex sketches composed of many strokes. The bottlenecks reside in the initialization of the algorithm and in the final mesh optimization, while the segmentation algorithm takes less than 10 seconds in most cases. Users can thus refine the segmentation multiple times with relatively short wait times and only compute the final mesh once satisfied.

  Fig. 16 and Fig. 17provide a gallery of results obtained by surfacing 3D sketches from varied sources, ranging from clean curve networks created with 2D sketch-based modeling interfaces [Bae et al. 2008; Xu et al. 2014], to imprecise and over-sketched stroke clouds created with 2D and VR ideation tools [Colò 2021; Kim and Bae 2016; Yu et al. 2021a]. All these results were obtained without user correction.We use a light gray material to shade proxy meshes computed automatically with VIPSS[START_REF] Huang | Variational implicit point set surfaces[END_REF], and a dark gray material for the ones created manually. Note how our method recovers the fine details and sharp edges depicted in the sketch even from very smooth, approximate proxy meshes. Our supplementary website displays each result with an interactive 3D viewer.We detail in Table1the performance of our algorithm on all results shown in the main paper. Timings vary from a few seconds on simple sketches up to a few minutes on complex sketches composed of many strokes. The bottlenecks reside in the initialization of the algorithm and in the final mesh optimization, while the segmentation algorithm takes less than 10 seconds in most cases. Users can thus refine the segmentation multiple times with relatively short wait times and only compute the final mesh once satisfied.

[

  Fig. 13. Our method tends to miss small yet important details (middle) in an effort to produce smooth, simple outputs robust to imprecise drawings. The method by Pan et al. [2015] does not face this challenge as they consider that all input curves are drawn precisely, and as such should be kept in the output. Similar to Pan et al. [2015], we surface the coffee cup separately from the coffee machine.

  Fig.14. Our method cannot take higher-level cues into account, and will thus fail to align with user intent if a critical sharp feature stroke is missing, as on the car headlight (a). Adding a stroke and re-running the method successfully recovers that feature (b).

Fig. 16 .

 16 Fig. 16. Results of our method applied to a variety of 3D sketches (a). We display proxy meshes computed with VIPSS in light gray, and the ones created manually in dark gray (b). For each sketch, we show the segmented patches with random colors, and the final surface in blue (c). Sketch robbins_car2 ©James Robbins, used with permission.

Fig. 17 .Fig. 18 .

 1718 Fig. 17. Results of our method applied to sketches of organic shapes (a). We display proxy meshes computed with VIPSS in light gray, and the ones created manually in dark gray (b). For each sketch, we show the segmented patches with random colors, and the final surface in blue (c). Sketch tolentino_shoe ©Arturo Tolentino, used with permission.

Fig. 19 .

 19 Fig. 19. To locate the boundary of an open surface, we ask the user to indicate which strokes correspond to boundaries in the input sketch (a, blue strokes).We perform a graph cut on the mesh to separate inside nodes (b, green) from outside nodes (b, red), while encouraging cuts to happen at edges crossed by boundary strokes (b, blue edges). The resulting labeling separates the inside of the surface (c, white) from the outside (c, red).

  While the Laplacian regularization in Equation 9 favors smooth, regular meshes, it has the adversarial effect of shrinking open meshes. We prevent this effect by adding to 𝐸 mesh a term that constrains boundary vertices of open meshes to stay close to their associated stroke points 𝐸 attract = ∑ 𝑣 ∈V boundary ∑ p∈ P 𝑣 ||𝑣 -p|| 2 (13)

  𝛼-expansion[Boykov et al. 2001] ∀𝛼 ∈ L to optimize previous labeling 𝑙 prev for energy (Equation1), given the current parametric surface models Θ L ←-{𝛼 = 𝑙 (𝑣) | 𝑣 ∈ V} // Remove unassigned labels Θ ←-𝛼 ∈ L 𝜽 𝛼 // And the corresponding inactive models (2) Optimize parametric surface models for 𝛼 ∈ L doV 𝛼 ←-{𝑣 ∈ V | 𝑙 (𝑣) = 𝛼 } //Vertices with label 𝛼 𝜽 𝛼 ←-Best fit for points in P 𝛼 = 𝑣 ∈V 𝛼 P 𝑣 end (3) Propose new models Propose new labels L new and models Θ new

	Algorithm 1: PEARL algorithm [Isack and Boykov 2012]
	applied to our problem
	(0) Initialize variables
	Propose 𝑚 vectors of parameters {𝜽 0 1 , • • • , 𝜽 0 𝑚 } as initial
	surface models
	Θ ←-𝑚 𝑖=1 𝜽 0 𝑖 // Concatenated surface parameters
	L ←-{1, . . . , 𝑚} // Set of existing labels
	𝑙 ←-𝑙 : 𝑣 ∈ V ↦ → 1 // Initialize with single label
	repeat
	𝑙 prev ←-𝑙
	(1) Optimize labeling
	𝑙 ←-

ACM Trans. Graph., Vol. 41, No. 4, Article 88. Publication date: July 2022.

ACM Trans. Graph., Vol. 41, No. 4, Article 88. Publication date: July 2022.Piecewise-Smooth Surface Fitting onto Unstructured 3D Sketches • 88:5

We use the C++ multi-label optimization library of[START_REF] Delong | Fast approximate energy minimization with label costs[END_REF] https://vision. cs.uwaterloo.ca/code/ ACM Trans. Graph., Vol. 41, No. 4, Article 88. Publication date: July

ACM Trans. Graph., Vol. 41, No. 4, Article 88. Publication date: July 2022. Piecewise-Smooth Surface Fitting onto Unstructured 3D Sketches • 88:7

ACKNOWLEDGMENTS

We thank the artists who kindly shared some of their artworks with us, James Robbins, Jacopo Colò, Arturo Tolentino. We thank the anonymous reviewers for their helpful comments. This work was supported by ERC Starting Grant D3 (ERC-2016-STG 714221), NSERC, Advokat Bent Thorbergs Fond (66.531) and by software and research donations from Adobe.

A ASSOCIATING 3D STROKES WITH THE PROXY MESH

We provide additional implementation details on how we associate stroke points 𝑝 ∈ P and proxy mesh vertices 𝑣 ∈ V.

We start by projecting each stroke point to its closest face on the proxy mesh. We then trace out stroke segments as geodesics lying on the mesh surface. Note that we do not need to trace precise geodesics on the mesh, since we are only interested in detecting which mesh edges are crossed by a given stroke segment. Therefore, we approximate the geodesics with shortest paths in the dual graph of the mesh (light blue path in Fig. 18a). This gives us the list of primal mesh edges that are crossed by the projected stroke segment.

We can then associate stroke points to these crossed edges (Fig. 18, red) by sampling additional points on the stroke segment (Fig. 18, blue). If 𝑛 edges are crossed while tracing a segment, we sample 𝑛 regularly-spaced points on that segment. Each of these points is then associated to the corresponding mesh edge.

Finally, we associate to each vertex 𝑣 𝑖 the stroke points associated with each edge of its 1-ring. If there are multiple stroke points associated to a single edge, we only associate the closest point to 𝑣 𝑖 .

The projection of stroke segments onto the proxy also serves to define the edge weights of the graph, since we assign a low weight to edges crossed by a projected segment -red edges in Fig. 18b.