Basándome en mi respuesta a otra pregunta sobre ajuste de círculos, la adaptación es muy sencilla.
En lugar las listas x
e y
que usaba el problema original para contener las coordenadas de los puntos, en esta ocasión usaremos arbol.x
y arbol.y
, siendo arbol
un dataframe pandas que contiene los datos de un árbol (y que tiene al menos dos columnas llamadas x
e y
).
El siguiente código muestra cómo operar (he ejecutado todos los comandos mostrados en celdas de un Jupyter Notebook):
Bajarse los datos proporcionados por la usuaria
!wget "https://www.dropbox.com/s/mc7njkif3ql78ts/arbol1.txt?dl=1" -o /dev/null -O arbol1.txt
Leerlos en un dataframe pandas y ver que tienen buen aspecto
# Leer los datos con pandas
import pandas as pd
df = pd.read_table("arbol1.txt", header=0)
df.head()
x y z0 z arbol
0 23.215 73.319 1.210 1068.3 1
1 23.209 73.332 1.210 1068.3 1
2 23.203 73.344 1.209 1068.3 1
3 23.208 73.321 1.397 1068.5 1
4 23.207 73.321 1.369 1068.4 1
Quedarse con los datos del "arbol 1" y pintarlos para ver que siguen una forma más o menos circular:
# Quedarse con un arbol, por ejemplo el 3
arbol = df[df.arbol == 1]
ax = arbol.plot(x="x", y="y", kind="scatter")
ax.set_aspect('equal')
Utilizar el código de la otra pregunta, adaptado a este caso, para calcular el centro y radio del círculo:
# Hacer el ajuste (basado en otra respuesta mía)
# https://es.stackoverflow.com/a/166354/7123
from scipy import optimize as opt
import numpy as np
def calc_R(x, y, xc, yc):
""" Calculate the distance of each 2D points from the center (xc, yc) """
return np.sqrt((x-xc)**2 + (y-yc)**2)
def f_2(c, x, y):
""" Calculate the algebraic distance between the 2D points and the mean circle centered at c=(xc, yc) """
Ri = calc_R(x, y, *c)
return Ri - Ri.mean()
# Coordenadas del baricentro (estimación inicial del centro del círculo)
x_m = arbol.x.mean()
y_m = arbol.y.mean()
barycenter = x_m, y_m
# Aplicamos un algoritmo para aproximar el centro de la circunferencia.
circle_center, ier = opt.leastsq(f_2, barycenter, args=(arbol.x, arbol.y))
xc, yc = circle_center
print("Centro del círculo:", xc, yc)
# El radio encontrado
r = calc_R(arbol.x, arbol.y, *circle_center).mean()
print("Radio del árbol:", r)
Centro del círculo: 22.98550108357253 73.229096837513
Radio del árbol: 0.2332831108455863
Finalmente podemos pintar encima de la gráfica anterior el círculo que hemos encontrado. Así:
from matplotlib.patches import Circle
from matplotlib.collections import PatchCollection
ax = arbol.plot(x="x", y="y", kind="scatter")
ax.plot(x=xc, y=yc, mark="*", color="red")
ax.add_collection(PatchCollection([Circle((xc, yc), r)], color="orange", alpha=0.2))
ax.set_aspect('equal')
ax.set_xlim((xc-r*1.1, xc+r*1.1))
ax.set_ylim((yc-r*1.1, yc+r*1.1))
.7z
que está vacío.