Visualization of workpieces during simulation of milling processes

[0001] This application claims the priority of German Patent Application, Serial No. 101 44 932.1, filed Sep. 12, 2001, pursuant to 35 U.S.C. 119(a)–(d), the subject matter of which is incorporated herein by reference.

[0002] The present invention relates to a method for simulating a workpiece contour by simulating a milling operation as well as a corresponding device, in particular a simulation computer.

[0003] With CNC-controlled processing machines, such as milling machines, a workpiece is typically encoded either directly or the workpiece is modeled first with a CAD system and then converted into an equivalent CNC parts program. The resulting CNC parts program and/or the CAD model then correspond to the idealized processing commands for the processing machine. The CNC program is loaded into a CNC controller, and the processing machine is controlled according to the CNC program.

[0004] If the workpiece produced according to the CNC program is within the desired manufacturing tolerances of an ideal workpiece, then this method poses no problems. However, if the manufactured workpiece does not satisfy the requirements, then the process should be optimized to determine which changes are necessary, for example in the CNC program, to produce a perfect workpiece.

[0005] It may be possible to change sequentially individual processing commands and/or individual operating parameters of the processing machine, and to produce a new workpiece and then to test the newly produced workpiece. This approach, however, is very cumbersome and also expensive as well as material- and time-consuming. This is particularly the case because the causes for the deviations between the actually produced workpiece and the intended workpiece are often not known.

[0006] For this reason, mechatronic systems, such as processing machines, are more frequently simulated. This is particularly true for milling machines and the workpieces produced with milling machines. In order to analyze the result of such simulation, a visualization environment is required and hence a process for realistically simulating the surface of a workpiece obtained by simulation.

[0007] In particular with simulation of milling processes, the visualization component is quite important. This helps to represent data and results, such as milling points, milling paths, and coordinate systems, as well as the workpiece contours calculated by the simulation system.

[0008] Today's visualization components are capable of three-dimensional rendering and can be used interactively. The displayed content can be orthographic and/or perspective 3-D projections which can be changed interactively the user. The user is thereby able to rotate, translate and zoom the displayed content. By selecting certain partial content, easy to use visualization components can also produce associated information, such as dimensions, spatial arrangement or correspondence with other partial content of the scene.

[0009] Together, this leads to a better understanding of the production process. Moreover, the surface quality of the workpiece to be machined can be readily determined and analyzed in the design stage, allowing optimization of the current parameter set that controls and drives the machine tool.

[0010] This makes it possible to produce and machine a “virtual workpiece”, so that an actual workpiece does not have to be produced. Moreover, the processing machine doesn't even have to exist.

[0011] To make quantitative statements concerning the expected workpiece contour, a high-quality visualization is hence essential. The virtual workpiece contour should preferably be rendered as a 3-D model. A representation of the milling points and paths using just points and lines is therefore not sufficient. The workpiece contour should be computed across an area from the milling data. This requires a method for producing a surface grid, which is supplied to the visualization component and subsequently rendered.

[0012] Known methods for visualizing a workpiece surface typically process only interpolated milling points, i.e., the travel of the tool center.

[0013] In the simplest case, a surface grid can be obtained by linearly interpolating the milling points. For example, conventional grid algorithms compute a triangular grid. The algorithm includes the interpolation as well as other characteristic properties, such as the accuracy of the milling path and the contour. Disadvantageously, this method uses only the milling points which typically describe the tool center points (TCP) to calculate the surface grid. With a cylindrical cutter, the TCP—which is sometimes also referred to as tool tip—represent those points where the longitudinal axis of the tool penetrates the workpiece surface. The tool center point is calculated as the point displaced by the tool radius correction (WRK). For spherical cutters, the TCP and the center of the sphere are typically identical.

[0014] This models a point-shaped cutter (cutter with a dimension zero). It will be understood, that this model does reflect a realistic milling processes. For example, if the cutter has a diameter of 12 mm, then the location of the actual milling cut is displaced by 6 mm from the position of the cutter center point.

[0015] Disadvantageously, the geometry of a cutter is not taken into account in the simulation. To obtain realistic workpiece contours, the geometry of the cutter should therefore be considered. Conventional methods that take the cutter geometry into account correspond more closely to the original, but cannot compete due to delay caused by their internal data structures and the long computing times.

[0016] The known algorithms for solving this problem are therefore either imprecise in their approximation of the contour or inadequate in their performance.

[0017] It would therefore be desirable to provide a method based on the initial milling data, which enables a sufficiently precise simulation of the expected workpiece contour with acceptable computing times. In addition, the cutter geometry should be taken into account without unduly impairing the computing time.

[0018] According to an aspect of the present invention, a method for simulating a workpiece contour by simulating a milling operation is provided, wherein a blank workpiece is modeled by a family of straight lines in three-dimensional space. The straight lines are all parallel to an axis, for example, the z-axis. The workpiece contour is determined by cutting the family of straight lines with a virtual milling edge following a predetermined milling path in the three-dimensional space based on the intersecting points of the straight lines of all swept volumes. A surface grid that describes the workpiece contour is formed from the obtained set of intersecting points.

[0019] For simulating a milling operation with three axes, only the intersecting point with the smallest coordinate in the direction of the straight line is considered for forming the surface grid. If a straight line has more than one intersecting point, then only that point is used that has the smallest value in the direction of the dimension parallel to the family of straight lines.

[0020] If in addition to the straight lines with intersection points, straight lines without an intersection point are also to be considered, then an initialization value is advantageously used for forming the surface grid.

[0021] Advantageously, a modeled blank workpiece can be subdivided into discrete sections in the direction of the other two dimensions (x, y), wherein the straight lines of the family of straight lines can be equidistantly arranged in the direction of the other two dimensions (x, y). The spacing between the straight lines can be changed interactively.

[0022] The computation can be more efficient if for modeling the blank workpiece, the family of straight lines is bounded by a rectangle that is parallel to the family of straight lines, with the edge length and spatial position of the rectangle determined by using the extreme coordinates of the predetermined milling path in the three-dimensional space, wherein all calculations are performed within the space defined by the rectangle.

[0023] The radius of a cutter can be taken into account by enlarging the determined edge lengths of the bounding axis-parallel rectangle in the positive and negative direction by the radius of the actual milling head.

[0024] Advantageously, the milling path can be represented by a path of the center of the cutter, which is determined as a linearly interpolated polygon path in the three-dimensional space from predetermined milling data, in particular from the parts program or from an interpolator.

[0025] According to another advantageous embodiment, a milling operation can be simulated with a cutter by simulating a volume of the family of straight lines that is swept during the machining operation from geometric partial volumes and by determining the intersecting points with the family of straight lines separately for each of these partial volumes.

[0026] Advantageously, a milling processing can be simulated with a spherical cutter by simulating a volume, that is swept during the machining operation of the family of straight lines, from partial volumes in the form of two spheres or hemispheres, which are connected by an additional partial volume of a cylinder aligned in the milling direction, wherein the spheres or the hemispheres as well as the connecting cylinder each have the cutter radius, and by determining the intersecting points separately with this family of straight lines for each of these partial volumes.

[0027] The computing time can be further optimized by not determining the intersecting points of the straight lines with the partial volumes of the cylinders if the milling path is oriented vertically.

[0028] Advantageously, only the volume underside of the swept volume of the family of straight lines needs to be taken into account for simulating a milling operation with a cylindrical cutter.

[0029] Preferably, a volume underside of the cylindrical cutter can be simulated as a circular disk and a displacement of this volume underside of the cylindrical cutter can be simulated by an elliptical cylinder that is bounded by two planes, wherein the circular disks and the connecting elliptical cylinders each have the cutter radius and wherein intersecting points with the family of straight lines are determined separately for each of these partial volumes.

[0030] A further improvement of the computing time can be achieved by locally restricting a set of the straight lines to be considered to those lines that are bounded by the axis-parallel rectangle, in that separate local axis-parallel sub-rectangles (Bn) are formed in this axis-parallel rectangle for corresponding partial sections of the milling path, with the dimensions of the sub-rectangles being limited to the corresponding swept cutter volume, and by taking into account only the straight lines located inside a sub-rectangle when determining relevant intersecting points.

[0031] When simulating a milling operation with a spherical cutter, the set of the straight lines to be considered can be locally limited to the family of straight lines located inside the axis-parallel rectangle by forming for each corresponding partial volume sphere separate local axis-parallel sub-rectangles inside this axis-parallel rectangle, hemisphere or cylinder. The dimensions of the sub-rectangles can be limited to the respective partial volume, with only the straight lines located inside a sub-rectangle being considered when determining relevant intersecting points.

