如何画出点集构成的光滑闭合曲面
mathematica吧
全部回复
仅看楼主
level 8
jiaoeyushushu 楼主
SeedRandom[2]
blot[smoothness_: 20, points_Integer: 10] :=
With[{fun = Exp[-smoothness #.#] &,
pts = RandomReal[1, {points, 2}]},
RegionPlot[
Total[fun[# - {x, y}] & /@ pts] > .5, {x, -.5, 1.5}, {y, -.5, 1.5},
Frame -> False, PlotStyle -> Black, BoundaryStyle -> Black]]
img = ImageAdjust@
DistanceTransform[
SelectComponents[
Binarize@ImagePad[ImageCrop@ColorNegate@Rasterize[blot[]], 10],
"Count", -1]];
pts = Catenate[
Select[Cases[
ListPointPlot3D[
13 ImageData@
ImageTake[
ColorConvert[Blur[img, 25], "Grayscale"], {1, -1,
10}, {1, -1, 10}]], {x_Real, y_, z_},
Infinity], #[[3]] > 2 &] /. {x_, y_,
z_} :> {{x, y, -z - 3}, {x, y, z}}];
ListPointPlot3D[pts, AspectRatio -> Automatic]
点集为pts,构成一个光滑闭合曲面,谢谢各位大佬了。
2020年02月19日 09点02分 1
level 8
jiaoeyushushu 楼主
前排@无影东瓜 求看[勉强][乖][太开心]
2020年02月19日 09点02分 2
吧务
level 10
看着这点集可以分成俩单值曲面,分别作插值然后Plot3D应该就行了吧
2020年02月19日 13点02分 4
关键是闭合[咦][勉强]
2020年02月20日 02点02分
@jiaoeyushushu 那在球坐标下插值行不?
2020年02月20日 05点02分
吧务
level 12
Plot3D截取,然后插值把中间空隙补上就行
SeedRandom[2]
blot[smoothness_: 20, points_Integer: 10] :=
With[{fun = Exp[-smoothness #.#] &,
pts = RandomReal[1, {points, 2}]},
RegionPlot[
Total[fun[# - {x, y}] & /@ pts] > .5, {x, -.5, 1.5}, {y, -.5, 1.5},
Frame -> False, PlotStyle -> Black, BoundaryStyle -> Black]]
img = ImageAdjust@
DistanceTransform[
SelectComponents[
Binarize@ImagePad[ImageCrop@ColorNegate@Rasterize[blot[]], 10],
"Count", -1]];
mat = MapIndexed[Append[#2, #1] &,
13 ImageData@
ImageTake[
ColorConvert[Blur[img, 25], "Grayscale"], {1, -1, 10}, {1, -1,
10}], {2}];
{nx, ny} = Dimensions[mat][[1 ;; 2]];
pts = Select[Catenate[mat], #[[3]] > 2 &];
funUpper = Interpolation[Catenate[mat], InterpolationOrder -> 2];
center = Mean[pts][[1 ;; 2]];
middlePoints = Table[Module[{\[Theta]rList, bound},
bound =
Cases[Normal@
ContourPlot[funUpper[x, y] == z, {x, 1, nx}, {y, 1, ny},
PlotPoints -> 50], Line[bounds_] :> bounds, \[Infinity]][[1]];
\[Theta]rList =
CoordinateTransform[
"Cartesian" -> "Polar", # - center & /@ bound][[
All, {-1, 1}]] // Sort;
\[Theta]rList =
Join[# - {2 \[Pi],
0} & /@ \[Theta]rList, \[Theta]rList, # + {2 \[Pi],
0} & /@ \[Theta]rList] // DeleteDuplicates // Sort;
Thread[{z, N@#,
Interpolation[\[Theta]rList, Method -> "Spline"][#]} &@
Range[-1.2 \[Pi], 1.2 \[Pi], \[Pi]/18]]
], {z, 2, 2.3, 0.1}
];
middlePoints =
Join[middlePoints, MapAt[-# - 3 &, middlePoints, {All, All, 1}]];
funMiddle =
Interpolation[Flatten[middlePoints, 1], InterpolationOrder -> 3];
With[{opts =
Sequence[Mesh -> None, PlotPoints -> 50, BoundaryStyle -> None]},
Show[
Plot3D[funUpper[x, y], {x, 1, nx}, {y, 1, ny},
RegionFunction -> Function[{x, y, z}, z > 2], opts],
Plot3D[-funUpper[x, y] - 3, {x, 1, nx}, {y, 1, ny},
RegionFunction -> Function[{x, y, z}, -z - 3 > 2], opts],
ParametricPlot3D[{center[[1]] +
funMiddle[z, \[Theta]] Cos[\[Theta]],
center[[2]] + funMiddle[z, \[Theta]] Sin[\[Theta]], z}, {z, -5.1,
2.1}, {\[Theta], -\[Pi], \[Pi]}, opts],
ListPointPlot3D[Join[pts, MapAt[-# - 3 &, pts, {All, 3}]],
PlotStyle -> PointSize[Large]],
BoxRatios -> Automatic, PlotRange -> All
]
]
2020年02月20日 10点02分 5
哇,谢谢了,我明天仔细看一下。[太开心],看来你写的也挺麻烦的,我之前以为要用B条样插值,结果没画出来。[吐舌][捂嘴笑]
2020年02月20日 14点02分
1