4
$\begingroup$

Context

enter image description here

We are trying to predict the formation of discs of intermediate mass black hole orbiting super massive black holes in galactic centres. This involves finding calorific curves in energy and angular momentum space.

Problem

I am trying to identify sets of curves from a (very large) sets of points in order to produce a lighter figure.

For instance I want to be able to redo such figure, but with less points (without loosing the poorly sampled lines).

enter image description here

The challenge is that these points are found by some complex optimisation routine and are over-numerous in places. So by lighter figure, I mean remove sets of point very close to each other while keeping more sparsely sampled points in order to draw curves through them.

Question

How can one achieve this simultaneous downsampling and curve matching procedure?


Attempt

Let me define a toy problem as follows: let me produce two curves:

data = Flatten[{Table[{x, Sin[x^2/15]}, {x, 0, 20, 0.01}],
    Table[{x, 2 + Cos[x]}, {x, 0, 20, 0.1}]}, 1];

and draw random points from those.

pts = RandomChoice[data, Length[data]];

ListPlot[pts]

enter image description here

Note that on purpose the sampling is not the same for the two curves.


I found on this site and in the documentation a couple of partial solutions to this problem, but none does both the curve identification and downsampling properly.

The first method relies on Nearest (from the documentation) but does neither job (interpolation or downsampling) particularly well:

p = Nearest[pts];
pts1 = p[#, 2][[2]] & /@ pts;
Graphics[{Line[{##}] & @@@ Transpose[{pts, pts1}]}, Axes -> True, 
 AspectRatio -> 2/3]

enter image description here

From this answer I tried also ListCurvePathPlot which seems like a good start but does not down sample.

pl = ListCurvePathPlot[pts, Mesh -> All, 
  MeshStyle -> Directive[ColorData[97][2], PointSize[Large]],
  AspectRatio -> 1]

enter image description here

Finally FindShortestTour doesnt do as good a job, but lets me downsample by brute force (which is clearly sub optimal).


ListLinePlot[pts[[FindShortestTour[pts][[2]][[;; ;; 20]]]], 
 PlotMarkers -> {Automatic, 10}, InterpolationOrder -> 1]

enter image description here


I am optimistic that this question has a Mathematica answer since it can be done by hand relatively easily (?)

$\endgroup$
4
  • 2
    $\begingroup$ I guess, the results of FindShortestTour should be post-processed with additional criteria. Jump from one fork to another is conjugated with significant change of the distance. So, you need go along the sorted points looking on the value of current distance and cut the array for segments. Further you can downsample them as you wish. $\endgroup$
    – Rom38
    Jul 7, 2021 at 9:35
  • $\begingroup$ @Rom38 agreed; may be a mixture between Nearest and FindShortestTour? Or Interpolation? It would be best to find as generic a solution to this problem as possible. $\endgroup$
    – chris
    Jul 7, 2021 at 9:45
  • 1
    $\begingroup$ In a sense, your example is "easy" as the two curves do not intersect. Intersection makes things much harder. $\endgroup$
    – A.G.
    Jul 11, 2021 at 18:37
  • $\begingroup$ @A.G. Absolutely. I have added another question (mathematica.stackexchange.com/q/250907/1089) to try and address this more general problem. $\endgroup$
    – chris
    Jul 12, 2021 at 15:45

2 Answers 2

1
$\begingroup$

Here is a way to downsample based on the distance between points. I use Manipulate to adjust the threshold that reduces the number of points. As the $x$ and $y$ axes are on different scales it is probably better to scale before taking euclidean distances, as is done in the second example.

SeedRandom[0];
data = Flatten[{Table[{x, Sin[x^2/15]}, {x, 0, 20, 0.01}], 
    Table[{x, 2 + Cos[x]}, {x, 0, 20, 0.1}]}, 1];
n = Length[data];
pts = RandomChoice[data, n];

(* Eliminate points based on euclidean distance *)
Manipulate[
 p1 = DeleteDuplicates[pts, 
   EuclideanDistance[#1, #2] <= minDistance1 &];
 ListPlot[p1, 
  PlotLabel -> 
   Row[{"Number of points: ", Length@p1, " out of ", n}]],
 {{minDistance1, .1}, 0, 0.25}]

(* Eliminate points based on scaled euclidean distance *)
xRange = With[{x = First /@ pts}, Max@x - Min@x];
yRange = With[{y = Last /@ pts}, Max@y - Min@y];
scaledDistance[{x1_, y1_}, {x2_, y2_}] := 
  Sqrt[(x2 - x1)^2/xRange^2 + (y2 - y1)^2/yRange^2 ];

Manipulate[
 p2 = DeleteDuplicates[pts, 
   scaledDistance[#1, #2] <= minDistance2 &];
 ListPlot[p2, 
  PlotLabel -> 
   Row[{"Number of points: ", Length@p2, " out of ", n}]],
 {{minDistance2, .02}, 0.01, 0.05}]

(* Now add interpolation *)
Manipulate[
 p3 = DeleteDuplicates[pts, 
   scaledDistance[#1, #2] <= minDistance3 &];
 Show[
  ListCurvePathPlot[p3, Mesh -> All, 
   MeshStyle -> Directive[ColorData[97][2], PointSize[Large]], 
   AspectRatio -> 1, 
   PlotLabel -> Row[{"Number of points: ", Length@p3, " out of ", n}]],
  ListPlot[p3, PlotStyle -> {Black, PointSize[0.01]}]],
 {{minDistance3, .02}, 0.01, 0.05}]

enter image description here

enter image description here

enter image description here

$\endgroup$
1
  • $\begingroup$ Thank you for your answer. $\endgroup$
    – chris
    Jul 7, 2021 at 21:06
1
$\begingroup$

Edit: Just noticed @Rom38 essentially mentioned the same approach in the comments above, credit goes to them.

This can be solved as two separate questions: namely separating the points into 'curves' and then down-sampling a list of points.

For the first question, I use FindShortestTour to find a single continuous curve, and then split it up according to EuclideanDistance jumps. I suspect this might need some fine-tuning for different use-cases, but seems to work well for the case of OP:

data = Flatten[{Table[{x, Sin[x^2/15]}, {x, 0, 20, 0.01}],Table[{x, 2 + Cos[x]}, {x, 0, 20, 0.1}]}, 1];
pts = RandomSample[data];

tour = Last[FindShortestTour[pts]];
tourPts = Extract[pts, List /@ tour];
peaks = Ordering[EuclideanDistance @@@ Partition[tourPts, 2, 1], -2];
{firstCurve, secondCurve} = TakeDrop[RotateLeft[tourPts, peaks[[1]]], Abs[Subtract @@ peaks]];

For the second question, there's many related questions on this site, e.g. here and here. I'll use the accepted solution from the first link for demonstration purposes.

np[f_][u_, dt_] := u + dt/Norm[f'[u]]
equallySpacedPts[pts_, dt_] := With[{bsf = BSplineFunction[pts]},
  bsf /@ Most[NestWhileList[np[bsf][#, dt] &, 0, # < 1 &]]]

equallySpacedPts[#, 0.25] & /@ {firstCurve, secondCurve} // ListPlot

enter image description here

$\endgroup$
3
  • $\begingroup$ Apologies, I had changed your pts definition as follows pts = RandomSample[data] to avoid duplicates and forgot to include in answer above. Note that if you want to only use some of the data, you can use something like pts = RandomSample[data, Floor[Length[data] 0.5]] for a duplicate-free random list of half the points. $\endgroup$ Jul 8, 2021 at 19:44
  • $\begingroup$ Thank you again .Would the method work with an unknown number of curves? Or a known one but not necessarily 2? $\endgroup$
    – chris
    Jul 8, 2021 at 19:48
  • $\begingroup$ The first part would need to be generalized slightly to identify multiple curves using EuclideanDistance jumps, but as long as FindShortestTour finds a reasonable path (e.g. the curves don't intersect), I reckon it should be easy to do. $\endgroup$ Jul 8, 2021 at 19:52

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.