[0032] The running or delay time can be further optimized when simulating a milling operation with a spherical cutter by determining the intersecting points with the family of straight lines for overlapping partial volumes only once, if the partial sections on the milling path are contiguous.

[0033] According to another aspect of the invention, an apparatus, for example, a suitably programmed simulation computer is provided, which includes application software for carrying out the method steps when the software is executed on the computer.

[0034] According to yet another aspect of the invention, a computer program product with software segments, which can be loaded directly into the internal memory of a digital computer for execution on the digital computer, is provided for carrying out the aforedescribed method steps.

[0035] In addition, a surface grid, such as a triangular grid, can be computed from the calculated contour points which can be subsequently supplied to a visualization component, such as a computer with a corresponding graphics display. A method for generating such surface grid from contour points is described, for example, in Maier, K.-H.: “Ein auf uniform mäandrisch vorliegenden 3D-Punktdaten basierender Vernetzungsalgorithmus”, [Grid algorithm based on uniform meander-shaped 3D point data] Technischer Bericht, Siemens A G, Nürnberg, 2001.

[0050] Throughout all the Figures, same or corresponding elements are generally indicated by same reference numerals.

[0051] The basis of the method of the invention are straight lines that are equally spaced in the x- and y-direction and are parallel to the vertical z-axis. Preferably, only those straight lines are permitted whose point-shaped orthogonal projections in the xy-plane are located inside an axis-parallel rectangle, which will be determined later. The family of straight lines hereby models a blank workpiece with discrete coordinate values in the x- and y-direction and extended to infinity in the positive and negative z-direction. The representation according to FIG. 1 is a two-dimensional representation of such uniform family of straight lines with a surrounding rectangle B in a top view parallel to the straight lines for an exemplary milling path. The observer looks along the z-direction. The straight lines are indicated by circles or crosses. A spherical cutter F is indicated by a large circle.

[0052] In addition to the family of straight line, the milling points FP, which are indicated by circles that are greater than those of the straight lines and smaller than the circle of the cutter F, as well as a resulting milling path FB are positioned as a linearly interpolated polygon in three-dimensional space. It should be noted that there need not be a spatial correspondence between the rectangle defining the family of straight lines and the milling points. However, the rectangle B is typically determined using the milling data.

[0053] In a realistic milling process, the cutter F moves according to the programmed cutter center path (milling path FB). The cutter hereby removes material from the blank workpiece. Only the workpiece to be fabricated remains as an end result of the machining process. The approach taken by the method of the invention is similar to this process. The virtual cutter F moves here likewise according to the corresponding milling path FB and cuts the blank workpiece which, as mentioned above, is modeled as a family of straight lines. The swept volume can be simulated by two spheres and a cylinder extending in the milling direction—shaped like a pill.

[0054] This is illustrated in FIG. 2 which is similar to FIG. 1. At this time, only the region of FIG. 2 along the milling edge FK will be of interest. FP1 indicates the starting point and FP2 the end point of a milling path. The milling edge is the line connecting FP1 with FP2. The volume removed by the cutter F at the starting point FP1 is indicated by a circle, likewise at the end point FP2. The intermediate swept volume along the line from FP1 to FP2 is indicated by two lines tangential to the periphery of both circles about FP1 to FP2. When the total area formed by these geometries is viewed in three-dimensional space, the resulting swept volume has the shape of a pill which is formed by two spheres and a cylinder extending in the milling direction.

[0055] The volume (pill) swept by the cutter F relative to a milling edge FK is intersected by the straight lines belonging to the family of straight lines either never, or once, or twice. The workpiece contour is discretely simulated for milling in three axes by the lower intersection points of all swept volumes (pills) with all straight lines and the initialization points of all straight lines that are not intersected. The volume swept by the cutter F relative to an edge produces intersecting points with the straight lines. FIG. 1 shows those swept straight lines that have intersecting points as points, such without intersecting point as crosses.

[0056] These intersecting points are candidates for the workpiece contour to be simulated. When machining in three axes, the intersecting points are part of the final workpiece contour exactly when no other intersecting points with a smaller z-coordinate exist.

[0057] For determining the intersecting points—which will be discussed later in more detail—the diameter of the employed tool (spherical cutter) can be defined. Another parameter <maxerror> enables an interactive control of the fineness of the discrete values within the rectangle B bounding the family of straight lines, i.e., the resolution of the computed workpiece contour. An exemplary unit for measuring dimensions can be a millimeter. The number of quantization steps mxsize in the x-direction and mysize in the y-direction are computed according to

• [0058] mxsize=ROUND(Rx/maxerror)
• [0059] mysize=ROUND(Ry/maxerror)
  wherein Rx and Ry represents a lengths of the sides of the surrounding rectangle R (so-called BoundingBox) and ROUND the rounding operator.

[0060] Another parameter <mr/or> can be used to determine the rectangle B bounding the family of straight lines. The user can select between 1 for mr (with radius) and 0 for or (without radius). If or is selected, i.e., if the dimension of the cutter F is neglected, then the rectangle B is defined as BoundingBox (minimal axis-parallel cuboid superset) of the xy-plane projection of the read-out 3-D data. Conversely, if mr is selected, then the radius is added to the x- and y-dimensions of the or-BoundingBox, thereby either increasing or decreasing the dimensions. It is recommended to select mr for workpieces that are substantially oriented in the z-direction, since the area of the BoundingBox is determined by a projection and can hence become infinitesimally small, i.e., the segment to be computed can shrink to a limit of zero.

[0061] The number of data points to be outputted can be defined by another parameter <wholematrix>. The user can here again select between 1 and 0. The first case, all computed z-intersections within the bounding rectangle and in addition the initializations of straight lines that are not intersected are outputted. Otherwise, only those computed data are outputted that are located on the workpiece contour. This is actually only a subset of all data, because the 3-D data describing the workpiece are generally not parallel to the axes, i.e., their xy-plane projection does not match the surrounding BoundingBox rectangle B.

[0062] The factual situation can be described as follows, in analogy to the representation of FIG. 1: the selection <wholematrix>=1 results in the output of the values for all straight lines of the family of straight lines, i.e., crosses and points, whereas with <wholematrix>=0 the output of crosses (straight lines without intersection) is suppressed.

[0063] The additional operation of an application according to the method of the present invention is indicated schematically in steps 1 to 9 of the flow diagram depicted in FIG. 3.

[0064] Step 1 only checks that the consistency of the simulation and that the simulation calls the proper parameters. In the second step 2, the input data are read, i.e., the milling data, for example based on an NC parts program. Step 3 relates to the calculation of the BoundingBox B. As described above, the method of the invention operates preferably with the help of a vertical family of straight lines bonded by a rectangle B, which is a representation of the blank workpiece by discrete values (see FIG. 1). This rectangle B, i.e., the two-dimensional BoundingBox, is determined as a function of the 3-D data points (cutter positions) and the radius of the (spherical) cutter F.

[0065] First, it will be described how the BoundingBox is determined, for example in the event that the user selects <mr/or>=0, i.e., without taking into account the radius of the cutter. For this purpose, the extreme coordinates xmin, xmax, ymin, and ymax of the data points set are determined. The edge length of the bounding rectangle B is then determined by

• • [0066] xmax−xmin and/or ymax−ymin.

[0067] The spatial location of the rectangle B will have to be determined and fixed. Since the rectangle is located in the xy-plane, it is sufficient to specify its corner points. For example, the corner point of the rectangle having the smallest x- and/or y-coordinates can be selected.

[0068] If the BoundingBox is to be changed later in the z-direction, for example, for adapting the blank workpiece, a z-component can also be included in the data structure and defined, like the x- and y-components, by setting the minimum z-coordinate and the lines of the 3-D cuboid and the z-direction zmax−zmin. A safety margin can be added to this value which is also corrected by the dimension of the radius in the negative z-direction, so as to eliminate rounding errors in the calculation.

[0069] The BoundingBox B bounding the family of straight lines is defined—in the case where <mr/or>=0 is selected—by the values

• [0070] B.x=xmin, B.wx=xmax−xmin
• [0071] B.y=ymin, B.wy=ymax−ymin
• [0072] B.z=zmin−2*radius, B.wz=zmax−zmin+10*radius
  wherein B.x, B.y and B.z denote the corners of the cuboid having the minimum coordinates and B.wx, B.wy and B.wz denote the length of the corresponding cuboid sides.

[0073] The BoundingBox B is actually present in the form of a 3-D cuboid which is also parallel to the coordinate axes. It should be mentioned again that for carrying out the method of the invention, the z-coordinate plays typically no role and needs to be kept only if the user later wishes to redefine the BoundingBox.

