Finden des Bereichs einer Überlappung zwischen Kurven (Python)

Question

Ist es möglich, die Fläche der Überlappung zweier Kurven zu berechnen? Ich habe hier zwei Antworten gefunden, aber sie sind in R geschrieben, mit denen ich nicht vertraut bin. Oder kämpfen, um sie in Python umzuwandeln.

Bereich zwischen den beiden Kurven und Überlappungsbereich zwischen zwei Kurven suchen

Zum Beispiel für einen bestimmten Datensatz mit definierten x, y-Punkten. (x1, y1, x2, y2)

Ich kann die Fläche jeder Kurve ermitteln mit:

np.trapz

Es ist jedoch eine Herausforderung, nur die Überlappung zu erhalten, und ich habe keine Lösung gefunden, die ich zeigen könnte. Alle Anleitungen oder mathematischen Formeln werden geschätzt.

Answer 1

Dies kann also mit dem Shapely-Modul in Python erfolgen.

Verbinden Sie zunächst die beiden Kurven, um ein sich selbst schneidendes Polygon zu erstellen (siehe Code unten).

unary_union()Wenn Sie dann die Funktion von Shapely verwenden, teilen Sie das komplexe Polygon in einfache Polygone auf, ermitteln die Fläche jedes Polygons und summieren sie, um die Gesamtfläche der beiden Kurven zu ermitteln.

