Determining closest lists near any selected point

Question

I'm plotting several lists of data (constant lines of density, temperature quality, etc. on a pressure-enthalpy diagram), none of which are fitted well by Interpolation or Fit. Since each table/line represents a constant density (in this case), I want to determine what density is at any selected point based off of the two nearest density lines (blue).

• Blue points are data points in rho lists.
• Orange points are all points from my function TotalClose which returns one point from each rho list.
• Green points are from my function PointClose.
• Largest black point is set with Manipulate controls.

The Locator is on a log scale so pt makes it linear for calculations:

Control[{{pt0,{2295,Log@145}},Locator}]

pt={First@pt0,Exp@Last@pt0};

In my attempt I found the point from each rho list closest to the point pt selected with the Locator:

TotalClose[list_]:=First@Nearest[list[[#]],pt]&/@Range@Length@list;

PointClose[list_]:=Nearest[TotalClose[list],pt,2];

Functions are evaluated at {rho1000,rho700,rho50,rho2,rho02}

The green points PointClose returned are not correct here. It needs to return points from the rho700 and rho50 lists. This would be easier if I could fit these lists to a function but they won't.

I was planning on using the lever rule to determine what the density would be at the selected point. I am stuck as to how to get which two rho lists are either immediately to the left and right or above and below pt.

I tried changing the DistanceFunction in Nearest with no luck:

 DistanceFunction->(((#1[[1]]-#2[[1]])^2+(Log@#1[[2]]-Log@#2[[2]])^2)^1/2&)

Here are the rho lists that make the blue constant density lines:

rho1000 = {{2674.112, 10000.}, {2442.182, 9000.}, {2211.362, 
    8000.}, {1980.958, 7000.}, {1750.036, 6000.}, {1517.27, 
    5000.}, {1280.666, 4000.}, {1036.931, 3000.}, {779.806, 
    2000.}, {493.584, 1000.}, {461.831, 900.}, {429.132, 
    800.}, {395.296, 700.}, {360.063, 600.}, {323.059, 
    500.}, {283.706, 400.}, {241.029, 300.}, {193.124, 200.}, {135.03,
     100.}, {120.986, 80.}, {105.405, 60.}, {87.389, 40.}, {64.549, 
    20.}, {48.572, 10.}, {44.443, 8.}, {39.658, 6.}, {34.028, 
    4.098}, {29.72, 2.984}, {25.44, 2.181}, {21.188, 1.692}, {16.964, 
    1.521}, {16.964, 0.05}};

rho700 = {{5000., 8801.33}, {4642.4, 8000.}, {4207.4, 7000.}, {3784.3,
     6000.}, {3371.92, 5000.}, {2968.32, 4000.}, {2570.4, 
    3000.}, {2173.14, 2000.}, {1767.73, 1000.}, {1726.18, 
    900.}, {1684.37, 800.}, {1642.24, 700.}, {1599.79, 
    600.}, {1556.95, 500.}, {1513.7, 400.}, {1469.97, 300.}, {1425.7, 
    200.}, {1380.8, 100.}, {1377.65, 93.070}};

rho50 = {{5000., 331.29}, {4673.64, 300.}, {3703.16, 200.}, {2815.85, 
    100.}, {2739.86, 91.84}, {2444, 71}, {2100, 50}, {1900, 
    37}, {1700, 26}, {1490, 15.5}, {1300, 9.5}, {1050, 4.1}, {820, 
    1.6}, {590, 0.5}, {425, 0.2}, {325, 0.1}, {250, 0.05}};

rho2 = {{5000., 13.07}, {4184.58, 10.}, {3693.09, 8.}, {3234.16, 
    6.}, {2802.61, 4.}, {2734.26, 3.68}, {2500, 3.3}, {2150, 
    2.5}, {1950, 2}, {1750, 1.5}, {1500, 1}, {1230, 0.6}, {1020, 
    0.38}, {815, 0.2}, {660, 0.115}, {450, 0.05}};

rho02 = {{5000., 1.31}, {4184.66, 1.}, {3692.28, 0.8}, {3232.73, 
    0.6}, {2803.20, 0.4}, {2626.43, 0.31}, {2500, 0.32}, {2250, 
    0.32}, {2000, 0.3}, {1840, 0.28}, {1660, 0.24}, {1480, 
    0.2}, {1250, 0.15}, {1000, 0.105}, {810, 0.08}, {550, 0.05}};

Answer 1

Update

I don't know what the "lever rule" is, even after visiting the Wikipedia page on it, but here's the routine wrapped up into a Module.

Clear[interpolateDensity]
interpolateDensity[pnt : {_, _}, fs : {{_InterpolatingFunction, _InterpolatingFunction} ..}, dens_List] /; Length@dens == Length@fs := Module[
  {mins}
  , mins = NMinimize[{#.# &@(pt - Through@#@t), 1. <= t <= #[[1]]["Domain"][[1, 2]]}, t] & /@ fs
  ; Dot @@ Transpose@#/Plus @@ #[[All, 1]] &@ Sort[Transpose[{mins[[All, 1]], dens}]][[1 ;; 2]]
 ]

This calculates the density by finding the nearest points on the nearest curves and using the relative distances as weights for the densities of those curves. Usage:

Clear[pt, fs, dens]
pt = {3604, 693};
fs = Interpolation /@ Transpose@# & /@ {rho02, rho2, rho50, rho700, rho1000};
dens = {0.2, 2, 50, 700, 1000};
interpolateDensity[pt, fs, dens]
(* 18.7113 *)

Alternatively, here is a routine that returns the list of the nearest points along with the densities corresponding to those curves (which you can then use to get the density via the lever rule?):

nearestPoints[pnt : {_, _}, fs : {{_InterpolatingFunction, _InterpolatingFunction} ..}, dens_List] /; Length@dens == Length@fs := Module[
  {mins}
  , mins = NMinimize[{#.# &@(pt - Through@#@t), 1. <= t <= #[[1]]["Domain"][[1, 2]]}, t] & /@ fs
  ; Sort[MapThread[{#2[[1]], Through@#1@t /. #2[[2]], #3} &, {fs, mins, dens}]][[1 ;; 2]][[All, 2 ;; 3]]
  ]

Usage:

nearestPoints[pt, fs, dens]
(* {{{3658.84, 195.126}, 50}, {{3606.93, 7.63443}, 2}} *)

Original Post

For starters, I will not re-paste in your data, but I will define a list of densities:

dens = {0.2, 2, 50, 700, 1000};

Here's something that I think will work. We first carefully Interpolate the data by interpolating the x and y lists separately to create a parametric interpolating function, like so:

fs = Interpolation /@ Transpose@# & /@ {rho02, rho2, rho50, rho700, rho1000};

These functions will be defined on the intervals 1 to Length@rho:

domains = Flatten /@ Through[Transpose[fs][[1]]@"Domain"]
(* {{1., 16.}, {1., 16.}, {1., 17.}, {1., 19.}, {1., 32.}} *)

We can plot the InterpolatingFunctions to see that they look nice:

p1 = Show[
 ListLogPlot[{rho02, rho2, rho50, rho700, rho1000}]
 , ParametricPlot[Evaluate[{#1, Log@#2} & @@ (Through@#@t)], Evaluate@Flatten@{t, #[[1]]["Domain"]}] & /@ fs
]

We now go through the process of minimizing the Euclidean distance between some point

pt = {RandomReal[{0, 5000}], RandomReal[{10, 10^3}]}
(* {3604.77,693.335} *)

and the curves using the code:

mins = NMinimize[{#.# &@(pt - Through@#@t), 1. <= t <= #[[1]]["Domain"][[1, 2]]}, t] & /@ fs;

We are minimizing the square of the Euclidean distance here, but putting the square root around #.# doesn't seem to spit out any errors, so it should work just as well. Obviously, if you think some other distance function is more relevant, then here's where to modify the code. For illustration purposes, we extract the points on the curves as

pts = MapThread[Through@#1@t /. Last@#2 &, {fs, mins}]
(* {{3605.06,0.763082},{3607.7,7.63772},{3659.63,195.214},{1921.07,1376.9},{649.59,1537.91}} *)

and plot everything together:

Show[
  p1
  , ListLogPlot[{pt}, PlotStyle -> {Black, PointSize[0.02]}]
  , ListLogPlot[pts, PlotStyle -> {Red, PointSize[0.02]}]
 ]

Finally, we choose the closest two points, keeping their distances and their densities:

weights = Sort[Transpose[{mins[[All, 1]], dens}]][[1 ;; 2]]
(* {{251135., 50}, {470189., 2}} *)

And finally, we linearly interpolate the densities by weighting the densities by the relative distances:

density = Dot @@ Transpose@#/Plus @@ #[[All, 1]] &@weights
(* 18.7116 *)

All of this can be automated into a neat little Module, and of course the code as written isn't the nicest to look at, but I think this works! Take it and run with it.

Answer 2

fs = Interpolation /@ {Mean /@ GatherBy[rho1000, First], rho700, rho50, rho2, rho02};



Show[
 ListLogPlot[{Mean /@ GatherBy[rho1000, First], rho700, rho50, rho2, 
   rho02}
  , PlotRange -> All
  , PlotTheme -> "Scientific"
  ]
 , Sequence @@ Table[
   With[{f = Evaluate[fs[[k]]]},
    LogPlot[f[t], Evaluate[{t, Sequence @@ First[f["Domain"]]}]]
    ], {k, Length[fs]}]
 ]

minpf[pnt_, func_] := Block[{r, x, xr, y},
  {r, xr} = NMinimize[
    {
     Norm[{pnt, {t, func[t]}}]
     , Between[t, func["Domain"]]
     }, t];
  x = Part[xr, 1, 2];
  {r, x, func[x]}
  ]

minpf[pnt_, func_] := Block[{r, x, xr, y},
  {r, xr} = NMinimize[
    {
     Norm[{pnt, {t, func[t]}}]
     , Between[t, func["Domain"]]
     }, t];
  x = Part[xr, 1, 2];
  {r, x, func[x]}
  ]

findpnt[pnt_] := Rest@First@MinimalBy[Map[minpf[pnt, #] &, fs], First]

findpnt[{100, 1}]
(* {16.964, 0.7855} *)

Answer 3

Following my OP comment, the first thing that worked. (Optimized? Heck, it barely works. I considered using the rho-sets as their labels, but that just seemed like nuking flies, so string labels.)

Module[{rhos, logRhos, minX, maxX, minY, maxY},
  rhos = {rho1000, rho700, rho50, rho2, rho02};
  logRhos = {#[[1]], Log[10, #[[2]]]} & /@ # & /@ rhos;
  {minX, maxX} = MinMax[(Join @@ logRhos)[[All, 1]]];
  {minY, maxY} = MinMax[(Join @@ logRhos)[[All, 2]]];
  {lrs1000, lrs700, lrs50, lrs2, lrs02} = (# - {minX, minY})/({maxX, maxY} - {minX, minY}) & /@ # & /@ logRhos;
  lrsPointToCurve = Association[Sequence[
    (# -> "lrs1000" &) /@ lrs1000,
    (# -> "lrs700" &) /@ lrs700,
    (# -> "lrs50" &) /@ lrs50,
    (# -> "lrs2" &) /@ lrs2,
    (# -> "lrs02" &) /@ lrs02]];
  logRhosScaled = Join[{lrs1000, lrs700, lrs50, lrs2, lrs02}];
  pointToFirstLRSPoint = Nearest[Join @@ logRhosScaled];
  knockout["lrs1000"] = Join[{lrs700, lrs50, lrs2, lrs02}];
  knockout["lrs700"] = Join[{lrs1000, lrs50, lrs2, lrs02}];
  knockout["lrs50"] = Join[{lrs1000, lrs700, lrs2, lrs02}];
  knockout["lrs2"] = Join[{lrs1000, lrs700, lrs50, lrs02}];
  knockout["lrs02"] = Join[{lrs1000, lrs700, lrs50, lrs2}];
  pointToSecondLRSPoint[pt_] := Block[{},
    Nearest[Join @@ knockout[
      lrsPointToCurve[pointToFirstLRSPoint[pt][[1]]]]][pt]];
]

ListPlot[logRhosScaled]
(* Point plot.  Happily, where they're expected to be. *)

VoronoiMesh[Join @@ logRhosScaled]
(* Pretty blue cellulation. *)

pointToFirstLRSPoint[{0.5, 0.5}]
lrsPointToCurve[pointToFirstLRSPoint[{0.5, 0.5}][[1]]]
pointToSecondLRSPoint[{0.5, 0.5}]
lrsPointToCurve[pointToSecondLRSPoint[{0.5, 0.5}][[1]]]
(* {{0.48706, 0.594656}} *)
(* lrs50 *)
(* {{0.559026, 0.359004}} *)
(* lrs2 *)

This strictly meets the request : go from an arbitrary point to the two nearest curves. The same sort of idea (use an Association[]) will let us go from logRhoScaled points to original points if that was actually intended.