[0074] If the user desires to take into account the cutter radius (<mr/or>=1), then the x- and y-coordinates of the BoundingBox B have to be adapted. This is done by increasing the length of the BoundingBox edges in the positive and negative direction by the radius, while the x-coordinates remain unchanged. Accordingly, for <mr/or>=1, the BoundingBox B is defined by

• • [0075] B.x=xmin−radius, B.wx=xmax−xmin+2*radius
  • [0076] B.y=ymin−radius, B.wy=ymax−ymin+2*radius
  • [0077] B.z=zmin−2*radius, B.wz=zmax−zmin+10*radius

[0078] Following the calculation of the BoundingBox B, the straight lines of the family of straight lines included in the computation of the z-intersections are determined. This is done with the help of the variables mxsize and mysize, which are defined by

• • [0079] mxsize=ROUND(B.wx/maxerror)
  • [0080] mysize=ROUND(B.wy/maxerror).

[0081] In other words, mxsize and mysize define the fineness of the discrete values of the blank workpiece in the x- and y-direction as a function of the computed BoundingBox and the resolution selected by the user.

[0082] The process steps 4 to 7 are optional and enable an interaction for a user-specific change of the BoundingBox. The determined parameters can be outputted for testing purposes in step 4. Steps 5, 6 and 7 make it possible to set a specially selected BoundingBox, wherein step 5 checks if the change has been requested. If a change has been requested in step 5, then the dimensions of the new BoundingBox B are inputted in step 6 and the new parameters are outputted in step 7—exactly like in step 4. Subsequent to the user-defined change of the BoundingBox, the variables mxsize and mysize are updated. Conversely, if no change has been requested in step 5, then the process jumps directly to step 8.

[0083] Step 8 is directed to the actual computation of the z-intersections, i.e., of the intersecting points with the straight lines of the family of straight lines. As described above, the workpiece contour is a result of discrete cutting of the parallel vertical straight lines bonded by the BoundingBox B. This is done by moving the cutter F (in the exemplary embodiment a sphere) according to the programmed cutter center point path FB and intersecting with the straight lines. When looking at an individual edge FK, i.e., the connecting line between consecutive milling points FP1 and FP2, then a volume is produced by the spatial movement of the spherical cutter. As mentioned above, the volume is advantageously simulated by two spheres and an enclosed cylinder. The spheres are hereby located at the two corresponding milling points, the cylinder is oriented in the direction of the milling edge, so that the swept volume resembles that of a pill.

[0084] Alternatively, the volume swept by the spherical cutter can also be modeled by two hemispheres and an interior cylinder. However, full spheres are used, since implicit functions are preferred over parametrically defined functions for calculating the intersections. Moreover, with hemispheres, a change in the orientation would require the milling edge to be recomputed which would negatively affect the computing time.

