VarianAPIs / VelocityEngine

Python and C# Clients for Velocity Scripting API
7 stars 5 forks source link

Rigid registration.getType() for deformable registration made with API #16

Closed mfizyczka closed 2 years ago

mfizyczka commented 2 years ago

When creating deformable registration in API I get the 'Rigid' type for it. It also seems Rigid while looking on Properties window in Velocity App. Vectors and grid clearly show that it is deformable registration.

aanghele commented 2 years ago

Can you please provide the block of code that produces this behavior so i can reproduce.

mfizyczka commented 2 years ago

Please find the code below. I also checked it in GUI and there is only one registration with rigid matrix (see attached PrtScr) GUI_registrationProperties .

registrationName = 'apiExtCont'
registration = regOps.createNewRegistration(registrationName)
orThrow(registration, regOps)
orThrow(engine.loadRegistration(registration.getVelocityId()))

overlapRegion = volOps.getVolumeOverlapRegionDICOM()

# RUN DEFORMABLE REGISTRATION
# deformableRegistrationSettitngs = velocity.DefaultBSplineDeformableRegistrationSettings
print('Performing deformable registration...')
deformableRegistrationSettitngs = createDeformableRegistrationSettingsObject(iRoiRegion = overlapRegion,
  iPrimaryStartLevel=-250, iPrimaryEndLevel=1700, iSecondaryStartLevel=-250, iSecondaryEndLevel=1700, 
  iNumberOfMultiResolutionLevels=5, iApplyBoundaryContinuityConstraints=True)
orThrow(regOps.performBsplineRegistrationDICOM(deformableRegistrationSettitngs), regOps)
print('done')

print(f'Registration type: {registration.getType()}')

Output:

Performing deformable registration...

done

Registration type: Rigid

I obtain the same result for iNumberOfMultiResolutionLevels = 3.

Function definition:

def createDeformableRegistrationSettingsObject(iRoiRegion = None,
  iPrimaryStartLevel=-250, iPrimaryEndLevel=1700, iSecondaryStartLevel=-250, iSecondaryEndLevel=1700,
  iPreProcessingMethod = 0, iNumberOfMultiResolutionLevels=3, iApplyBoundaryContinuityConstraints = False, iApplyTopologicalRegularizer = False,
  iGradientMagnitudeTolerance = 0.000000000000000000005, iGridCellSize = (5.0, 10.0, 15.0), iGridCellSizeType = 'n', 
  iMaximumNumberOfIterations = 30, iMaximumNumberOfConsecutiveOptimizerAttempts = 10, iMetricValuePercentageDifference = 0.0,
  iMinimumStepLength = 0.000001, iMaximumStepLength = 100.0, iNumberOfHistogramBins = 50, iRelaxationFactor = 0.9, iSamplesDenominator = 5,
  iTopologicalRegularizerDistanceLimitingCoefficient = velocity.VectorR3d(0.0)):  

  # Define Deformable Registration Settings
  bsplineSettings = velocity.BSplineDeformableRegistrationSettingsStructure()

  # if roiRegion is None set it to overlap:
  if iRoiRegion is not None: 
    bsplineSettings.roiStart[0] = iRoiRegion.start[0]
    bsplineSettings.roiStart[1] = iRoiRegion.start[1]
    bsplineSettings.roiStart[2] = iRoiRegion.start[2]
    bsplineSettings.roiEnd[0] = iRoiRegion.end[0]
    bsplineSettings.roiEnd[1] = iRoiRegion.end[1]
    bsplineSettings.roiEnd[2] = iRoiRegion.end[2]

  # setting intensity range as in HPTC protocol (for CBCT?); values in HU (200HU-1700HU selects bones)
  # notebook p.77
  bsplineSettings.primaryStartLevel = iPrimaryStartLevel
  bsplineSettings.primaryEndLevel = iPrimaryEndLevel
  bsplineSettings.secondaryStartLevel = iSecondaryStartLevel
  bsplineSettings.secondaryEndLevel = iSecondaryEndLevel

  bsplineSettings.preprocessingMethod = iPreProcessingMethod
  bsplineSettings.numberOfMultiResolutionLevels = iNumberOfMultiResolutionLevels # so each vector setting should be length 3
  bsplineSettings.applyBoundaryContinuityConstraints = velocity.BoolList([iApplyBoundaryContinuityConstraints]*iNumberOfMultiResolutionLevels)
  bsplineSettings.applyTopologicalRegularizer = velocity.BoolList([iApplyTopologicalRegularizer]*iNumberOfMultiResolutionLevels)
  bsplineSettings.gradientMagnitudeTolerance = velocity.DoubleList([iGradientMagnitudeTolerance]*iNumberOfMultiResolutionLevels)
  bsplineSettings.gridCellSize = velocity.VectorR3dList(( velocity.VectorR3d(iGridCellSize[0]), velocity.VectorR3d(iGridCellSize[1]), velocity.VectorR3d(iGridCellSize[2]) ))
  bsplineSettings.gridCellSizeType = velocity.CharList([ord(iGridCellSizeType)]*iNumberOfMultiResolutionLevels)
  bsplineSettings.maximumNumberOfIterations = velocity.IntList([iMaximumNumberOfIterations]*iNumberOfMultiResolutionLevels)
  bsplineSettings.maximumNumberOfConsecutiveOptimizerAttempts = velocity.IntList([iMaximumNumberOfConsecutiveOptimizerAttempts]*iNumberOfMultiResolutionLevels)
  bsplineSettings.metricValuePercentageDifference = velocity.DoubleList([iMetricValuePercentageDifference]*iNumberOfMultiResolutionLevels)
  bsplineSettings.minimumStepLength = velocity.DoubleList([iMinimumStepLength]*iNumberOfMultiResolutionLevels)
  bsplineSettings.maximumStepLength = velocity.DoubleList([iMaximumStepLength]*iNumberOfMultiResolutionLevels)
  bsplineSettings.numberOfHistogramBins = velocity.IntList([iNumberOfHistogramBins]*iNumberOfMultiResolutionLevels)
  bsplineSettings.relaxationFactor = velocity.DoubleList([iRelaxationFactor]*iNumberOfMultiResolutionLevels)

  bsplineSettings.samplesDenominator = velocity.IntList([iSamplesDenominator]*iNumberOfMultiResolutionLevels) # only 1/5 of pixels in ROI will be considered

  r3dZeroes = velocity.VectorR3d(0.0)
  bsplineSettings.topologicalRegularizerDistanceLimitingCoefficient = velocity.VectorR3dList([iTopologicalRegularizerDistanceLimitingCoefficient]*iNumberOfMultiResolutionLevels)

  return bsplineSettings
