Method roi_utils.get_mask() does not work for HFP patient set up

[Golpette]

Hi,

I tried to use this tool in a radiotherapy simulation in which the patient was in head first prone position (HFP). The mhd image I am using has the header:

ObjectType = Image
NDims = 3
BinaryData = True
BinaryDataByteOrderMSB = False
CompressedData = False
TransformMatrix = -1 0 0 0 -1 0 0 0 1
Offset = 250 250 -100
CenterOfRotation = 0 0 0
AnatomicalOrientation = LPI
ElementSpacing = 0.97656200000000004 0.97656200000000004 2.5
DimSize = 512 512 164
ElementType = MET_SHORT
ElementDataFile = ct.raw

Seems the issue is that the get_mask method is expecting the Origin to represent the centre of the voxel with the minimum coordinates, which is not the case for this HFP image.

I'm not comfortable enough with ITK / mhd / dicom to suggest a robust fix to this. But my hack solution for now is to determine a "direction" for each axis using GetDirection() (corresponding to the mhd TransformMatrix). I'm then using these +/-1 values to determine if the Origin field is representing a max or min position of the image. So I have:

(1) Replaced lines 517-518 of roi_utils.py to find the centre of the minimal voxel and the direction of the axes:

        orig = roimask.GetOrigin()  
        space = roimask.GetSpacing()   

        # directions +/- 1 for each axis       
        directions = list( roimask.GetDirection()*(1,1,1) )         

        # Find minimum corner (voxel centre)
        orig_min_corner = []
        for o,dim,s,direction in zip(orig,dims,space,directions):
            if direction>0:
                orig_min_corner.append(o)
            else:
                orig_min_corner.append(o-(dim-1)*s)

(2) Give this orig_min_corner to check "contained" in line 527:

for o,s,d,rmin,rmax in zip(orig_min_corner,space,dims,self.bb.mincorner,self.bb.maxcorner):

(3) Replaced 537-538 in case this min/max orientation issue ever happens with the z-axis too:

        if directions[2]>0:
            zmin = orig[2] - 0.5*space[2]
            zmax = orig[2] + (dims[2]-0.5)*space[2]
        else:
            zmin = orig[2] + (0.5-dims[2])*space[2]
            zmax = orig[2] + 0.5*space[2]

(4) Replaced 557-558 to give the xpoints and ypoints the correct ordering, using the 'directions' list above, before feeding them into the meshgrid:

        xpoints=np.linspace(orig[0],orig[0]+directions[0]*space[0]*dims[0],dims[0],False)
        ypoints=np.linspace(orig[1],orig[1]+directions[1]*space[1]*dims[1],dims[1],False)

The above seems to have solved my problem but I'm not comfortable enough with ITK/mhd/dicom images and their properties, or the rest of roi_utils.py, to know if this hack will generate problems elsewhere.

[djboersma]

Hi @Golpette , thanks for the report. I must admit that I only tested this code with HFS images and I agree that we should make sure that it also works with other anatomical orientations. I think I saw in some of the other gate tools that @tbaudier managed to ensure that they work with other orientations. The roi utils code is more messy than the other tools, unfortunately. We'll have look at your proposed fix. I'm currently very busy with other things, so I cannot promise that this is fixed very soon, but hopefully within a few weeks this should at least partly be addressed.

[Golpette]

Great, I'll also update the thread if I manage to do anything extra. I think using GetDirection() as above is fine if the matrix is always diagonal with values of +/-1 (which it should be for any standard image for radiotherapy). I initially wanted to just take the AnatomicalOrientation tag from the mhd file ("RAI", "LPI" etc) and generate the directions list from that, but as far as I can see ITK totally ignores that information.