[0085] If the intersections of all straight lines to be used in the computation relative to this milling edge FK are considered, then one can have two intersecting points, a single intersecting point or no intersecting point at all. The intersecting point (so-called z-intersection with a smaller z-coordinate than the current actual z-intersection for the particular straight line is then selected as a new workpiece contour point. In the case of two intersecting points, only the intersecting points with the smaller z-coordinate is considered.

[0086] Finally, in step 9, the method of the invention produces the discretely modeled workpiece contour by looking at all intersections of all edge volumes (pills) with all corresponding straight lines.

[0087] Step 8 is formally implemented with the functions InitializeZcoordinates, UpdateZcoordinates and EvalZcoordinates and/or intcyl and SphereLineIntersection. The workflow follows the process depicted in FIG. 4, wherein the area with the gray background is executed sequentially for all milling edges.

[0088] First, a two-dimensional field (matrix) with the dimension (mxsize x mysize) is created using the function InitializeZcoordinates and initialized with the value of the largest z-coordinate of the calculated or selected BoundingBox, i.e., with the value B.z×B.wz. The function UpdateZcoordinates processes all straight lines with respect to a milling edge, i.e., with respect to a single “pill.” The function UpdateZcoordinates calls the function EvalZcoordinates, which takes on the calculation of its intersecting points (for example, a test for an intersecting point) with the two spheres and the interior cylinder by using the additional functions intcyl and SphereLineIntersection. The function intcyl can be implemented very efficiently since the particular straight lines extend in the vertical direction.

[0089] The results of this step (intersection of all straight lines with a pill) represents the workpiece contour, which is created by the corresponding milling edge FK, i.e., by milling from one point FP1 to the following point FP2. The coordinates of the respective current straight line aktLine are then given by the variables

• • [0090] xstep=B.wx/mxsize
  • [0091] ystep=B.wy/mysize
    according to the following procedure:
  • [0092] aktLine_x=B.x+i*xstep, i=0, . . . , mxsize,
  • [0093] aktLine_y=B.y+j*ystep, j=0, . . . , mysize.

[0094] This ensures that all the dimensions of the BoundingBox B with respect to the x- and y-coordinates are taken into account. The discretely computed workpiece contour is produced by processing all milling edges, i.e., by processing all pills.

[0095] The function SphereLineIntersection computes the intersection of the straight line and the sphere, which will be described in more detail below. Initially, the sign of the residue

• • [0096] res=radius*radius−(x−m.x)*(x−m.x)−(y−m.y)*(y−m.y)
    is checked, whereby radius describes the radius of the milling sphere, x and y the x- and y-coordinates, respectively, of the considered straight line and m.x and m.y the center of the sphere (milling point coordinates). If the sign is negative, then the straight line does not intersect the sphere. Otherwise, the z-intersection
  • [0097] *zvalue=m.z−sqrt(res)
    is returned by the function SphereLineIntersection.

[0098] For forming the intersection between the straight line and cylinder, which is much more difficult, the function “intersect a ray with a cylinder” described in Heckbert, P. S. (ed.): Graphics Gems IV (The Graphic Gems Series), Harcourt Publishers Ltd, a subsidiary of Harcourt International Ltd., 1994, ISBN 0-12-336156-7, can be used within the function intcyl. This function calculates the potential intersecting point of the actual straight line with the cylinder located in the direction of the milling path, whereby the cylinder is initially considered to extend to infinity. For testing if a obtained intersecting point is located within the cylinder portion bonded by the two milling points, a “normal test” is performed. If the result is positive, then the obtained intersecting point is taken as the new z-intersection, i.e., its z-coordinate is compared with the current z-coordinate, possibly replacing it.

[0099] The method which was briefly described above will now be explained in more detail. The method is illustrated in FIGS. 5 and 6 wherein the notation is slightly different from the previous notation. Hereby denotes:

[0100]

d, xg, λ	Scalars (e.g., length of a vector)
a, r	Vectors
Q, Φ	Matrices
BR, BC	Points, i.e., elements of an affine vector space
g1, h	One-dimensional geometric shapes (e.g., straight line)
E1, K	Two-dimension geometric shapes (e.g., planes, spheres)
|d |, |r |	(Euklidian) norm of a vector (in particular, absolute value of
	scalar).

[0101] FIGS. 5 and 6 show geometrically an intersection between a straight line and a cylinder, as explained more particularly by the following mathematical description.

[0102] The predefined circular cylinder is defined by the point B_c and the direction r (normalized), the straight line which in the example is vertical with the point B_R and the direction r (normalized).

[0103] Initially, it is checked with a “quick elimination test” if the objects cylinder and straight line intersect each other, i.e., if an intersecting point can be computed at all. For this purpose, according to

[0104] $d=\frac{r\times a}{r\times a},$
a vector d which is orthogonal to the directions r and a is computed, wherein x represents the cross (vector) product of the vectors. If the distance between the two inclined straight lines defined by r and a is greater than the radius of the cylinders, then the elimination test gives a negative result—i.e., the straight line and the cylinder do not intersect.

[0105] Otherwise, two intersecting points H_in and H_out (which do not necessarily have to be different) exist which are located symmetric to the point P (see FIG. 6 with H_in=H₁, H_out=H₂). These satisfy the equation
H_in=B_R+t_inr, H_out=B_R+t_outr,
whereby the parameters t_in and t_out are given by
t_in=t−s, t_out=t+s.

[0106] The following relationship can be defined for the parameters t and s:

[0107] $t=-\frac{\left\{\left(B_R-B_C\right)\times a\right\}·d}{r\times a}$ $\mathrm{and}$ $s=\frac{\sqrt{r^2-d^2}}{r·o},o=\frac{d\times a}{d\times a}.$

[0108] The second relationship is the result of the Pythagorean relation between k, d and r (see FIG. 6) and the observation
sr·o=±k.

[0109] If the straight line extends parallel to the cylinder, i.e., that d in the above form cannot be computed, then the “intersecting points” have to be derived differently than described above. However, this does not apply to the method of the invention and will therefore not be discussed further.

[0110] Finally, after the intersecting point with the smaller z-coordinate has been calculated, it has to be checked if this intersecting point is located inside the finite cylinder under consideration. The relevant cylinder section is defined by the two boundary planes that are located in the milling points and are orthogonal to the milling edge. This is illustrated schematically in FIG. 7, which shows how a normal test is performed for an interior intersecting point H_in.

[0111] The vector H_in−P1 is first scalar multiplied with the milling edge vector P2−P1. If the result is negative, then H_in−P1 does not point in the direction P2−P1, i.e., the intersecting point H_in is not located in the relevant cylinder section. Otherwise, the vector H_in−P2 can be multiplied scalar with the milling edge vector P1−P2 and the sign checked. If the sign is positive, then H_in−P2 points in the direction P1−P2, and the intersecting point H_in is located in the relevant cylinder section.

[0112] However, instead of comparing the last step, it is more efficient to compare the length of the projection of the vector Hin−P1 onto P2−P1 with the length of the milling edge vector P2−P1. If the projection is shorter (or has the same length), then the intersecting point H_in is located in the relevant cylinder section.

[0113] It will now be shown that the formula

[0114] $t=-\frac{\left\{\left(B_R-B_C\right)\times a\right\}·d}{r\times a}$
correctly describes the intersection between a straight line and a circular cylinder with reference to the diagram of FIG. 8.

[0115] For this purpose, the plane defined by B_c and the vectors a and d is intersected by the straight line which is defined by B_r and the direction vector r. The parameter value t which is to be calculated can subsequently be provided as the distance between the obtained intersecting point and the point B_r.

[0116] Given the plane which is defined by B_c and the direction vectors a and d. This plane has the form
E: X=B_C+μa+νd, μ,ν∈.

[0117] The straight line shall be defined by B_r and the direction vector r as follows.
g: X=B_R+λr, λ∈.

[0118] For intersecting the two geometric shapes, the two formulas are set equal, giving
B_R+λr=B_C+μa+νd
or equivalently
B_R−B_C=μa+νd−λr.

[0119] By forming on both sides of equations the cross product with a, one obtains
(B_R−B_C)×a=(μa+νd−λr)×a
or
(B_R−B_C)×a=(μa)×a+(νd)×a−(λr)×a.

[0120] Based on the definition of the cross product,
(μa)×a=μ(a×a)=0,
so that
(B_R−B_C)×a=(νd)×a−(λr)×a

[0121] After multiplying the two sides of the equation by d, one obtains
{(B_R−B_C)×a}·d={(νd)×a}d−{(λr)×a}·d.

[0122] Since d is parallel to the plane spanned by d and a (coplanar vectors), it follows
{(νd)×a}·d=ν{d×a}·d=0
and therefore
{(B_R−B_C)×a}·d=−{(λr)×a}·d=−λ(r×a)·d.

[0123] By dividing both sides of equations by |r×a|, one obtains

[0124] $-\frac{\left\{\left(B_R-B_C\right)\times a\right\}·d}{r\times a}=\lambda \frac{\left(r\times a\right)·d}{r\times a}=\lambda d·d=\lambda .$

[0125] The intersecting point P of the computed plane with the straight line is then obtained by inserting the calculated □-value into the equation for the straight line

[0126] $P=g\left(\lambda \right)=B_R-\frac{\left\{\left(B_R-B_C\right)\times a\right\}·d}{r\times a}r.$

[0127] Finally, the distance t between P and B_R is calculated as

[0128] $t=B_R-P=B_R-B_R+\frac{\left\{\left(B_R-B_C\right)\times a\right\}·d}{r\times a}r=\frac{\left\{\left(B_R-B_C\right)\times a\right\}·d}{r\times a}=-\frac{\left\{\left(B_R-B_C\right)\times a\right\}·d}{r\times a},$
whereby in the last conversion it has to be observed that B_R−B_C, a and d form a left-handed system.

[0129] The spacing between the inclined straight lines can be determined as follows. Given two inclined straight lines g₁ and g₂ with the origins B_R and B_C, respectively, and the direction vectors r and a, i.e.,
g₁: X=B_R+λr, g₂: X=B_C+λa, λ,μ∈.

[0130] There exists exactly one straight line h which orthogonally intersect the straight lines g₁ and g₂. The distance d between the two intersecting points represents the separation between the inclined straight lines g₁ and g₂.

[0131] Calculation of d can be realized by defining, for example, two planes E₁ and E₂, wherein E₁ includes the straight line g₁ and E₂ includes the straight line g₂. If, in addition, the direction vector of h is normal to the planes, then the two planes are uniquely defined (see diagram in FIG. 8).

[0132] Based on the definition for enclosed angles, it holds that (see diagram in FIG. 9)

[0133] $\cos \left(\beta \right)=\cos \left(\angle \left(B_R-B_C,d\right)\right)=\frac{\left(B_R-B_C\right)·d}{B_R-B_C·d},d=\frac{r\times a}{r\times a}.$

[0134] Using the relationship for right triangles,

[0135] $\cos \left(\beta \right)=\frac{d}{B_R-B_C},$
one obtains (d is normalized)

[0136] $\frac{d}{B_R-B_C}=\frac{\left(B_R-B_C\right)·d}{B_R-B_C·d}=\frac{\left(B_R-B_C\right)·d}{B_R-B_C}$
and therefore
d=|(B_R−B_C)·d|.

[0137] An intersection between a straight line and a sphere can also be easily determined. Given the straight line g by

[0138] $g\text{:}X=\left(\begin{array}{c}{x}_g\\ y_g\\ 0\end{array}\right)+t\left(\begin{array}{c}0\\ 0\\ 1\end{array}\right),t\in ℛ$
and the sphere K by
K: (x−x_m)²+(y−y_m)²+(z−z_m)²=r².

[0139] The intersection of the two geometric shapes is formed by inserting the equation for the straight line into the equation for the sphere. One obtains
(x_g−x_m)²+(y_g−y_m)²+(t−z_m)²=r²,
or equivalently
t=z_m±√{square root over (r²−(x_g−x_m)²−(y_g−y_m)²)}{square root over (r²−(x_g−x_m)²−(y_g−y_m)²)}.

[0140] The process flow of the method of the invention described so far, whereby individual milling volumes (pills) are tested for forming an intersection with all available straight lines, and hence the lead-time and/or delay time can be improved considerably. In particular, when the dimensions of the cutter F are small in comparison to the calculated or selected BoundingBox B, a number of straight lines exist which give a negative result when testing for an intersection with a currently active pill. The method can therefore operate significantly more efficiently if the set of the straight lines to be tested is drastically restricted.

[0141] This can be implemented with local BoundingBoxes Bn and hence with local z-intersection updates. This advantageous method will hereinafter be described in more detail with reference to the already discussed FIG. 2.

[0142] Not all straight lines have to be taken into account when calculating the intersecting points between a straight line and the volume (pill) swept by the cutter F relative to the actual milling edge FK. It is sufficient to consider only those straight lines having xy-coordinates located inside the xy-projection of the pill Bounding Box.

[0143] These straight lines which are indicated in FIG. 2 by points, have necessarily the xy-coordinates located in the interval [lx, ux]×[ly, uy], wherein the indicated variables are defined by

• • [0144] lx=MIN(p1.x,p2.x)−radius
  • [0145] ux=MAX(p1.x,p2.x)+radius
  • [0146] ly=MIN(p1.y,p2.y)−radius
  • [0147] uy=MAX(p1.y,p2.y)+radius.

[0148] The indices required for selecting the individual straight lines can then be calculated according to

• • [0149] li=(int)ROUND((lx−B.x)/B.wx*mxsize)
  • [0150] ui=(int)ROUND((ux−B.x)/B.wx*mxsize)
  • [0151] lj=(int)ROUND((ly−B.y)/B.wy*mysize)
  • [0152] uj=(int)ROUND((uy−B.y)/B.wy*mysize).

[0153] As seen in FIG. 2, since the BoundingBox of the pill volumes can exceed the BoundingBox of the straight line, the algorithm only accesses the x- and y-values in the interval

• • [0154] [MAX(0,li), MIN(mxsize,ui)] and
  • [0155] [MAX(0,lj), MIN(mysize,uj)]

[0156] Alternatively, separate BoundingBoxes Bn can be calculated for the individual partial volumes cylinder and sphere. However, this is unnecessary for large dimensions of the cutter and small distances between milling points.

[0157] In the following, a different approach for calculating cylinder intersections will be described. The algorithm “intersect a ray with a cylinder” used for calculating the intersecting point between a straight line and a cylinder is described in Heckbert, P. S., referenced above for intersections of cylinders in any spatial arrangement with straight lines in any spatial arrangement.

[0158] However, the method of the invention presents a special situation in that the straight lines are always parallel to the z-axis, although the cylinders can be arranged arbitrarily in space. The source variables raybase and raycos of the function intcyl described therein can then be simplified, reducing the delay times of the method of the invention due to a faster computation of the scalar products.

[0159] The delay time can be further improved if an intersection between the cylinder and the straight line is not computed for a vertical milling edge. This is true because the minimal z-intersection of straight lines remains unchanged (if the straight line does not intersect a cylinder) or the intersecting point is already derived from the initially tested spheres located at the starting point and end point of the milling edge. If it is known which of the two spheres is lower in the z-direction, then it is sufficient to test only this sphere.

[0160] In the following, another modified embodiment for forming the intersecting volumes will be described. As already mentioned above, the z-intersections are computed relative to the edges. The volume swept relative to each milling edge is herein tested with the available straight lines. This test includes for each milling edge, in addition to the cylinder intersection test already described in detail, forming the intersection between the straight line with the spheres at the starting point FP1 and end point FP2 of the milling edge FK. This method implies, however, that for treating the subsequent milling edge, a redundant sphere-straight line-test is performed, since the sphere at the end point of the previously treated edge is identical to the starting point sphere of the following edge. If this swept volume is viewed as a sphere-cylinder-pair, then the intersecting points for approximately 50% of the spheres need not be computed.

[0161] Aside from the exemplary case of a spherical cutter, other types of cutters can also be simulated with the method of the invention. Practical milling processes can employ, for example, end mill cutters (with and without rounded corners), angular cutters (with/without rounded corners), surface cutters, disk cutters and tapered cutters.

[0162] In the following, a second embodiment using a cylindrical cutter will be described. Initially, the volume swept by the cutter-cylinder during the milling process (referred to as so-called sweep volume) is described analytically and intersected with a vertical family of straight lines. In this case, not the entire milling path, but only sections are considered (e.g., in a milling path represented by an NC parts program with several data sets). The volume swept by the cutter relative to a linear set (milling edge) is then formally simulated. It is sufficient to consider the volume underside, since the z-axis represents the tool travel direction in a three-axes milling process.

[0163] Accordingly, the intersecting points between the straight lines and the underside of the cutter cylinder are first calculated at the points P1 and P2 (circular disks, see FIG. 10). Displacement of the underside of the cutter cylinder in P1 in the direction P1P2, typically not in the horizontal direction, produces an elliptical cylinder which also has to be intersected with the straight lines.

[0164] In this case, the corresponding lower intersecting points H_u have to be calculated, whereby only those intersecting points are relevant which are spatially located between the two planes through the points P1 and P2 that are parallel to the xy-plane (a (finite) elliptical cylinder bonded by parallel planes). The intersection between the straight line and the circular disks should not be neglected. If the elliptical cylinder is intersected both above and below a circular disk, then the lower intersecting point is not in the relevant region, whereas the z-coordinate of the upper intersecting point is too large. In this case, the intersecting point between the straight line and the circular disk would have to be processed exactly.

[0165] Forming the intersection between straight lines and circular disks is non-critical, as discussed in detail above. It only has to be determined if the starting point of the straight line is inside the circular disks projected onto the xy-plane. If this is the case, then the z-coordinate of the intersecting point is the z-coordinate of the point P1 or P2, depending which circular disk is considered.

[0166] Forming the intersection between the straight lines and the elliptical cylinder is more complex and numerically less stable. It is hereby not sufficient to consider only the location of the points—the actual coordinates of the intersecting points have to be computed. This part of the problem can be solved by representing the marginal surfaces of the elliptical cylinders as a Null set of an implicit function, or more specifically as a quadrics. Quadrics (quadratic form) represents a vector set having elements x with the functional relationship
x^TAx+2x^Tb+c=0
with a symmetric matrix A and vectors b and c. The intersecting points are obtained by inserting the parameterized straight lines as well as by solving the resulting systems of second-order equations.

[0167] In principle, the intersecting points can be calculated in two different ways: either the calculations are performed in a standard coordinate system or the intersecting points are computed in a new coordinate system that is independent of the standard coordinate system.

[0168] If the first case (solution trial 1), i.e., when calculating in the standard coordinate system, the elliptical cylinder and the straight line have to be represented in the standard coordinate system (i.e., standard basis, origin at 0, the coordinate system (O, e1, e2, e3) in FIG. 11 wherein similar to the representation of FIG. 10, the elliptical cylinder is shown with the lower intersecting point H_u and the upper intersecting point H_o). The quadrics of the elliptical cylinder however, is not known a priori and needs first to be computed.

[0169] In the other case (solution trial 2), the elliptical cylinder and the straight line are represented in the new coordinate system and the intersecting points are computed in this coordinate system. The coordinate system is preferably selected so that the base of the coordinate system with the normalized major axes matches the ellipse generating the elliptical cylinder as well as a normalized orthogonal vector (in the direction P1P2). As a result, the elliptical cylinder in this coordinate system now has an normal form, i.e., mixed terms are absent in its implicit representation (“major axes transformation”, see in FIG. 11 the coordinate system [P1, {circumflex over (b)}, â, {circumflex over (n)}], whereby {circumflex over (b)}, â and {circumflex over (n)} refer to the vectors normalized with respect to b, a and n). The order of the vectors is arbitrary. In the following, the shorter major axes vector {circumflex over (b)} represents the first base vector of the new coordinate system, with the point P1 selected as the origin of the coordinate system. The computed intersecting points are therefore not represented in the standard coordinate system and have to be transformed back in a subsequent step by an inverse coordinate transformation (change of the coordinate system).

[0170] In the following, the two approaches for forming the intersection between straight lines and elliptical cylinders are once more schematically described. The relationship with the coordinate system is indicated by superscript numbers in parentheses: for example, z⁽³⁾ indicates the z-coordinate of a vector represented in the standard coordinate system.

[0171] Solution Trial 1 (Computing an Intersection in a Standard Coordinate System):

• [0172] a1) computing the formally implicit representation of the elliptical cylinder in the standard coordinate system;
• [0173] b1) intersecting the straight line (represented in the standard coordinate system) and the elliptical cylinder, i.e., inserting the parameter representation of the straight lines in the quadrics equations;
• [0174] c1) solving the quadratic equation for the parameter t=z^(S) of the straight line; and
• [0175] d1) computing the intersecting point having the smallest z^(S)-coordinate by inserting the parameter values computed in c).

