# Get best fit sphere of humeral head

**URL:** <https://discourse.slicer.org/t/get-best-fit-sphere-of-humeral-head/22975>\
**Category:** Support\
**Tags:** segmentation, python\
**Created:** [April 15, 2022, 11:09am UTC](https://discourse.slicer.org/t/get-best-fit-sphere-of-humeral-head/22975 "2022-04-15T11:09:30Z")\
**Posts on this page:** 11\
**Page:** 1

<div class="post-metadata">

**Author:** ![osints](https://avatars.discourse-cdn.com/v4/letter/o/4af34b/32.png) [@osints](https://discourse.slicer.org/u/osints)\
**Post date:** [April 15, 2022, 11:09am UTC](https://discourse.slicer.org/t/get-best-fit-sphere-of-humeral-head/22975/1 "2022-04-15T11:09:30Z")

</div>

Hi,  
I am tyring to makes best fit sphere of humeral head.  
So i makes 3 land marks and then import above code to python interactor.  
but there is no sphere and some red text is shown…  
plase some help…  
thank you

 ![image](https://us1.discourse-cdn.com/flex002/uploads/slicer/original/3X/3/9/39abd41525623548dd32836cf746cbe89940c6e5.png)  
 ![image](https://us1.discourse-cdn.com/flex002/uploads/slicer/original/3X/3/2/32cbd449d2009566db530bd7ac5fe1a9416f12c3.png)

---

<div class="post-metadata">

**Author:** ![jamesobutler](https://sea2.discourse-cdn.com/flex002/user_avatar/discourse.slicer.org/jamesobutler/32/7511_2.png) [@jamesobutler](https://discourse.slicer.org/u/jamesobutler)\
**Post date:** [April 15, 2022, 11:38am UTC](https://discourse.slicer.org/t/get-best-fit-sphere-of-humeral-head/22975/4 "2022-04-15T11:38:34Z")

</div>

`slicer.util.arrayfromMarkupsControlPoints(markups)` has an undefined “markups” object. This variable must be a reference to your point list.

`markups=slicer.util.getNode("F")` is example code of getting a point list with a current name of “F” which is a common default name for it.

---

<div class="post-metadata">

**Author:** ![osints](https://avatars.discourse-cdn.com/v4/letter/o/4af34b/32.png) [@osints](https://discourse.slicer.org/u/osints)\
**Post date:** [April 16, 2022, 11:16am UTC](https://discourse.slicer.org/t/get-best-fit-sphere-of-humeral-head/22975/5 "2022-04-16T11:16:16Z")

</div>

Thank you…  
That`s works.  
I have some more questions.

1. Can i get center point of that sphere?
2. Is there any method that one plane move to crossing the point which made by question 1. (Center ponints of sphere)
3. Is there any method that calculate the distance between above 2 planes.

 ![image](https://us1.discourse-cdn.com/flex002/uploads/slicer/original/3X/e/c/ec1ab4cc637e23071b6867dd9f839b09b859ad9c.jpeg)

---

<div class="post-metadata">

**Author:** ![lassoan](https://sea2.discourse-cdn.com/flex002/user_avatar/discourse.slicer.org/lassoan/32/13_2.png) [@lassoan](https://discourse.slicer.org/u/lassoan)\
**Post date:** [April 16, 2022, 2:19pm UTC](https://discourse.slicer.org/t/get-best-fit-sphere-of-humeral-head/22975/6 "2022-04-16T14:19:29Z")

</div>

> [@osints](#):
>
> Can i get center point of that sphere?

Yes, sure. You have set it yourself in the [script](https://slicer.readthedocs.io/en/latest/developer_guide/script_repository.html#specify-a-sphere-by-multiple-control-points): `sphere.SetCenter(centerX, centerY, centerZ)`

> [@osints](#):
>
> Is there any method that one plane move to crossing the point which made by question 1. (Center ponints of sphere)

You can create a plane at any position. Add a new plane node, set its position to the sphere centerpoint, and set a normal.

> [@osints](#):
>
> Is there any method that calculate the distance between above 2 planes.

You can use basic computational geometry methods for computing distances, angles between points, lines, planes, project points, lines, determine intersections, etc. There are a couple of [examples in the script repository](https://slicer.readthedocs.io/en/latest/developer_guide/script_repository.html#measure-angle-between-two-markup-planes) to see how to get inputs and show outputs in Slicer, and you can learn the rest by typing your question into Google search (for example `python distance between two planes`) and reading through the search results. You will get both detailed mathematical explanation and complete Python implementations (typically few lines of numpy code). If you really get stuck then you can ask here, but it could be useful for you to spend some with reading and experimenting.

---

<div class="post-metadata">

**Author:** ![osints](https://avatars.discourse-cdn.com/v4/letter/o/4af34b/32.png) [@osints](https://discourse.slicer.org/u/osints)\
**Post date:** [April 17, 2022, 2:54am UTC](https://discourse.slicer.org/t/get-best-fit-sphere-of-humeral-head/22975/7 "2022-04-17T02:54:48Z")

</div>

thank you…  
I have some more question…  
I can`t find by google search…

1. Can i get center point of that sphere? → I`m trying to find center point of sphere aleady made by humeral head fiduciel not set a center point before making sphere.

2. I want make plane which are exactly parallel to one plane (Aleady made by GT foot print) and crossing center points of sphere.

3. can i subtract plane by overlapped with humerus bone defect?  
I want measure the depth of defect. image below  
 ![image](https://us1.discourse-cdn.com/flex002/uploads/slicer/original/3X/2/b/2b236a91ac75676c1f193cedc57bef267b64470e.jpeg)

---

<div class="post-metadata">

**Author:** ![lassoan](https://sea2.discourse-cdn.com/flex002/user_avatar/discourse.slicer.org/lassoan/32/13_2.png) [@lassoan](https://discourse.slicer.org/u/lassoan)\
**Post date:** [April 17, 2022, 3:29am UTC](https://discourse.slicer.org/t/get-best-fit-sphere-of-humeral-head/22975/8 "2022-04-17T03:29:52Z")

</div>

> [@osints](#):
>
> Can i get center point of that sphere? → I`m trying to find center point of sphere aleady made by humeral head fiduciel not set a center point before making sphere.

The [script](https://slicer.readthedocs.io/en/latest/developer_guide/script_repository.html#specify-a-sphere-by-multiple-control-points) computes the centerpoint position and stores it in `centerX, centerY, centerZ`.

> [@osints](#):
>
> I want make plane which are exactly parallel to one plane (Aleady made by GT foot print) and crossing center points of sphere.

The plane orientation is defined by the plane normal vector. You can create parallel planes by setting the same plane normal vector.

> [@osints](#):
>
> can i subtract plane by overlapped with humerus bone defect?  
> I want measure the depth of defect. image below

You can use Dynamic Modeler to cut a model with planes.

For making measurements you probably don’t need to actually cut a model, but it is enough to view a cross-section in a slice view. You can enable display of slice intersections and then rotate the slice intersections; or you can enable [reformat](https://slicer.readthedocs.io/en/latest/user_guide/user_interface.html#slice-view) widget to translate/rotate a slice plane in 3D.

---

<div class="post-metadata">

**Author:** ![gsirca](https://avatars.discourse-cdn.com/v4/letter/g/a587f6/32.png) [@gsirca](https://discourse.slicer.org/u/gsirca)\
**Post date:** [July 28, 2022, 6:10pm UTC](https://discourse.slicer.org/t/get-best-fit-sphere-of-humeral-head/22975/9 "2022-07-28T18:10:25Z")

</div>

Hello,  
I’m experiencing the same problem:

```auto
# Get markup node from scene
>>> pointListNode = slicer.util.getNode("F")
>>> 
>>> from scipy.optimize import least_squares
>>> import numpy
>>> 
>>> def fit_sphere_least_squares(x_values, y_values, z_values, initial_parameters, bounds=((-numpy.inf, -numpy.inf, -numpy.inf, -numpy.inf),(numpy.inf, numpy.inf, numpy.inf, numpy.inf))):
... """
... Source: https://github.com/thompson318/scikit-surgery-sphere-fitting/blob/master/sksurgeryspherefitting/algorithms/sphere_fitting.py
... Uses scipy's least squares optimisor to fit a sphere to a set
... of 3D Points
... :return: x: an array containing the four fitted parameters
... :return: ier: int An integer flag. If it is equal to 1, 2, 3 or 4, the
... solution was found.
... :param: (x,y,z) three arrays of equal length containing the x, y, and z
... coordinates.
... :param: an array containing four initial values (centre, and radius)
... """
... return least_squares(_calculate_residual_sphere, initial_parameters, bounds=bounds, method="trf", jac="3-point", args=(x_values, y_values, z_values))
... 
>>> def _calculate_residual_sphere(parameters, x_values, y_values, z_values):
... """
... Source: https://github.com/thompson318/scikit-surgery-sphere-fitting/blob/master/sksurgeryspherefitting/algorithms/sphere_fitting.py
... Calculates the residual error for an x,y,z coordinates, fitted
... to a sphere with centre and radius defined by the parameters tuple
... :return: The residual error
... :param: A tuple of the parameters to be optimised, should contain [x_centre, y_centre, z_centre, radius]
... :param: arrays containing the x,y, and z coordinates.
... """
... #extract the parameters
... x_centre, y_centre, z_centre, radius = parameters
... #use numpy's sqrt function here, which works by element on arrays
... distance_from_centre = numpy.sqrt((x_values - x_centre)**2 + (y_values - y_centre)**2 + (z_values - z_centre)**2)
... return distance_from_centre - radius
... 
>>> # Fit a sphere to the markups fidicual points
>>> markupsPositions = slicer.util.arrayFromMarkupsControlPoints("F")
Traceback (most recent call last):
  File "<console>", line 1, in <module>
  File "D:\gsirc\AppData\Local\NA-MIC\Slicer 4.11.20210226\bin\Python\slicer\util.py", line 1808, in arrayFromMarkupsControlPoints
    numberOfControlPoints = markupsNode.GetNumberOfControlPoints()
AttributeError: 'str' object has no attribute 'GetNumberOfControlPoints'
>>> import numpy as np
>>> # initial guess
>>> center0 = np.mean(markupsPositions, 0)
Traceback (most recent call last):
  File "<console>", line 1, in <module>
NameError: name 'markupsPositions' is not defined
>>> radius0 = np.linalg.norm(np.amin(markupsPositions,0)-np.amax(markupsPositions,0))/2.0
Traceback (most recent call last):
  File "<console>", line 1, in <module>
NameError: name 'markupsPositions' is not defined
>>> fittingResult = fit_sphere_least_squares(markupsPositions[:,0], markupsPositions[:,1], markupsPositions[:,2], [center0[0], center0[1], center0[2], radius0])
Traceback (most recent call last):
  File "<console>", line 1, in <module>
NameError: name 'markupsPositions' is not defined
>>> [centerX, centerY, centerZ, radius] = fittingResult["x"]
Traceback (most recent call last):
  File "<console>", line 1, in <module>
NameError: name 'fittingResult' is not defined
>>> 
>>> # Create a sphere using the fitted parameters
>>> sphere = vtk.vtkSphereSource()
>>> sphere.SetPhiResolution(30)
>>> sphere.SetThetaResolution(30)
>>> sphere.SetCenter(centerX, centerY, centerZ)
Traceback (most recent call last):
  File "<console>", line 1, in <module>
NameError: name 'centerX' is not defined
>>> sphere.SetRadius(radius)
Traceback (most recent call last):
  File "<console>", line 1, in <module>
NameError: name 'radius' is not defined
>>> sphere.Update()
>>> 
>>> # Add the sphere to the scene
>>> modelsLogic = slicer.modules.models.logic()
>>> model = modelsLogic.AddModel(sphere.GetOutput())
>>> model.GetDisplayNode().SetSliceIntersectionVisibility(True)
>>> model.GetDisplayNode().SetSliceIntersectionThickness(3)
>>> model.GetDisplayNode().SetColor(1,1,0)
>>> 
>>> 

```

Can anybody provide some help?

---

<div class="post-metadata">

**Author:** ![lassoan](https://sea2.discourse-cdn.com/flex002/user_avatar/discourse.slicer.org/lassoan/32/13_2.png) [@lassoan](https://discourse.slicer.org/u/lassoan)\
**Post date:** [July 28, 2022, 6:16pm UTC](https://discourse.slicer.org/t/get-best-fit-sphere-of-humeral-head/22975/10 "2022-07-28T18:16:33Z")

</div>

> [@gsirca](#):
>
> `markupsPositions = slicer.util.arrayFromMarkupsControlPoints("F")`

`F` is the a name of a node. `slicer.util.arrayFromMarkupsControlPoints` requires a node object as input. You can get the node object from a node name by calling `getNode()`:

`markupsPositions = slicer.util.arrayFromMarkupsControlPoints(slicer.util.getNode("F"))`

---

<div class="post-metadata">

**Author:** ![gsirca](https://avatars.discourse-cdn.com/v4/letter/g/a587f6/32.png) [@gsirca](https://discourse.slicer.org/u/gsirca)\
**Post date:** [July 30, 2022, 2:19pm UTC](https://discourse.slicer.org/t/get-best-fit-sphere-of-humeral-head/22975/11 "2022-07-30T14:19:20Z")

</div>

> [@lassoan](#):
>
> slicer.util.getNode(“F”)

Thank you for the help. It worked smoothly.

---

<div class="post-metadata">

**Author:** ![Bkeohane](https://avatars.discourse-cdn.com/v4/letter/b/c68b51/32.png) [@Bkeohane](https://discourse.slicer.org/u/Bkeohane)\
**Post date:** [December 7, 2022, 1:48pm UTC](https://discourse.slicer.org/t/get-best-fit-sphere-of-humeral-head/22975/13 "2022-12-07T13:48:23Z")

</div>

Hi, I succesfully extratced the coordinates and radius and I assume, they are in mm, right?

---

<div class="post-metadata">

**Author:** ![lassoan](https://sea2.discourse-cdn.com/flex002/user_avatar/discourse.slicer.org/lassoan/32/13_2.png) [@lassoan](https://discourse.slicer.org/u/lassoan)\
**Post date:** [December 7, 2022, 2:58pm UTC](https://discourse.slicer.org/t/get-best-fit-sphere-of-humeral-head/22975/14 "2022-12-07T14:58:09Z")

</div>

Yes, length unit is millimeters by default. It is always recommended to validate your entire workflow by testing on some object of known size.