PyESAPI icon indicating copy to clipboard operation
PyESAPI copied to clipboard

Calculate overlapping volume (OAR and PTV)

Open Shea1986 opened this issue 1 year ago • 0 comments

Hi, may I ask how to calculate the overlapping volume in pyesapi using this method (https://github.com/mtparagon5/ESAPI-Projects/blob/d497d7818013296e5fdffbe6cdafcb601e953ac2/Plugins/_CheckOverlap.cs)?

I see you have the ovh.py but it requires write access to patient objects and above algorithm from Matthew's github is the alternative way using C#. I converted C# into python by writing the following functions but unable to make out how to use VVector method in pyesapi - can you help?

def volume_overlap(structure1, structure2):
    # Initialize variables for overlap calculation
    volume_intersection = 0
    intersection_count = 0
    
    # Get bounds for both structures
    structure1_bounds = structure1.MeshGeometry.Bounds
    structure2_bounds = structure2.MeshGeometry.Bounds
    
    # Combine bounds of both structures
    combined_rect_bounds = {
        "min_x": min(structure1_bounds.X, structure2_bounds.X),
        "max_x": max(structure1_bounds.X + structure1_bounds.SizeX, structure2_bounds.X + structure2_bounds.SizeX),
        "min_y": min(structure1_bounds.Y, structure2_bounds.Y),
        "max_y": max(structure1_bounds.Y + structure1_bounds.SizeY, structure2_bounds.Y + structure2_bounds.SizeY),
        "min_z": min(structure1_bounds.Z, structure2_bounds.Z),
        "max_z": max(structure1_bounds.Z + structure1_bounds.SizeZ, structure2_bounds.Z + structure2_bounds.SizeZ),
    }

    # Set the resolution and range in each direction
    start_z = math.floor(combined_rect_bounds["min_z"] - 1)
    end_z = start_z + round(combined_rect_bounds["max_z"] - combined_rect_bounds["min_z"] + 2)
    start_x = math.floor(combined_rect_bounds["min_x"] - 1)
    end_x = start_x + round(combined_rect_bounds["max_x"] - combined_rect_bounds["min_x"] + 2)
    start_y = math.floor(combined_rect_bounds["min_y"] - 1)
    end_y = start_y + round(combined_rect_bounds["max_y"] - combined_rect_bounds["min_y"] + 2)

    # Check if one structure fully contains the other
    if structure1_bounds.Contains(structure2_bounds):
        volume_intersection = structure2.Volume
    elif structure2_bounds.Contains(structure1_bounds):
        volume_intersection = structure1.Volume
    else:
        # Iterate through the bounding box for overlap
        for z in np.arange(start_z, end_z, 0.5):
            for y in np.arange(start_y, end_y, 1.0):
                for x in np.arange(start_x, end_x, 1.0):
                    # Create a VVector point using the Eclipse API's VVector
                    point = VVector(x, y, z)
                    # Check if the point is inside both structures
                    if structure1.IsPointInsideSegment(point) and structure2.IsPointInsideSegment(point):
                        intersection_count += 1

        # Calculate volume intersection based on voxel counts
        volume_intersection = intersection_count * 0.001 * 0.5
    
    return volume_intersection

def percent_overlap(structure, volume_intersection):
    # Calculate percentage overlap
    percent_overlap = volume_intersection / structure.Volume
    return min(percent_overlap, 1)

Shea1986 avatar Oct 01 '24 21:10 Shea1986