[0176] Solution Trial 2 (Computing an Intersection in the Coordinate System [P1, {circumflex over (b)}, â, {circumflex over (n)}]):

• [0177] a2) computing the representation of straight line and the coordinate system [P1, {circumflex over (b)}, â, {circumflex over (n)}];
• [0178] b2) intersecting the straight line and the elliptical cylinder (each represented in the coordinate system [P1, {circumflex over (b)}, â, {circumflex over (n)}]), i.e., inserting the parameter representation of the straight lines in the quadrics equations;
• [0179] c2) solving the quadratic equation for the parameter;
• [0180] d2) computing the two (not necessarily different) intersecting points by inserting the parameter values computed in c2) into the equations for the straight line; and
• [0181] e2) representing the intersecting points computed in d2) in the standard coordinate system (inverse transformation) and selecting the point having the minimal z^(S)-coordinate.

[0182] In the following, the mathematical solutions for solving the aforedescribed tasks will be summarized. It is assumed that a vertically extending straight line is to be intersected with the lower sweep volume relative to a single milling edge, i.e., with circular disks and elliptical cylinders. In the aforedescribed approach, the process can be optimized by taking into account only the intersection with the cutter cylinder in the point P1 (starting point of the milling edge) as well as with the sweep volume, without taking into account the intersection with the cutter cylinder in the point P2. This is possible because in the treatment of the subsequent milling edge, the cutter cylinder at the end point of the previous milling edge can be regarded as being identical to the cutter cylinder at the starting point of the following milling edge.

[0183] Assuming, that the vertical straight line has the form

[0184] $g\text{:}X=\left(\begin{array}{c}{x}_g\\ y_g\\ 0\end{array}\right)+\lambda \left(\begin{array}{c}0\\ 0\\ 1\end{array}\right),\lambda \in 3,$
with the radius of the cutter cylinder denoted with r.