jakecobb commented 2 years ago

Registration algorithm results are only in memory until saved as shown in the pet_ct_registrations.py example. You need to call RegistrationOperations::saveRegistration to persist the changes. This returns an updated handle to the registration record that will reflect the type and name change (name due to versioning):

orThrow(regOps.performBsplineRegistrationDICOM(deformableRegistrationSettitngs), regOps)
print('done')

registration = regOps.saveRegistration()
orThrow(registration, regOps)
print(f'Registration type: {registration.getType()}')
mfizyczka commented 2 years ago

I am saving the results, and I checked this:

# RUN DEFORMABLE REGISTRATION
print('Performing deformable registration...')
deformableRegistrationSettitngs = createDeformableRegistrationSettingsObject(iRoiRegion = overlapRegion,
  iPrimaryStartLevel=-250, iPrimaryEndLevel=1700, iSecondaryStartLevel=-250, iSecondaryEndLevel=1700, 
  iNumberOfMultiResolutionLevels=3, iApplyBoundaryContinuityConstraints=False)
orThrow(regOps.performBsplineRegistrationDICOM(deformableRegistrationSettitngs), regOps)
print('done')

print(f'Registration type: {registration.getType()}')

# changes are just in memory, save changes to the databases
orThrow(regOps.saveRegistration(), regOps)

print(f'Registration type: {registration.getType()}')

Output:

Performing deformable registration...

done
Registration type: Rigid

Registration type: Rigid
jakecobb commented 2 years ago

The patientdata::Registration object only has the database info from the time it was fetched, it is not live in the sense of refreshing automatically. If you assign the return value with:

registration = regOps.saveRegistration()

Are you still seeing the type remain Rigid?

mfizyczka commented 2 years ago

Now this works! Indeed I was missing the assignment. I was just about to write it here.

Why I don't see this is GUI?

jakecobb commented 2 years ago

After saveRegistration you should be able to see the difference in the GUI. If you don't try closing and re-opening the patient in the GUI client, possibly it didn't refresh if you had the patient open already in the GUI when the script ran.