Code unten gezeigt:

    import numpy as np
    from shapely.geometry import LineString
    from shapely.ops import unary_union, polygonize
    
    avg_coords = [(0.0, 0.0), (4.872117, 2.29658), (5.268545, 2.4639225), (5.664686, 2.6485724), (6.059776, 2.8966842), (6.695151, 3.0986626), (7.728006, 3.4045217), (8.522297, 3.652668), (9.157002, 3.895031), (10.191483, 4.1028132), (10.827622, 4.258638), (11.38593, 4.2933016), (11.86478, 4.3048816), (12.344586, 4.258769), (12.984073, 4.2126703), (13.942729, 4.1781383), (14.58212, 4.137809), (15.542498, 3.99943), (16.502588, 3.878359), (17.182951, 3.7745714), (18.262657, 3.6621647), (19.102558, 3.567045), (20.061789, 3.497897), (21.139917, 3.4806826), (22.097425, 3.5153809), (23.65388, 3.5414772), (24.851482, 3.541581), (26.04966, 3.507069), (27.72702, 3.463945), (28.925198, 3.429433), (29.883854, 3.3949006), (31.08246, 3.3344274), (31.92107, 3.317192), (33.716183, 3.3952322), (35.63192, 3.4213595), (37.427895, 3.4474766), (39.343628, 3.473604), (41.49874, 3.508406), (43.773468, 3.5518723), (46.287716, 3.595359), (49.28115, 3.6302335), (52.633293, 3.6997545), (54.30922, 3.7431688), (55.8651, 3.8038807), (58.738773, 3.8387446), (60.893887, 3.8735466), (63.647655, 3.9170544), (66.760704, 3.960593), (68.79663, 3.9607692), (70.23332, 3.986855), (72.867905, 3.995737), (75.38245, 4.0219164), (77.778656, 3.9615464), (79.337975, 3.8145657), (80.41826, 3.6675436), (80.899734, 3.5204697), (81.62059, 3.38207), (82.34045, 3.3042476), (83.30039, 3.1918304), (84.38039, 3.062116), (84.50359, 2.854434), (83.906364, 2.7591898), (83.669716, 2.586092), (83.43435, 2.3351095), (83.19727, 2.1879735), (82.84229, 1.9283267), (82.48516, 1.7984879), (81.65014, 1.5993768), (80.454544, 1.4781193), (79.13962, 1.3308897), (77.944595, 1.1750168), (76.39001, 1.0364205), (74.59633, 0.87184185), (71.60447, 0.741775), (70.04903, 0.6551017), (58.3, 0.0)]
    model_coords = [(0.0, 0.0), (0.6699889, 0.18807), (1.339894, 0.37499), (2.009583, 0.55966), (2.67915, 0.74106), (3.348189, 0.91826), (4.016881, 1.0904), (4.685107, 1.2567), (5.359344, 1.418), (6.026172, 1.5706), (6.685472, 1.714), (7.350604, 1.8508), (8.021434, 1.9803), (8.684451, 2.0996), (9.346408, 2.2099), (10.0066, 2.311), (10.66665, 2.4028), (11.32436, 2.4853), (11.98068, 2.5585), (12.6356, 2.6225), (13.29005, 2.6775), (13.93507, 2.7232), (14.58554, 2.7609), (15.23346, 2.7903), (15.87982, 2.8116), (16.52556, 2.8254), (17.16867, 2.832), (17.80914, 2.8317), (18.44891, 2.825), (19.08598, 2.8124), (19.72132, 2.7944), (20.35491, 2.7713), (20.98673, 2.7438), (21.61675, 2.7121), (22.24398, 2.677), (22.86939, 2.6387), (23.49297, 2.5978), (24.1147, 2.5548), (24.73458, 2.51), (25.3526, 2.464), (25.96874, 2.4171), (26.58301, 2.3697), (27.1954, 2.3223), (27.80491, 2.2751), (28.41354, 2.2285), (29.02028, 2.1829), (29.62512, 2.1384), (30.22809, 2.0954), (30.82917, 2.0541), (31.42837, 2.0147), (32.02669, 1.9775), (32.62215, 1.9425), (33.21674, 1.9099), (33.80945, 1.8799), (34.40032, 1.8525), (34.98933, 1.8277), (35.5765, 1.8058), (36.16283, 1.7865), (36.74733, 1.7701), (37.33002, 1.7564), (37.91187, 1.7455), (38.49092, 1.7372), (39.06917, 1.7316), (39.64661, 1.7285), (40.22127, 1.7279), (40.79514, 1.7297), (41.36723, 1.7337), (41.93759, 1.7399), (42.50707, 1.748), (43.07386, 1.7581), (43.63995, 1.7699), (44.20512, 1.7832), (44.76772, 1.7981), (45.3295, 1.8143), (45.88948, 1.8318), (46.44767, 1.8504), (47.00525, 1.8703), (47.55994, 1.8911), (48.11392, 1.9129), (48.6661, 1.9356), (49.21658, 1.959), (49.76518, 1.9832), (50.31305, 2.0079), (50.85824, 2.033), (51.40252, 2.0586), (51.94501, 2.0845), (52.48579, 2.1107), (53.02467, 2.1369), (53.56185, 2.1632), (54.09715, 2.1895), (54.63171, 2.2156), (55.1634, 2.2416), (55.69329, 2.2674), (56.22236, 2.2928), (56.74855, 2.3179), (57.27392, 2.3426), (57.7964, 2.3668), (58.31709, 2.3905), (58.83687, 2.4136), (59.35905, 2.4365), (59.87414, 2.4585), (60.38831, 2.4798), (60.8996, 2.5006), (61.40888, 2.5207), (61.91636, 2.5401), (62.42194, 2.5589), (62.92551, 2.577), (63.42729, 2.5945), (63.92607, 2.6113), (64.42384, 2.6275), (64.91873, 2.643), (65.4127, 2.658), (65.90369, 2.6724), (66.39266, 2.6862), (66.87964, 2.6995), (67.36373, 2.7123), (67.84679, 2.7246), (68.32689, 2.7364), (68.80595, 2.7478), (69.28194, 2.7588), (69.756, 2.7695), (70.22709, 2.7798), (70.69707, 2.7898), (71.16405, 2.7995), (71.62902, 2.809), (72.0919, 2.8183), (72.55277, 2.8273), (73.01067, 2.8362), (73.46734, 2.845), (73.92112, 2.8536), (74.37269, 2.8622), (74.82127, 2.8706), (75.26884, 2.8791), (75.71322, 2.8875), (76.15559, 2.8958), (76.59488, 2.9042), (77.03304, 2.9126), (77.46812, 2.921), (77.90111, 2.9294), (78.33199, 2.9379), (78.75986, 2.9464), (79.18652, 2.955), (79.60912, 2.9637), (80.03049, 2.9724), (80.44985, 2.9811), (80.86613, 2.99), (81.2802, 2.9989), (81.69118, 3.0078), (82.10006, 3.0168), (82.50674, 3.0259), (82.91132, 3.035), (83.31379, 3.0441), (83.71307, 3.0533), (84.10925, 3.0625), (84.50421, 3.0717), (84.8961, 3.0809), (85.28577, 3.0901), (85.67334, 3.0993), (86.05771, 3.1085), (86.43989, 3.1176), (86.81896, 3.1267), (87.19585, 3.1358), (87.57063, 3.1448), (87.94319, 3.1537), (88.31257, 3.1626), (88.67973, 3.1713), (89.04372, 3.18), (89.40659, 3.1886), (89.7652, 3.197), (90.12457, 3.2053), (90.47256, 3.2135), (90.82946, 3.2216), (91.17545, 3.2295), (91.52045, 3.2373), (91.86441, 3.2449), (92.20641, 3.2524), (92.54739, 3.2597), (92.88728, 3.2669), (93.21538, 3.2739), (93.55325, 3.2807), (93.87924, 3.2874), (94.20424, 3.2939), (94.52822, 3.3002), (94.85012, 3.3064), (95.16219, 3.3123), (95.48208, 3.3182), (95.79107, 3.3238), (96.09807, 3.3293), (96.40505, 3.3346), (96.71003, 3.3397), (97.01401, 3.3447), (97.31592, 3.3496), (97.60799, 3.3542), (97.90789, 3.3587), (98.19686, 3.3631), (98.48386, 3.3673), (98.77085, 3.3714), (99.05574, 3.3753), (99.32983, 3.3791), (99.6127, 3.3828), (99.8837, 3.3863), (100.1538, 3.3897), (100.4326, 3.393), (100.6897, 3.3961), (100.9566, 3.3991), (101.2215, 3.402), (101.4756, 3.4048), (101.7375, 3.4075), (101.9885, 3.4101), (102.2385, 3.4126), (102.4875, 3.4149), (102.7354, 3.4172), (102.9714, 3.4194), (103.2163, 3.4214), (103.4493, 3.4234), (103.6823, 3.4253), (103.9133, 3.4271), (104.1433, 3.4288), (104.3712, 3.4304), (104.5882, 3.4319), (104.8141, 3.4333), (105.0291, 3.4346), (105.2421, 3.4358), (105.4541, 3.437), (105.6651, 3.438), (105.8751, 3.439), (106.083, 3.4399), (106.28, 3.4407), (106.4759, 3.4414), (106.6699, 3.442), (106.8629, 3.4425), (107.0549, 3.443), (107.2458, 3.4433), (107.4249, 3.4435), (107.6128, 3.4437), (107.7897, 3.4438), (107.9647, 3.4437), (108.1387, 3.4436), (108.3116, 3.4433), (108.4737, 3.443), (108.6436, 3.4426), (108.8027, 3.4421), (108.9706, 3.4414), (109.1265, 3.4407), (109.2814, 3.4399), (109.4255, 3.439), (109.5784, 3.4379), (109.7195, 3.4368), (109.8694, 3.4356), (110.0084, 3.4342), (110.1454, 3.4328), (110.2813, 3.4313), (110.4162, 3.4296), (110.5403, 3.4279), (110.6722, 3.426), (110.7932, 3.424), (110.9132, 3.422), (111.0322, 3.4198), (111.1492, 3.4175), (111.2651, 3.4151), (111.3701, 3.4127), (111.483, 3.4101), (111.585, 3.4074), (111.686, 3.4046), (111.786, 3.4017), (111.884, 3.3987), (111.9809, 3.3956), (112.0669, 3.3924), (112.1608, 3.3891), (112.2448, 3.3857), (112.3268, 3.3822), (112.4078, 3.3786), (112.4867, 3.3749), (112.5548, 3.3711), (112.6317, 3.3672), (112.6978, 3.3632), (112.7726, 3.3591), (112.8356, 3.3549), (112.8975, 3.3506), (112.9476, 3.3462), (113.0076, 3.3417), (113.0655, 3.3372), (113.1125, 3.3325), (113.1584, 3.3278), (113.2024, 3.3229), (113.2464, 3.318), (113.2884, 3.313), (113.3283, 3.3079), (113.3584, 3.3027), (113.3963, 3.2974), (113.4233, 3.292), (113.4492, 3.2865), (113.4742, 3.281), (113.4972, 3.2753), (113.5201, 3.2696), (113.5312, 3.2638), (113.5501, 3.2579), (113.5591, 3.2519), (113.5661, 3.2459), (113.5721, 3.2397), (113.577, 3.2335), (113.5809, 3.2272), (113.573, 3.2208), (113.5749, 3.2143), (113.5649, 3.2077), (113.5539, 3.2011), (113.5409, 3.1944), (113.5278, 3.1876), (113.5128, 3.1807), (113.4967, 3.1737), (113.4697, 3.1667), (113.4418, 3.1596), (113.4227, 3.1524), (113.3917, 3.145), (113.3597, 3.1375), (113.3266, 3.1298), (113.2827, 3.1218), (113.2475, 3.1136), (113.2016, 3.1051), (113.1635, 3.0964), (113.1155, 3.0873), (113.0655, 3.0779), (113.0144, 3.0683), (112.9525, 3.0583), (112.8994, 3.048), (112.8345, 3.0373), (112.7793, 3.0264), (112.7123, 3.0152), (112.6453, 3.0037), (112.5763, 2.9919), (112.5063, 2.9798), (112.4352, 2.9674), (112.3533, 2.9548), (112.2801, 2.9419), (112.1952, 2.9287), (112.1102, 2.9153), (112.034, 2.9017), (111.9361, 2.8879), (111.8481, 2.8739), (111.7581, 2.8597), (111.667, 2.8453), (111.5661, 2.8307), (111.473, 2.816), (111.3689, 2.801), (111.2639, 2.786), (111.1579, 2.7708), (111.0509, 2.7555), (110.9428, 2.74), (110.8239, 2.7245), (110.7138, 2.7088), (110.5928, 2.6931), (110.4709, 2.6772), (110.3578, 2.6613), (110.2338, 2.6453), (110.1087, 2.6292), (109.9826, 2.613), (109.8457, 2.5968), (109.7176, 2.5805), (109.5787, 2.5642), (109.4496, 2.5478), (109.3086, 2.5314), (109.1666, 2.5149), (109.0236, 2.4984), (108.8806, 2.4819), (108.7355, 2.4653), (108.5905, 2.4488), (108.4434, 2.4322), (108.2865, 2.4155), (108.1384, 2.3989), (107.9794, 2.3822), (107.8195, 2.3655), (107.6684, 2.3488), (107.5063, 2.3321), (107.3374, 2.3156), (107.1744, 2.2989), (107.0104, 2.2822), (106.8442, 2.2654), (106.6683, 2.2487), (106.5012, 2.232), (106.3242, 2.2152), (106.1452, 2.1985), (105.9662, 2.1818), (105.7862, 2.165), (105.6052, 2.1483), (105.4232, 2.1316), (105.2402, 2.1149), (105.0572, 2.0981), (104.8721, 2.0814), (104.6772, 2.0647), (104.492, 2.048), (104.295, 2.0313), (104.098, 2.0147), (103.9, 1.998), (103.701, 1.9813), (103.502, 1.9647), (103.301, 1.948), (103.1, 1.9314), (102.899, 1.9148), (102.6959, 1.8982), (102.483, 1.8816), (102.2789, 1.865), (102.0649, 1.8484), (101.8588, 1.8318), (101.6428, 1.8153), (101.4268, 1.7988), (101.2098, 1.7822), (100.9918, 1.7657), (100.7728, 1.7492), (100.5538, 1.7328), (100.3338, 1.7163), (100.1128, 1.6999), (99.89169, 1.6834), (99.65978, 1.667), (99.43769, 1.6506), (99.20477, 1.6343), (98.98066, 1.6179), (98.74665, 1.6016), (98.51164, 1.5852), (98.27574, 1.5689), (98.04964, 1.5527), (97.81264, 1.5364), (97.57562, 1.5202), (97.33752, 1.5039), (97.08962, 1.4877), (96.8506, 1.4716), (96.61061, 1.4554), (96.37051, 1.4393), (96.12058, 1.4232), (95.87949, 1.4071), (95.62759, 1.391), (95.38547, 1.375), (95.13258, 1.359), (94.88946, 1.343), (94.63548, 1.3271), (94.38145, 1.3111), (94.12645, 1.2952), (93.87144, 1.2793), (93.61545, 1.2635), (93.35946, 1.2477), (93.10343, 1.2319), (92.84642, 1.2161), (92.58843, 1.2004), (92.33042, 1.1846), (92.07232, 1.169), (91.8034, 1.1533), (91.54331, 1.1377), (91.2744, 1.1221), (91.0133, 1.1065), (90.7434, 1.091), (90.48229, 1.0755), (90.21139, 1.0601), (89.9493, 1.0446), (89.67728, 1.0292), (89.40428, 1.0139), (89.13137, 0.99855), (88.86826, 0.98325), (88.59427, 0.96799), (88.32026, 0.95277), (88.04527, 0.93758), (87.77126, 0.92242), (87.4972, 0.90731), (87.21732, 0.89222), (86.94719, 0.87718), (86.66711, 0.86217), (86.3773, 0.8472), (86.10719, 0.83227), (85.82721, 0.81738), (85.5472, 0.80252), (85.26721, 0.7877), (84.9872, 0.77292), (84.7071, 0.75819), (84.41721, 0.74349), (84.1371, 0.72883), (83.84721, 0.71421), (83.5671, 0.69963), (83.27721, 0.68509), (82.99711, 0.6706), (82.70711, 0.65615), (82.41721, 0.64173), (82.1371, 0.62736), (81.8471, 0.61304), (81.55722, 0.59875), (81.27709, 0.58451), (80.98712, 0.57031), (80.697, 0.55616), (80.39711, 0.54205), (80.10722, 0.52798), (79.8271, 0.51396), (79.53701, 0.49999), (79.23711, 0.48605), (78.9471, 0.47217), (78.65701, 0.45833), (78.3571, 0.44453), (78.06712, 0.43078), (77.77701, 0.41708), (77.4771, 0.40343), (77.18701, 0.38982), (76.8871, 0.37626), (76.59711, 0.36274), (76.30701, 0.34928), (76.0071, 0.33586), (75.7169, 0.32249), (75.4071, 0.30917), (75.11701, 0.29589), (74.8171, 0.28267), (74.52701, 0.26949), (74.22711, 0.25636), (73.937, 0.24329), (73.63691, 0.23026), (73.3271, 0.21728), (73.03699, 0.20436), (72.73712, 0.19148), (72.4469, 0.17865), (72.13712, 0.16588), (71.84701, 0.15315), (71.547, 0.14048), (71.24701, 0.12786), (70.947, 0.11528), (70.64701, 0.10277), (70.3471, 0.090298), (70.05691, 0.077883), (69.74712, 0.06552), (69.457, 0.05321), (69.1569, 0.040952), (68.84709, 0.028747), (68.557, 0.016595), (68.25701, 0.0)]
    
    polygon_points = [] #creates a empty list where we will append the points to create the polygon
    
    for xyvalue in avg_coords:
        polygon_points.append([xyvalue[0],xyvalue[1]]) #append all xy points for curve 1
    
    for xyvalue in model_coords[::-1]:
        polygon_points.append([xyvalue[0],xyvalue[1]]) #append all xy points for curve 2 in the reverse order (from last point to first point)
    
    for xyvalue in avg_coords[0:1]:
        polygon_points.append([xyvalue[0],xyvalue[1]]) #append the first point in curve 1 again, to it "closes" the polygon
    
    avg_poly = [] 
    model_poly = []
    
    for xyvalue in avg_coords:
        avg_poly.append([xyvalue[0],xyvalue[1]]) 
    
    for xyvalue in model_coords:
        model_poly.append([xyvalue[0],xyvalue[1]]) 
    
    
    line_non_simple = LineString(polygon_points)
    mls = unary_union(line_non_simple)
    
    Area_cal =[]
    
    for polygon in polygonize(mls):
        Area_cal.append(polygon.area)
        print(polygon.area)# print area of each section 
        Area_poly = (np.asarray(Area_cal).sum())
        
    print(Area_poly)#print combined area