[0185] According to the aforedescribed approach, it is initially checked if the straight line intersects the circular disk defined by the underside of the cutter cylinder in the point P1. Thereafter, the lower intersecting point of the straight line with the elliptical cylinder has to be computed either according to solution trial 1 or solution trial 2. Those straight lines that do not intersect with the elliptical cylinder (as described above with reference to the spherical cutter) can be ready eliminated using the aforedescribed “quick elimination test.”

[0186] Intersection of Straight Line and Circular Disk

[0187] If the straight line intersects the cutter cylinder located at point P1, then this indicates that the straight line intersects with the circular-disk-shaped underside of the cutter cylinder, since the straight line is vertical. The resulting z-intersection coordinate then matches the z-coordinate of the point P1. Finally, it needs to be checked only if the starting point of the straight line is located inside the projection of the circular disk onto the xy-plane. The test condition can be written formally as follows
√{square root over ((x_g−x_m)²+(y_g−y_m)²)}{square root over ((x_g−x_m)²+(y_g−y_m)²)}≦r,
wherein x_m and y_m represent the x- and y-coordinate of the center P1 of the circular disk.

[0188] Intersection of Straight Line and Elliptical Cylinder

[0189] When intersecting a straight line and an elliptical cylinder, the situations for horizontal and non-horizontal milling edges have to be treated differently. For horizontal edges, the shorter half-axis b of the generating ellipse (see FIG. 11) has a length of zero, so that the ellipse degenerates into a line and the elliptical cylinder degenerates into a rectangle. Initially, the case of a non-horizontal milling edge will be described.

[0190] Approach with Non-horizontal Milling Edge

[0191] Computation of the new coordinate system [P1, {circumflex over (b)}, â, {circumflex over (n)}]:

[0192] The coordinate system to be calculated has its origin at P1. In addition, it will be assumed that the base vectors of the coordinate system point in the direction of the major axes of the ellipse generating the cylinder and/or in the direction of the milling edge extending perpendicular thereto. By normalizing the milling edge vector P1P2, the base vector {circumflex over (n)} of the coordinate system to be determined is already found. By also computing and normalizing the longer major axis vector a or the shorter major axis vector b of the ellipse, an additional base vector of the system is obtained. The straight base vector is obtained by forming the cross product of the above calculated vectors. Since the major axes vectors of the generating ellipse are not known a priori, at least one of these vectors has to be initially computed. The variant (i) described hereinafter computes the longer major axis vector a and therefrom the shorter major axis vector by forming the cross product with P1P2. The variant (ii) performs this operation in the reverse order.

[0193] Computation of the Major Axes of the Cylinder—Variant (i)

[0194] When using the variant (i), the longer major axis base vector a is computed first. For this purpose, the projection of the milling edge P1P2 is rotated in the xy-plane by 90° in the mathematically negative direction. This already provides the direction a_richt of the longer major axis vector (the absolute value of a is equal to the radios r of the cutter cylinder). The coordinate system base vector â is obtained by normalizing a_richt. The direction of the shorter major axis of the ellipse can be computed from the two computed vectors P1P2 and a_richt by forming the cross product
b_richt=a_richt×P1P2,
wherein one notices that the z-coordinate of a_richt becomes zero. The base vector {circumflex over (b)} is highly obtained by normalizing b_richt.

[0195] Computation of the Major Axes of the Cylinder—Variant (ii)

[0196] Unlike variant (i), variant (ii) first computes the shorter major axis vector b. It is hereby assumed that P1P2 is not horizontal, i.e., the vector b is not the null vector. The following figure illustrates the situation.

[0197] The major axis vector b is computed with the two auxiliary points S and R, which has to be computed first. The point S is assumed to be located horizontally from the point P1 in the direction of the xy-plane projection of P1P2 (see FIG. 12). S can hence be computed formally according to

[0198] $S=\mathrm{P1}+\frac{r}{p}p,$
wherein p is the projection of P1P2 onto the xy-plane. The point R, on the other hand, is obtained by projecting the point S onto the plane through the point P1, to which P1P2 is orthogonal. The formula for R is therefore
R=S+[(P1−S)·{circumflex over (n)}]·{circumflex over (n)}.

[0199] The major axis vector b is then
b=R−P1,
with a coordinate system base vector {circumflex over (b)} obtained by normalizing b.

[0200] The longer major axis vector a of the ellipse can then be computed by forming the cross product according to
a=P1P2×b,
which yields the coordinate system base vector â by normalizing a.

[0201] Due to the significantly longer derivation, computing the new coordinate system according to variant (ii) can be expected to be less efficient.

[0202] Determining the Quadrics Representation of the Elliptical Cylinder

[0203] Four-dimensional coordinates (homogeneous coordinates) can be used to represent the set of the boundary points of the elliptical cylinder, i.e., its quadrics. In general, a quadrics has the following representation

[0204] $f\left(x,y,z\right)=\left(\begin{array}{cccc}x& y& z& 1\end{array}\right)·\left(\begin{array}{cccc}a& b& c& d\\ b& e& f& g\\ c& f& g& h\\ d& g& h& i\end{array}\right)·\left(\begin{array}{c}x\\ y\\ z\\ 1\end{array}\right)=X^T\mathrm{AX}=0$
with a symmetric matrix A.

[0205] With the aforedescribed approach 1, i.e., the computation of intersection in a standard coordinate system, the above equation can advantageously be represented first in the coordinate system of the elliptical cylinder (see FIG. 11) as

[0206] $\left(\begin{array}{cccc}{x}^{\left(\mathrm{EZ}\right)}& y^{\left(\mathrm{EZ}\right)}& z^{\left(\mathrm{EZ}\right)}& 1\end{array}\right)·\left(\begin{array}{cccc}{q}_{11}& q_{12}& q_{13}& q_{14}\\ q_{12}& q_{22}& q_{23}& q_{24}\\ q_{13}& q_{23}& q_{33}& q_{34}\\ q_{14}& q_{24}& q_{34}& q_{44}\end{array}\right)·\left(\begin{array}{c}{x}^{\left(\mathrm{EZ}\right)}\\ y^{\left(\mathrm{EZ}\right)}\\ z^{\left(\mathrm{EZ}\right)}\\ 1\end{array}\right)={\left(X^{\left(\mathrm{EZ}\right)}\right)}^T{\mathrm{QX}}^{\left(\mathrm{EZ}\right)}=0,$
whereby the superscripts (EZ) represents the association with the coordinate system of the elliptical cylinder.

[0207] After construction of the coordinate system [P1, {circumflex over (b)}, â, {circumflex over (n)}], the elliptical cylinder coordinate system has a standard form, i.e.,

[0208] $Q=\left(\begin{array}{cccc}{q}_{11}& q_{12}& q_{13}& q_{14}\\ q_{12}& q_{22}& q_{23}& q_{24}\\ q_{13}& q_{23}& q_{33}& q_{34}\\ q_{14}& q_{24}& q_{34}& q_{44}\end{array}\right)=\left(\begin{array}{cccc}\frac{1}{{b}^2}& 0& 0& 0\\ 0& \frac{1}{{a}^2}& 0& 0\\ 0& 0& 0& 0\\ 0& 0& 0& -1\end{array}\right).$

[0209] As already mentioned above, the coordinate system has to be changed for calculating the implicit representation of the elliptical cylinder, i.e., for representing the points of the elliptical cylinder in the standard coordinate system. It will now be described how points in one coordinate system can be represented in a second coordinate system. In the following example, the point P is represented relative to the coordinate system 2, with its representation in the coordinate system 1 to be computed. The relationship to the corresponding systems will be indicated by superscripted numbers. For example, P⁽¹⁾ indicates the homogeneous representation of the point P in the coordinate system 1. The point P is represented in coordinate system 2 using homogeneous coordinate as follows

[0210] $P^{\left(2\right)}=\left(\begin{array}{c}1\\ 2\\ 1\end{array}\right).$

[0211] If the representation of the point in the coordinate system 1 is of interest, then one has to change into that coordinate system. This is done by multiplying the transformation matrix with the vector P⁽²⁾.

[0212] In general,
P⁽ⁱ⁾=M_i←j·P^(j),
wherein M_i←j represents the transformation which converts the representation of a point in the coordinate system j into its representation in the system i. In the given example, the above equation is then
P⁽¹⁾=M_1←2·P⁽²⁾,
i.e., the matrix M_1←2 has to be calculated.

[0213] In general,
F=M_i←j⁻¹=M_j←i,
whereby F represents the transformation that maps the base vectors of the coordinate system i onto those of the coordinate system j, as represented in the coordinates of system i.

[0214] Applied to the example, it follows

