A.M. Hoornweg 126 Posted May 20 Hello all, does anyone know a nice algorithm that will connect the points in a polygon in a rounded fashion? Its intended use is a contour map. The polygon must go through all vertices but I'd like to avoid sharp corners for a more natural look; it's OK if the lines between the vertices are slightly curved. The few solutions I found on the internet tend to bend without touching the vertices. Share this post Link to post
ULIK 6 Posted May 20 Have you already tried it with a combination of splines? They would touch the vertices. Share this post Link to post
Sherlock 472 Posted May 20 From the Bézier curve you'll eventually have to move on to (Bézier) splines to get from a polygon to a bezigon. And not sure about VCL, but FMX has some Bézier procedures. Edith says: Naturally VCL has it, because MS has it: https://docwiki.embarcadero.com/Libraries/Sydney/en/Vcl.Graphics.TCanvas.PolyBezier Share this post Link to post
vfbb 216 Posted May 20 If you don't mind third-party solutions, Skia has an effect that applies rounded edges automatically and works for Console, Vcl and FMX: https://github.com/skia4delphi/skia4delphi/blob/eed4afbf8a34137a9bfa308bcb5ef87cee84abcb/Samples/Demo/FMX/Source/Sample.Form.PathsAndEffects.pas#L193 Result: Share this post Link to post
A.M. Hoornweg 126 Posted May 20 (edited) It's not trivial at all. I've downloaded a random contour map from the internet and attach it here. This is a map of a terrain, consisting of polygons that describe points of equal altitude above sea level. So there's a polygon describing an altitude of 350 metres, another one that corresponds to 360 metres, another one for 370 metres etcetera. The problem: The polygons may never intersect (any point of the topography has only one altitude). If an automated routine is to be used to create smooth connections between points of a polygon, it is necessary to avoid "wild deflections" that would make these lines cross adjacent polygons. A human being would have no problem drawing such connections by hand, intuitively. But an algorithm ? Edited May 20 by A.M. Hoornweg Share this post Link to post
A.M. Hoornweg 126 Posted May 20 3 minutes ago, vfbb said: If you don't mind third-party solutions, Skia has an effect that applies rounded edges automatically and works for Console, Vcl and FMX: https://github.com/skia4delphi/skia4delphi/blob/eed4afbf8a34137a9bfa308bcb5ef87cee84abcb/Samples/Demo/FMX/Source/Sample.Form.PathsAndEffects.pas#L193 Result: Such rounded edges "miss" the points of the polygon, they deviate before hitting the point. Share this post Link to post
Sherlock 472 Posted May 20 Ok, so the points need to be hit, but the lines in between may be rounded. That might not be covered by splines to well. 😞 Well strike that, of course the Points will be hit. It is "just" a matter of choosing the controlpoints well to get a decent shape. But therein lies the rub I guess. Share this post Link to post
vfbb 216 Posted May 20 2 hours ago, A.M. Hoornweg said: Such rounded edges "miss" the points of the polygon, they deviate before hitting the point. Here's an example: uses Skia; function MakeCubicSplineInterpolation(const APoints: TArray<TPointF>): ISkPath; var LPathBuilder: ISkPathBuilder; LSegments: Integer; I: Integer; mx: Single; my: Single; LScratches: array of record a, b, c, r, p: TPointF; end; begin LPathBuilder := TSkPathBuilder.Create; if Length(APoints) < 2 then Exit(LPathBuilder.Detach); if Length(APoints) = 2 then begin LPathBuilder.MoveTo(APoints[0]); LPathBuilder.LineTo(APoints[1]); Exit(LPathBuilder.Detach); end; LSegments := Length(APoints) - 1; SetLength(LScratches, LSegments); LScratches[0].a := PointF(0, 0); LScratches[0].b := PointF(2, 2); LScratches[0].c := PointF(1, 1); LScratches[0].r := PointF(APoints[0].X + 2 * APoints[1].X, APoints[0].Y + 2 * APoints[1].Y); for I := 1 to LSegments - 2 do begin LScratches[I].a := PointF(1, 1); LScratches[I].b := PointF(4, 4); LScratches[I].c := PointF(1, 1); LScratches[I].r := PointF(4 * APoints[i].X + 2 * APoints[I + 1].X, 4 * APoints[I].Y + 2 * APoints[I + 1].Y); end; LScratches[LSegments - 1].a := PointF(2, 2); LScratches[LSegments - 1].b := PointF(7, 7); LScratches[LSegments - 1].c := PointF(0, 0); LScratches[LSegments - 1].r := PointF(8 * APoints[LSegments - 1].X + APoints[LSegments].X, 8 * APoints[LSegments - 1].Y + APoints[LSegments].Y); for I := 1 to LSegments - 1 do begin mx := LScratches[I].a.X / LScratches[I - 1].b.X; my := LScratches[I].a.Y / LScratches[I - 1].b.Y; LScratches[I].b := LScratches[I].b - PointF(mx * LScratches[I - 1].c.X, my * LScratches[I - 1].c.Y); LScratches[I].r := LScratches[I].r - PointF(mx * LScratches[I - 1].r.X, my * LScratches[I - 1].r.Y); end; LScratches[LSegments - 1].p := PointF(LScratches[LSegments - 1].r.X / LScratches[LSegments - 1].b.X, LScratches[LSegments - 1].r.Y / LScratches[LSegments - 1].b.Y); for I := Length(APoints) - 3 downto 0 do begin LScratches[I].p := PointF((LScratches[I].r.X - LScratches[I].c.X * LScratches[I + 1].p.X) / LScratches[I].b.X, (LScratches[I].r.Y - LScratches[I].c.Y * LScratches[I + 1].p.Y) / LScratches[I].b.Y); end; LPathBuilder.MoveTo(APoints[0]); for I := 0 to LSegments - 2 do begin LPathBuilder.CubicTo(LScratches[I].p, PointF(2 * APoints[I + 1].X - LScratches[I + 1].p.X, 2 * APoints[I + 1].Y - LScratches[I + 1].p.Y), APoints[I + 1]); end; LPathBuilder.CubicTo(LScratches[LSegments - 1].p, PointF(0.5 * (APoints[LSegments].X + LScratches[LSegments - 1].p.X), 0.5 * (APoints[LSegments].Y + LScratches[LSegments - 1].p.Y)), APoints[LSegments]); Result := LPathBuilder.Detach; end; procedure TForm1.SkPaintBox1Draw(ASender: TObject; const ACanvas: ISkCanvas; const ADest: TRectF; const AOpacity: Single); var LPaint: ISkPaint; LMyPoints: TArray<TPointF>; begin LMyPoints := [PointF(62, 511), PointF(162, 605), PointF(262, 610), PointF(362, 402), PointF(462, 959), PointF(562, 58), PointF(662, 272), PointF(762, 99), PointF(862, 759), PointF(962, 945)]; LPaint := TSkPaint.Create(TSkPaintStyle.Stroke); LPaint.Color := TAlphaColors.Red; LPaint.AntiAlias := True; LPaint.StrokeWidth := 3; LPaint.StrokeCap := TSkStrokeCap.Round; ACanvas.DrawPath(MakeCubicSplineInterpolation(LMyPoints), LPaint); LPaint.StrokeWidth := 10; LPaint.Color := TAlphaColors.Black; ACanvas.DrawPoints(TSkDrawPointsMode.Points, LMyPoints, LPaint); end; Result: Note: You don't need to use Skia, it was just a facilitator for the example. 6 Share this post Link to post
A.M. Hoornweg 126 Posted May 21 @vfbb: thanks for the example, I'll look into it! Share this post Link to post
Rollo62 363 Posted May 22 So you need to fit all datapoints, but want to avoid large overshootings, then this paper is maybe also interesting for you. https://towardsdatascience.com/numerical-interpolation-natural-cubic-spline-52c1157b98ac Share this post Link to post
Anders Melander 1129 Posted May 22 11 hours ago, Rollo62 said: So you need to fit all datapoints, but want to avoid large overshootings, then this paper is maybe also interesting for you. https://towardsdatascience.com/numerical-interpolation-natural-cubic-spline-52c1157b98ac Quote We will use the top-down approach and make sure you visualize while you’re reading to understand it better. I counted one photo of a Chinese paper lamp, two irrelevant meme pics, two general conceptual illustrations and 10 pages of math with no illustrations 😕 Share this post Link to post
Anders Melander 1129 Posted May 22 I just remembered that Graphics32 has two examples which demonstrates interpolation: https://github.com/graphics32/graphics32/tree/master/Examples/Drawing/CubicSpline https://github.com/graphics32/graphics32/tree/master/Examples/Drawing/Curves Both of these just uses Graphics32 for output. The curve generation is independent. Built into Graphics32 there's also the TCanvas32.CurveTo method which does cubic Bézier interpolation (4 control points) and the TCanvas32.ConicTo method which does quadratic Bézier interpolation (3 control points). 2 Share this post Link to post
Rollo62 363 Posted May 23 Ok, here is something to play around with natural spline interpolation vs. other kinds. https://tools.timodenk.com/cubic-spline-interpolation Maybe that has enough gamification factor to fulfil the high expectations here. Share this post Link to post
angusj 90 Posted May 23 (edited) Here's my GetSmoothPath() routine. It requires no specific graphics library to use, just a few extra functions (also included below). This function generates an array of control points that's very easily converted into a flattened cubic bezier path using just about any 2D graphics library. (nb: The code below has been written with simplicity as the focus rather than performance.) uses SysUtils, Math; type TPointD = record X, Y: double; end; TPathD = array of TPointD; TArrayOfDouble = array of double; function DistanceSqrd(const pt1, pt2: TPointD): double; begin result := Sqr(pt1.X - pt2.X) + Sqr(pt1.Y - pt2.Y); end; function Distance(const pt1, pt2: TPointD): double; begin Result := Sqrt(DistanceSqrd(pt1, pt2)); end; function OffsetPoint(const pt: TPointD; dx, dy: double): TPointD; begin result.x := pt.x + dx; result.y := pt.y + dy; end; function GetAvgUnitVector(const vec1, vec2: TPointD): TPointD; var inverseHypot: Double; begin Result.X := (vec1.X + vec2.X) * 0.5; Result.y := (vec1.Y + vec2.Y) * 0.5; inverseHypot := 1 / Hypot(Result.X, Result.Y); Result.X := Result.X * inverseHypot; Result.Y := Result.Y * inverseHypot; end; procedure MakeSymmetric(var val1, val2: double); begin val1 := (val1 + val2) * 0.5; val2 := val1; end; function GetUnitVector(const pt1, pt2: TPointD): TPointD; var dx, dy, inverseHypot: Double; begin if (pt1.x = pt2.x) and (pt1.y = pt2.y) then begin Result.X := 0; Result.Y := 0; Exit; end; dx := (pt2.X - pt1.X); dy := (pt2.Y - pt1.Y); inverseHypot := 1 / Hypot(dx, dy); dx := dx * inverseHypot; dy := dy * inverseHypot; Result.X := dx; Result.Y := dy; end; // GetSmoothPath - returns cubic bezier control points // parameters: 1. path for smoothing // 2. whether or not the smoothed path will closed // 3. percent smoothness (0..100) // 4. maximum dist control pts from path pts (0 = no limit) // 5. symmetric vs asymmmetric control pts function GetSmoothPath(const path: TPathD; pathIsClosed: Boolean; percentOffset, maxCtrlOffset: double; symmetric: Boolean): TPathD; var i, len, prev: integer; vec: TPointD; pl: TArrayOfDouble; unitVecs: TPathD; d, d1,d2: double; begin Result := nil; len := Length(path); if len < 3 then Exit; d := Max(0, Min(100, percentOffset))/200; if maxCtrlOffset <= 0 then maxCtrlOffset := MaxDouble; SetLength(Result, len *3 +1); prev := len-1; SetLength(pl, len); SetLength(unitVecs, len); for i := 0 to len -1 do begin pl[i] := Distance(path[prev], path[i]); unitVecs[i] := GetUnitVector(path[prev], path[i]); prev := i; end; Result[len*3] := path[0]; for i := 0 to len -1 do begin if i = len -1 then begin vec := GetAvgUnitVector(unitVecs[i], unitVecs[0]); d2 := pl[0]*d; end else begin vec := GetAvgUnitVector(unitVecs[i], unitVecs[i+1]); d2 := pl[i+1]*d; end; d1 := pl[i]*d; if symmetric then MakeSymmetric(d1, d2); if i = 0 then Result[len*3-1] := OffsetPoint(path[i], -vec.X * Min(maxCtrlOffset, d1), -vec.Y * Min(maxCtrlOffset, d1)) else Result[i*3-1] := OffsetPoint(path[i], -vec.X * Min(maxCtrlOffset, d1), -vec.Y * Min(maxCtrlOffset, d1)); Result[i*3] := path[i]; Result[i*3+1] := OffsetPoint(path[i], vec.X * Min(maxCtrlOffset, d2), vec.Y * Min(maxCtrlOffset, d2)); end; if not pathIsClosed then begin Result[1] := Result[0]; dec(len); Result[len*3-1] := Result[len*3]; SetLength(Result, Len*3 +1); end; end; And here's what it produces ... the path to smooth (black), the cubic bezier control path produced by GetSmoothPath() (blue) and the flattened cubic bezier path (2D graphics library of you choice required) (red). var TPathD path; begin path := MakePath([190,120, 260,270, 560,120, 190,490]); path := GetSmoothPath(path, true, 20, 0, false); path := ThirdParty2DGraphicsLibrary.FlattenCBezier(path); end; var TPathD path; begin path := MakePath([190,120, 260,270, 560,120, 190,490]); path := GetSmoothPath(path, true, 80, 0, false); path := ThirdParty2DGraphicsLibrary.FlattenCBezier(path); end; On 5/20/2022 at 8:45 PM, A.M. Hoornweg said: The problem: The polygons may never intersect ( Edit: The best way to avoid intersections is to make sure you have enough data points before generating your curves. Edited May 24 by angusj Note to explain the code focuses on simplicity not performance 3 Share this post Link to post
A.M. Hoornweg 126 Posted May 23 10 hours ago, Anders Melander said: I just remembered that Graphics32 has two examples which demonstrates interpolation: https://github.com/graphics32/graphics32/tree/master/Examples/Drawing/Curves Built into Graphics32 there's also the TCanvas32.CurveTo method which does cubic Bézier interpolation (4 control points) and the TCanvas32.ConicTo method which does quadratic Bézier interpolation (3 control points). Very interesting, thanks for the link! 1 Share this post Link to post
angusj 90 Posted July 17 (edited) On 5/23/2022 at 5:18 PM, angusj said: Here's my GetSmoothPath() routine. I've simplified this and uploaded it here: https://github.com/AngusJohnson/Image32/blob/0bb2e258f66c37f9eb263909480c473abf654f74/source/Img32.Extra.pas#L2044 Edited July 18 by angusj Edit: updated link 1 Share this post Link to post