Essays/Hidden Line Removal
Hidden Line Removal
Hidden Line Removal is used in wireframe graphics, to delete lines that are hidden by nearer polygons.
Here is a J verb to perform hidden line removal. The x operand is a table of triangles in 3D (shape nx3x3), and the y operand is an edge consisting of 2 points in 3D (shape 2x3). The result is a table containing the portions of the edge y that remain visible. There may be 0, 1, or more such portions.
This code assumes that the (x,y) values are the screen positions of the lines and triangles. The z values are for occlusion calculation only and will be either the orthographic z-value or the perspective pseudo-z; in either case, the triangle is planar in xyz.
Frontfacing triangles should have CCW winding in (x,y). Higher z values occlude lower z values.
usedtocull =: 1 : 'u#]' ifany =: ^:(*@#@]) NB. y is a line, start and end points 2x3 NB. x is a nx3x3 array of triangles, frontfaces in CCW order NB. Result is a table of visible 2-D segments for the line, with Z component removed NB. (2 3 3 $0 0 1 1 0 1 0 1 1 1 0 1 1 1 1 0 1 1) visiblesegments 2 3 $ 0 0 0 2 2 3 visiblesegments =: 4 : 0"3 2 line2d =. 2 {."1 line =. y NB. Discard triangles that do not intersect the bounding square for the line NB. Calculate bounding square x,y, min,:max NB. Exclusion if min>max or (-max)>(-min) NB. If any exclusion, discard the tri qtri =. line ([: -. ((+./@,"2)@:(>"2&:(1 _1&*"2)) |.)~&:((<./ ,: >./)"2@:(2&{."1))) usedtocull x NB. Calculate volume to each endpoint (+ volume means the face is in front NB. of the point, having higher z, i. e. line is occludable) nx2 epvol =. qtri -/ . *@:(-"1)"2 1/ line NB. Keep faces that are in front of at least one endpoint (or have 0 volume, which NB. is most likely the face that contributed the line) ftri =. (+./"1 epvol > 0) # qtri epvol =. (+./"1 epvol > 0) # epvol NB. Create an interval of applicability for each face. This is [0,crosspoint) NB. or (crosspoint,0] for interpenetrations, or [0,1] for faces in front of NB. the edge NB. 0,1 if both + NB. 0,cross if +,- NB. cross,1 if -,+ applintvl =. 0 1 +"1 (0&<. % -/"1) epvol NB. For each triangle, calculate 2-D area from each point to each edge. nx3x2 NB. The sign will be + if the point is on the inside of the edge areas =. 0 2 1 |: ftri (-/ . *)@:(0 1 2&(}.@|."0 _)"2)@:(-"1"2 1"2)&:(2&{."1) line NB. Convert the areas to intervals in (0,1) giving the part of the line that NB. is inside the edge of the triangle nx3x2 edgeintvls =. 0 1 +"1 (0&<. % -/"1) areas NB. Take the intersection of the intervals for each triangle, to give the NB. interval of occlusion for the triangle NB. Intersect the occlusion interval with the applicability interval to give NB. the total occluded interval for each triangle nx2 occludeintvl =. (>./&.:(1 _1&*"1))"2 edgeintvls ,"2 1 applintvl NB. Union the occluded intervals for all triangles to get the total occluded area. NB. We go through pairs in ascending order of start point. Result at each point NB. is, if overlap, lower start,larger end; if no overlap, keep the second intvl NB. Then discard intervals that share a start point occludeintvl =. \:~ </"1 usedtocull occludeintvl NB. Discard empty or invalid intervals occludeintvl =. ({./.~ {."1) ({.@[ , >./&{:)^:({:@[ >: {.@])~/\. ifany occludeintvl NB. Complement the occluded areas to produce final result: list of visible segments, in 2D ({. line) +"1 (-~/line) */~ (0.001 < -~/"1) usedtocull _2 ]\ 0 , (,occludeintvl) , 1 )
This verb can be used to create a hidden-surface-removed plot of a parametric surface, as in the following example.
load 'plot' endtoend =: 1 : ';@:(<@u)' createlines =: 3 : 0 NB. Create the polygons for the object, in uv space uvgrid =. (,. {."_1) (, {.) (2p1 * (%~ i.) 60) ,"0/ (2p1 * (%~ i.) 20) NB. Apply the function to get mesh in world space worldgrid =. ((cos@[ * 2 + cos@]) , (sin@[ * 2 + sin@]) , (3 * sin@[ * cos@]))/"1 uvgrid NB. Calculate viewing matrix (world-to-view) NB. Apply azimuth first, then elevation. nxmx3 azmtx =. ((cos,sin,0:) , (-@sin,cos,0:) ,: (0,0,1:)) 2p1*214%360 elmtx =. ((cos,0,-@sin) , (0,1,0:) ,: (sin,0,cos)) 2p1*220%360 viewgrid =. (elmtx mp azmtx)&mp&.|:"2 worldgrid NB. Create triangles (CCW winding). nxmx2x3x3 trimesh =. 2 2 3 ((0 1 3,:3 2 0) { ,/);._3 viewgrid NB. Cull backfaces. nx3x3 fronttris =. (#~ 0 <: -/ . *@:(}. -"1 {.)"2@:(0 1&{"1)) ,/^:3 trimesh NB. Create table of lines. nx2x3 trilines =. ~. /:~"2 ,/ (2 ]\ (,{.))"2 fronttris NB. Produce the visible portion of each line vislines =. fronttris visiblesegments endtoend trilines NB. Truncate to 2D 2 {."1 vislines ) 'type line;color 0 0 0' plot ~. /:~"1 ,/^:(0 >. 2 -~ #@$) j./"1 createlines''
About 80% of the time in this example is taken calculating qtri in visiblesegments, which suggests that a spatial index on the triangles would provide a performance improvement.
Contributed by Henry Rich