[0215] $F=\left(\begin{array}{ccc}1& 0& t_1\\ 0& 1& t_2\\ 0& 0& 1\end{array}\right)=\left(\begin{array}{ccc}1& 0& 3\\ 0& 1& 2\\ 0& 0& 1\end{array}\right)=M_{2\leftarrow 1}^{-1},$
since the transformation F only represents the translation (t₁, t₂). As a result

[0216] $P^{\left(1\right)}=M_{1\leftarrow 2}·P^{\left(2\right)}=M_{2\leftarrow 1}^{-1}·P^{\left(2\right)}=\left(\begin{array}{ccc}1& 0& 3\\ 0& 1& 2\\ 0& 0& 1\end{array}\right)·\left(\begin{array}{c}1\\ 2\\ 1\end{array}\right)=\left(\begin{array}{c}4\\ 4\\ 1\end{array}\right).$

[0217] From the relationship
X^(EZ)=M_EZ←SX^(S)
the following equation is obtained directly for the elliptical cylinder
(X^(EZ))^TQX^(EZ)=(M_EZ←SX^(S))^TQM_EZ←SX^(S)=(X^(S))^TM_EZ←S^TQM_EZ←SX^(S)=(X^(S))^T{tilde over (Q)}X^(S)=0
with the symmetric matrix {tilde over (Q)}=M_EZ←S^TQM_EZ←S, i.e., the desired representation of the elliptical cylinder in the standard coordinate system.

[0218] Finally, the matrix M_EZ←S has to be determined. Taking into account the formal relationship
F=M_EZ←S⁻¹F⁻¹=M_EZ←S,
wherein F is the transformation matrix which transforms the standard base vectors into the base vectors of the coordinate system EZ. This transformation is composed of a rotation R followed by a translation T about P1 (see FIG. 11). As a result,

[0219] $M_{\mathrm{EZ}\leftarrow S}=F^{-1}={\left(\mathrm{TR}\right)}^{-1}=R^{-1}{T}^{-1}=R^T{T}^{-1}=\left(\begin{array}{cccc}& & & 0\\ & {\Phi }^T& & 0\\ & & & 0\\ 0& 0& 0& 1\end{array}\right)·\left(\begin{array}{cccc}1& 0& 0& \\ 0& 1& 0& -p_1\\ 0& 0& 1& \\ 0& 0& 0& 1\end{array}\right)=\left(\begin{array}{cccc}& & & \\ & {\Phi }^T& & -{\Phi }^T{p}_1\\ & & & \\ 0& 0& 0& 1\end{array}\right)$
wherein the (3×3)-rotation matrix Φ after the construction has the following form:

[0220] $\Phi =\left(\begin{array}{ccc}\hat{b}& \hat{a}& \hat{n}\end{array}\right),$
resulting in the following representation of the quadrics in the standard coordinate system:

[0221] ${\left(X^{\left(S\right)}\right)}^T\tilde{Q}{X}^{\left(S\right)}={\left(X^{\left(S\right)}\right)}^T{\left(\begin{array}{cccc}& & & \\ & {\Phi }^T& & -{\Phi }^T{p}_1\\ & & & \\ 0& 0& 0& 1\end{array}\right)}^T·\left(\begin{array}{cccc}\frac{1}{{b}^2}& 0& 0& 0\\ 0& \frac{1}{{a}^2}& 0& 0\\ 0& 0& 0& 0\\ 0& 0& 0& -1\end{array}\right)·{\left(\begin{array}{cccc}& & & \\ & {\Phi }^T& & -{\Phi }^T{p}_1\\ & & & \\ 0& 0& 0& 1\end{array}\right)}^T{X}^{\left(S\right)}=0.$

[0222] It should finally be mentioned that an approach according to solution trial 2, i.e., computing the intersection in the coordinate system of the elliptical cylinder, is performed with the quadrics Q represented in standard form. Accordingly, no additional manipulation is required.

[0223] Computation of the z-intersection with a Smaller z-coordinate

[0224] The intersecting points can be computed in the corresponding coordinate systems, taking advantage of the aforedescribed results, by solving only a quadratic equation. According to the aforedescribed solution trials, a selection can be made between computing the intersecting point in the standard coordinate system or in the coordinate system of the elliptical cylinder.

[0225] Computing the Intersection in the Coordinate System 1 (Solution Trail 1)

[0226] The elliptical cylinder has the following form in the standard coordinate system

[0227] ${\left(X^{\left(S\right)}\right)}^T\tilde{Q}{X}^{\left(S\right)}={\left(X^{\left(S\right)}\right)}^T\left(\begin{array}{cccc}{\tilde{q}}_{11}& {\tilde{q}}_{12}& {\tilde{q}}_{13}& {\tilde{q}}_{14}\\ {\tilde{q}}_{12}& {\tilde{q}}_{22}& {\tilde{q}}_{23}& {\tilde{q}}_{24}\\ {\tilde{q}}_{13}& {\tilde{q}}_{23}& {\tilde{q}}_{33}& {\tilde{q}}_{34}\\ {\tilde{q}}_{14}& {\tilde{q}}_{24}& {\tilde{q}}_{34}& {\tilde{q}}_{44}\end{array}\right)X^{\left(S\right)}=0,$
while the straight line has the form

[0228] $g^{\left(S\right)}\text{:}{X}^{\left(S\right)}=\left(\begin{array}{c}{x}_g^{\left(S\right)}\\ y_g^{\left(S\right)}\\ 0\end{array}\right)+\lambda \left(\begin{array}{c}0\\ 0\\ 1\end{array}\right),\lambda \in 3.$

[0229] By inserting the equation for the straight line into the representation of the elliptical cylinder, one obtains with x=x_g^(S) and y=y_g^(S)
{tilde over (q)}₃₃λ²+2({tilde over (q)}₁₃x+{tilde over (q)}₂₃y+{tilde over (q)}₃₄)λ+({tilde over (q)}₁₁x²+2{tilde over (q)}₁₂xy+2{tilde over (q)}₁₄x+{tilde over (q)}₂₂y²+2{tilde over (q)}₂₄y+{tilde over (q)}₄₄)=0
and hence

[0230] ${\lambda }_{+,-}=\frac{\begin{array}{c}-\left({\tilde{q}}_{13}x+{\tilde{q}}_{23}y+{\tilde{q}}_{34}\right)\pm \\ \sqrt{\begin{array}{c}{\left({\tilde{q}}_{13}x+{\tilde{q}}_{23}y+{\tilde{q}}_{34}\right)}^2-\\ {\tilde{q}}_{33}\left({\tilde{q}}_{11}{x}^2+2{\tilde{q}}_{12}\mathrm{xy}+2{\tilde{q}}_{14}x+{\tilde{q}}_{22}{y}^2+2{\tilde{q}}_{24}y+{\tilde{q}}_{44}\right)\end{array}}\end{array}}{{\tilde{q}}_{33}}.$

[0231] The lower intersecting point, i.e., the intersecting point with a smaller z^(S)-coordinate, can then be represented as

[0232] $S=\left(\begin{array}{c}{x}_g^{\left(S\right)}\\ y_g^{\left(S\right)}\\ {\lambda }_{-}\end{array}\right).$

[0233] Computing the Intersection in the Coordinate System EZ (Solution Trial 2)

[0234] For computing the intersecting points in the coordinate system of the elliptical cylinder, the straight line has to be first represented in its coordinate system. The transformed end point has then the form

[0235] $\left(\begin{array}{c}{x}_g^{\left(\mathrm{EZ}\right)}\\ y_g^{\left(\mathrm{EZ}\right)}\\ z_g^{\left(\mathrm{EZ}\right)}\\ 1\end{array}\right)=F^{-1}\left(\begin{array}{c}{x}_g^{\left(S\right)}\\ y_g^{\left(S\right)}\\ 0\\ 1\end{array}\right),$
yielding for the transformed the direction vector of the straight line

[0236] $r^{\left(\mathrm{EZ}\right)}=\left(\begin{array}{c}{r}_x^{\left(\mathrm{EZ}\right)}\\ r_y^{\left(\mathrm{EZ}\right)}\\ r_z^{\left(\mathrm{EZ}\right)}\end{array}\right)=\Phi ·\left(\begin{array}{c}0\\ 0\\ 1\end{array}\right),$
so that the straight line as a whole has the form

[0237] $g^{\left(\mathrm{EZ}\right)}\text{:}{X}^{\left(\mathrm{EZ}\right)}=\left(\begin{array}{c}{x}_g^{\left(\mathrm{EZ}\right)}\\ y_g^{\left(\mathrm{EZ}\right)}\\ z_g^{\left(\mathrm{EZ}\right)}\end{array}\right)+\mu r^{\left(\mathrm{EZ}\right)},\mu \in 3.$