Answer 2

Stellen Sie Ihre Überlappungsbereiche nach Möglichkeit als Polygone dar. Von dort aus kann die Polygonfläche durch eine bemerkenswert präzise Formel berechnet werden, wie auf der Website von Paul Bourke erläutert .

Angenommen, (x [i], y [i]), i = 0, ..., N sind die Polygonscheitelpunkte mit (x [0], y [0]) = (x [N], y [N. ]), so dass das Polygon geschlossen ist und konsistent alle im Uhrzeigersinn oder alle gegen den Uhrzeigersinn. Dann ist die Gegend

area = |0.5 * sum_i (x[i] * y[i+1] - x[i+1] * y[i])|

wobei die Summe über i = 0, ..., N-1 geht. Dies gilt auch für nicht konvexe Polygone. Diese Formel ist im Wesentlichen das gleiche Prinzip, wie ein Planimeter die Fläche einer beliebigen zweidimensionalen Form misst, ein Sonderfall des Greenschen Theorems .

Answer 3

Wenn Ihre Funktionen tatsächlich "Funktion" sind, was bedeutet, dass keine vertikalen Linien die Funktionen mehr als einmal schneiden, ist das Finden der Überlappungen die Frage des Findens von Nullen.

import numpy as np
import matplotlib.pyplot as plt

dx = 0.01
x = np.arange(-2, 2, dx)
f1 = np.sin(4*x)
f2 = np.cos(4*x)

plt.plot(x, f1)
plt.plot(x, f2)

eps = 1e-1; # threshold of intersection points.
df = f1 - f2
idx_zeros = np.where(abs(df) <= eps)[0]

area = 0
for i in range(len(idx_zeros) - 1):
    idx_left = idx_zeros[i]
    idx_rite = idx_zeros[i+1]
    area += abs(np.trapz(df[idx_left:idx_rite])) * dx

• Ich habe angenommen, dass Bereiche als positiv angesehen werden.
• Der analytische Wert für das von mir verwendete Beispiel ist

ausreichend nahe am berechneten Wert ( area=2.819). Natürlich können Sie dies verbessern, wenn Ihre Gitter feiner und der Schwellenwert epskleiner sind.