%!PS % An improved Bezier Cubic Spline length finder and LINEAR(!) subdivider % ======================================================================= % by Don Lancaster % Copyright c 2006 by Don Lancaster & Synergetics, Box 809, Thatcher, AZ, 85552 % (928) 428-4073 Email: don@tinaja.com Website: http://www.tinaja.com % Consulting services available http://www.tinaja.com/info01.html % All commercial rights and all electronic media rights ~fully~ reserved. % Linking usually welcome. Reposting expressly forbidden. Version 1.1 %%%%%%%%%%% % internal gonzo utilities provided. %%%%%%%%%%% % Plus the next line if Gonzo is to be activated... % gonzo begin ps.util.1 begin nuisance begin % NOTE THAT ALL PS FILENAME STRINGS !!!DEMAND!!! DOUBLE REVERSE SLASHES. % A collection of utilities that let you find the length of a full or partial cubic % spline using actual along-the-curve measurement for best possible accuracy. % A "S-list" method provides accurate and constant spacing. Spacing is uniform, % compared to "constant-T" methods which are closer together along the "more bent" % portions of the cuve. % Methods presented are accurate and conceptually simple but computationally intensive. % Run modified text as PostScript ASCII textfile sent to Acrobat Distiller. %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % % % TEMPORARY COPIES OF GONZO ROUTINES FOR DEMOS BELOW % % % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Use of "real" gonzo is strongly encouraged. % Free download from http://www.tinaja.com/post01.asp#gonzo /xrpt{gsave aload pop /trips exch def /dist exch def /rproc exch def trips { gsave rproc grestore dist 0 translate } repeat grestore} def /yrpt{gsave aload pop /trips exch def /dist exch def /rproc exch def trips { gsave rproc grestore 0 dist translate } repeat grestore} def /fat5 true def % 3 pixels wide every fifth line? /fatter10 true def % 5 pixels wide every tenth line? /fatterborder true def % 5 pixel gray outline border? /setgrid { gsave /blocksize exch def translate blocksize dup scale} def /thingridlines 0 store /ggray {60 0 {sub abs 0.21 lt {1}{0} ifelse} setscreen gsave 0 60 0 {pop pop 1 add 1} setscreen grestore dup 25 eq {pop 0.8}{100 eq {0.72}{0.765} ifelse} ifelse setgray 0.3 1 0.5 setrgbcolor % modded - not sure why } def /showgrid {gsave ggray /vblocks exch def /hblocks exch def thingridlines setlinewidth [{0 0 moveto 0 vblocks rlineto stroke} 1 hblocks 1 add] xrpt [{0 0 moveto hblocks 0 rlineto stroke} 1 vblocks 1 add] yrpt fatterborder { gsave newpath 0 0.96 blocksize div dtransform round idtransform setlinewidth pop 2 setlinecap 0 0 moveto hblocks 0 rlineto 0 vblocks rlineto hblocks neg 0 rlineto closepath stroke grestore} if fat5 { gsave newpath 0 0.48 blocksize div dtransform round idtransform setlinewidth pop mark {5 0 moveto 0 vblocks rlineto stroke} 5 hblocks 5 div cvi] xrpt mark {0 5 moveto hblocks 0 rlineto stroke} 5 vblocks 5 div cvi] yrpt grestore} if fatter10 { gsave newpath 0 0.96 blocksize div dtransform round idtransform setlinewidth pop mark {10 0 moveto 0 vblocks rlineto stroke} 10 hblocks 10 div cvi] xrpt mark {0 10 moveto hblocks 0 rlineto stroke} 10 vblocks 10 div cvi] yrpt grestore} if grestore} def /line1 {.06 dup setlinewidth 5 mul /erase exch def} def /line2 {.12 dup setlinewidth 5 mul /erase exch def} def /line3 {.18 dup setlinewidth 5 mul /erase exch def} def /m {moveto} store /dot { currentpoint newpath 0.150 0 360 arc fill } def /mdot { m dot} def /mergestr {2 copy length exch length add string dup dup 4 3 roll 4 index length exch putinterval 3 1 roll exch 0 exch putinterval} def %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % % % CUBIC SPLINE SERVICE UTILITIES % % % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %% code to convert x0y0...x3y3 influence points to A-H equation constants. % See http://www.tinaja.com/glib/bezmath.pdf for a derivation... /findAH { /A x3 x2 3 mul sub x1 3 mul add x0 sub store /E y3 y2 3 mul sub y1 3 mul add y0 sub store /B x2 3 mul x1 6 mul sub x0 3 mul add store /F y2 3 mul y1 6 mul sub y0 3 mul add store /C x1 3 mul x0 3 mul sub store /G y1 3 mul y0 3 mul sub store /D x0 store /H y0 store } store % code to find current x and y given A-H and t and using "cubless" method % input tt, exit with curx or cury defined but not on stack /findx {/ttt exch store A ttt mul B add ttt mul C add ttt mul D add /curx exch store} store /findy {/ttt exch store E ttt mul F add ttt mul G add ttt mul H add /cury exch store} store % code to find curent s incremental chord length given curx cury prevx prevy... % the usual sum-of squares of the x and y increment. Returns to curs but not on stack. /finds {curx prevx sub dup mul cury prevy sub dup mul add sqrt /curs exch store} store % key to the routines are three arrays of identical length... % The tlist is an array of Bezier t values, initially constant spaced. % The slist is the RUNNING TOTAL of the chord lengths, ending with length. % The nslist is a NORMALIZED 0-1 version of the slist. % code to make a t list of numchords+1 samples. In the case of 100 samples, % would be [0 0.01 0.02 .... 1.00] /numchords 100 store % recommended default /makecurtlist {/curtlist mark 0 numchords {dup 1 numchords div add} repeat ] store } store % code to make a slist from a tlist. Each entry is SUMMED distance to the tlist value; % final value is Bezier length. /makecurslist { /prevx 0 store % initialize previous values /prevy 0 store /curslist mark 0 % start slist array 0 1 numchords {curtlist exch get /tt exch store % grab current tt findx tt findy % find curx and cury finds % current chord curs add dup % integrate running result and dup for next /prevx curx store % save new chord end /prevy cury store } for pop % undo final integration! ] store % save slist } store % /findsplinelen reads the last entry in the slist, equal to the total spline length. /findsplinelen {curslist dup length 1 sub get /splinelen exch store} store % /findtofs finds an interpolated t value for a given s distance. The fractional % distance through a s array entry is calculated and converted to dt/ds. Once t % is known, the corresponding x and y positions are easily found. % findbezlen must be previously run. /findtofs { /ss exch store % must be in absolute range of 0 to length. 1 1 curslist length 1 sub % scan each s value { /ii exch store ss curslist ii get % view each accumulated s in turn till too big le {exit} if } for % then exit with high s location /ii ii 1 sub store % back up one so s is in box curslist ii 1 add get curslist ii get sub % current s delta ss curslist ii get sub % excess into interval exch div % current interpolated s fraction curtlist ii 1 add get % t delta curtlist ii get sub mul % fractional t interpolation distance curtlist ii get add % final t for selected s } store % Yeah, the above "count till overflow" is somewhat crude, but it seems more % than fast enough. You can speed it up by a binary pretree (0, 0.25, 0.5 0.75 etc) % Or else by "inverting" fot a linear curslist and nonlinear curtlist. %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % % % HIGH LEVEL BEZIER LENGTH & SUBDIVISION PROCS % % % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % /findbezlen calculates the total Bezier spline length % and returns it as /splinelen. Requires predefined x0y0...x3y3 % and predefined numchords. (100 recommended). /findbezlen {findAH % calculate cubic coefficients makecurtlist % create linear t array makecurslist % create nonlinear s array findsplinelen % report length as splinelen } store % /findxyofs finds the x and y position for a given sublength position. % The sublength position is presently entered NORMALIZED in the 0 to 1.0 % range. findbezlen must be run previously. /findxyofs {splinelen mul % denormalize findtofs % find the interpolated t for the LINEAR s dup findx % get x(t) findy % get y(t) curx cury % and return to stack } store %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % % % DEMOS - REMOVE OR ALTER BEFORE REUSE % % % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%% %%%% (A) Find the length of a Bezier cubic spline %%%%%% %%%%%%%%%%%%% % Spline coefficients or defaults need entered here... /testspline { /x0 0 store /y0 0 store /x1 3 store /y1 20 store /x2 16 store /y2 4 store /x3 20 store /y3 20 store } store testspline % enter values /numchords 100 store % new entry here if not default=100 findbezlen % find length of spline % and report to log file... (\n\nThe length of the test spline is ) splinelen 20 string cvs mergestr (.\n) mergestr print flush %%%%%%%%%%%%%% % Accuracy check of a spline of known length of 28.28427... %%%%%%%%%%%%%% /x0 0 store /y0 0 store /x1 0.001 store /y1 0.001 store /x2 19.999 store /y2 19.999 store /x3 20 store /y3 20 store findbezlen % find length of spline (\nThe length of the 28.28427 known length spline is ) splinelen 20 string cvs mergestr (.\n) mergestr print flush %%%%%%%%%%%%%%% % Accuracy check for increasing numchords. Values >> 10,000 may trip PS math limits. %%%%%%%%%%%%%%% testspline % enter values /numchords 10 store % new entry here if not default=100 findbezlen % find length of spline % and report to log file... (\nThe length of the 10 sample test spline is ) splinelen 20 string cvs mergestr (.\n) mergestr print flush testspline % enter values /numchords 100 store % new entry here if not default=100 findbezlen % find length of spline % and report to log file... (The length of the 100 sample test spline is ) splinelen 20 string cvs mergestr (.\n) mergestr print flush testspline % enter values /numchords 1000 store % new entry here if not default=100 findbezlen % find length of spline % and report to log file... (The length of the 1000 sample test spline is ) splinelen 20 string cvs mergestr (.\n) mergestr print flush testspline % enter values /numchords 10000 store % new entry here if not default=100 findbezlen % find length of spline % and report to log file... (The length of the 10000 sample test spline is ) splinelen 20 string cvs mergestr (.\n\n) mergestr print flush %%%%%%%%%%%%%%% % Find x and y values a given distance into a spline %%%%%%%%%%%%%%% testspline % enter values /numchords 100 store % new entry here if not default=100 findbezlen % find length of spline 0.271 findxyofs (\n\nGoing 0.271 of the distance into testspline results in an x value of ) curx 20 string cvs mergestr ( and a y value of ) mergestr cury 20 string cvs mergestr (.\n) mergestr print % Plotted in red in next example below. %%%%%%%%%%%%%%% % Subdivide a cubic spline into 20 EQUAL chord length parts %%%%%%%%%%%%%%% % requires full or partial gonzo utilities 100 100 10 setgrid % draw a pretty graph 20 20 showgrid testspline % enter values /numchords 100 store % new entry here if not default=100 findbezlen % find length of spline x0 y0 moveto x1 y1 x2 y2 x3 y3 curveto line1 stroke 0 0.05 1.001 { findxyofs m dot} for 1 0 0 setrgbcolor 2.67118 7.9798 m dot % plot previous example result showpage % eof