[0238] It should be noted that only the end point of the straight line and not its direction vector is translated.

[0239] After the transformed equation for the straight line is inserted into the expression (X^(EZ))^TQX^(EZ)=0, one obtains

[0240] ${\left(\frac{1}{{b}^2}{\left(r_x^{\left(\mathrm{EZ}\right)}\right)}^2+\frac{1}{{a}^2}{\left(r_y^{\left(\mathrm{EZ}\right)}\right)}^2\right)}^2{\mu }^2+\left(\frac{2x_g^{\left(\mathrm{EZ}\right)}{r}_x^{\left(\mathrm{EZ}\right)}}{{b}^2}+\frac{2y_g^{\left(\mathrm{EZ}\right)}{r}_y^{\left(\mathrm{EZ}\right)}}{{a}^2}\right)\mu +\frac{{\left(x_g^{\left(\mathrm{EZ}\right)}\right)}^2}{{b}^2}+\frac{{\left(y_g^{\left(\mathrm{EZ}\right)}\right)}^2}{{a}^2}-1=0.$

[0241] The two intersecting points

[0242] $S_{-}^{\left(\mathrm{EZ}\right)}=\left(\begin{array}{c}{x}_g^{\left(\mathrm{EZ}\right)}\\ y_g^{\left(\mathrm{EZ}\right)}\\ y_g^{\left(\mathrm{EZ}\right)}\end{array}\right)+{\mu }_{-}\left(\begin{array}{c}{r}_x^{\left(\mathrm{EZ}\right)}\\ r_y^{\left(\mathrm{EZ}\right)}\\ r_z^{\left(\mathrm{EZ}\right)}\end{array}\right),S_{+}^{\left(\mathrm{EZ}\right)}=\left(\begin{array}{c}{x}_g^{\left(\mathrm{EZ}\right)}\\ y_g^{\left(\mathrm{EZ}\right)}\\ y_g^{\left(\mathrm{EZ}\right)}\end{array}\right)+{\mu }_{+}\left(\begin{array}{c}{r}_x^{\left(\mathrm{EZ}\right)}\\ r_y^{\left(\mathrm{EZ}\right)}\\ r_z^{\left(\mathrm{EZ}\right)}\end{array}\right)$
can be computed with the help of the parameters μ_+,−, which are obtained by solving the quadratic equation (if a solution exists). The intersecting points are available at this time only in the coordinate system of the elliptical cylinder. Accordingly, it cannot be decided yet which of the intersecting points can be considered to represent the lower intersecting point in the standard coordinate system based on its smaller z^(S)-coordinate. The inverse transformation of S₋^(EZ) and S₊^(EZ) is again performed using the mapping matrix F

[0243] $S_{-}^{\left(S\right)}=F·\left(\begin{array}{c}{S}_{-}^{\left(\mathrm{EZ}\right)}\\ 1\end{array}\right),S_{+}^{\left(S\right)}=F·\left(\begin{array}{c}{S}_{+}^{\left(\mathrm{EZ}\right)}\\ 1\end{array}\right).$

[0244] One of the two points can now be disregarded based on its larger z^(S)-coordinate.

[0245] Testing the Location of the Intersecting Point within the Finite Elliptical Cylinder

[0246] According to solution trial 1 or solution trial 2, the last calculated lower intersecting point is located on the infinitely extending elliptical cylinders, based on its implicit representation. When the cylindrical cutter moves from point P1 in the direction P1P2 to point P2, only those intersecting points are relevant which are located between the planes E1 and E2 (see FIG. 13).

[0247] This can be typically checked afterwards with a normal test. It is initially checked if

[0248] $\cos \alpha =\frac{\mathrm{P1S}·\mathrm{P1P2}}{\mathrm{P1S}·\mathrm{P1P2}}=\frac{p}{\mathrm{P1S}}\iff p=\frac{\mathrm{P1S}·\mathrm{P1P2}}{\mathrm{P1P2}}\iff \mathrm{P1S}·\mathrm{P1P2}=p·\mathrm{P1P2}<{\mathrm{P1P2}}^2,$
wherein α refers to the angle enclosed between P1S and the projection p. The following relevance condition can then be deduced
0<P1S·P1P2<|P1P2|².

[0249] The inequality on the left states that the intersecting point is located above the lower plane E₁, whereas the inequality on the right states that the projection of the intersecting point S onto the line P1P2 is shorter than the distance P1P2. This ensures that the intersecting point is located within the planes E₁ and E₂. However, in the given example, the planes E₁ and E₂ are parallel to the xy-plane, so that the aforedescribed general approach for the normal test can be simplified considerably.

[0250] It only has to be tested if the z^(S)-coordinate of the intersecting point is located between the z^(S)-coordinate of the points P1 and P2. Accordingly, if formulated for an approach according to solution trial 1, i.e., the computation of the intersecting points in the standard coordinate system, one obtains
min(P1_z^(S), P2_z^(S))≦λ_≦max(P1_z^(S), P2_z^(S)).

[0251] Approach with a Horizontal Milling Edge

[0252] If the milling edge extends the horizontal direction, then the boundary surface of the elliptical cylinder cannot be represented as a quadrics, because in this case its shorter major axes b becomes the null vector so that the elliptical cylinder degenerates to a rectangle. However, since in the representation of the quadrics a division by the length b is required, this approach does not lead to the desired result. However, the z^(S)-coordinate of an existing intersecting point is identical to the z^(S)-coordinate of the point P1. Accordingly, it has to be tested only if the corresponding straight line intersects with the sweep triangle (see FIG. 14).

[0253] For this purpose, it is initially tested, for example, if the end point of the straight line is located between the bounding straight lines g₁ and the g₂, and h₁ and h₂, respectively. Because the latter are typically not parallel to the x^(S)- and/or y^(S)-axis, the aforedescribed general normal test has to be applied. The z^(S)-coordinate can hereby be disregarded, since it is a planar problem. The following relationships has to be tested:
0≦P1S·P1P2≦|P1P2|²
for the pair of straight lines g₁ and g₂, and
0≦Z1S·Z1Z2≦|Z1Z2|²=4r²
for the pair of straight lines h₁ and h₂.

[0254] In summary, the straight line then intersects the “elliptical cylinder” and the point is relevant for further processing.

[0255] Quick Elimination Test

[0256] Using the same quick elimination test already described above with reference to the embodiment of a spherical cutter, the method can likewise be further optimized for a cylindrical cutter. In this case, all vertical straight lines that never intersect with the elliptical cylinder are disregarded. The elliptical cylinder is embedded in a circular cylinder with the same major axes P1P2 and the radius r. If the straight line misses the circular cylinder, then also no intersecting point with the elliptical cylinder exists, making further computational steps obsolete.

[0257] For this purpose, according to

[0258] $d=\frac{r\times \mathrm{P1P2}}{r\times \mathrm{P1P2}}$
a normalized vector d which is orthogonal to r (direction vector of the straight line, see FIG. 11) and P1P2 is computed. If the distance between the two (typically) inclined straight lines defined by r and P1P2,
d=|((x_gy_g0)^T−P1)·d|
is greater than the radius r of the circular cylinder, then the straight line and the cylinder do not intersect.

[0259] For adapting the exemplary method as described above for a spherical cutter and a cylindrical cutter to other types of cutters, the milling volume swept for each milling edge can be composed of suitable partial volumes in analogy to the two aforedescribed embodiments and expressed implicitly with analytical formulas. For forming the intersection with the vertical straight lines, typically a method for determining the zero points is applied. Accordingly, a method according to the invention adapted for more general cutter types can typically not compete with the variant spherical cutter with respect to the delay time characteristic.

[0260] A more generalized milling process with five axes instead of three axes has to take into account also the spatial location of the milling tool and its orientation. An approach which only stores the intersecting point with the smallest z-coordinate, can then no longer be used. For example, with a cutter oriented in the positive z-direction, the upper intersecting point, i.e. the intersecting point with the maximum z-coordinate, has to be recorded. In particular when milling undercuts or lateral pockets, several intersecting points for each vertical straight line are typical, so that a more general process must be contemplated.

[0261] While the invention has been illustrated and described in connection with currently preferred embodiments shown and described in detail, it is not intended to be limited to the details shown since various modifications and structural changes may be made without departing in any way from the spirit of the present invention. The embodiments were chosen and described in order to best explain the principles of the invention and practical application to thereby enable a person skilled in the art to best utilize the invention and various embodiments with various modifications as are suited to the particular use contemplated.

[0262] What is claimed as new and desired to be protected by Letters Patent is set forth in the appended claims and their equivalents: