[Insight-users] [Insight-users]itkMutualInformationImageToImageMetric.txx

[duboisda@tele.ucl.ac.be]

Some optimisation methods using the MutualInformationMetric need the
derivative of the MI. The slow part of the derivative is the line
   derivative += ( derivB - (*aditer) ) * weight;

You can accelerate it by 3 if you replace

    for( aiter = m_SampleA.begin(), aditer = sampleADerivatives.begin();
      aiter != aend; ++aiter, ++aditer )
        {

	valueFixed = ( (*biter).FixedImageValue - (*aiter).FixedImageValue ) /
	  m_FixedImageStandardDeviation;
	valueFixed = m_KernelFunction->Evaluate( valueFixed );

	valueMoving = ( (*biter).MovingImageValue - (*aiter).MovingImageValue ) /
	  m_MovingImageStandardDeviation;
	valueMoving = m_KernelFunction->Evaluate( valueMoving );

	weightMoving = valueMoving / dDenominatorMoving;
  	weightJoint = valueMoving * valueFixed / dDenominatorJoint;

	weight = ( weightMoving - weightJoint );
	weight *= (*biter).MovingImageValue - (*aiter).MovingImageValue;

  	derivative += ( derivB - (*aditer) ) * weight;

      } // end of sample A loop
    } // end of sample B loop

by

    totalweight = 0;

    for( aiter = m_SampleA.begin(), aditer = sampleADerivatives.begin();
      aiter != aend; ++aiter, ++aditer )
        {

	valueFixed = ( (*biter).FixedImageValue - (*aiter).FixedImageValue ) /
	  m_FixedImageStandardDeviation;
	valueFixed = m_KernelFunction->Evaluate( valueFixed );

	valueMoving = ( (*biter).MovingImageValue - (*aiter).MovingImageValue ) /
	  m_MovingImageStandardDeviation;
	valueMoving = m_KernelFunction->Evaluate( valueMoving );

	weightMoving = valueMoving / dDenominatorMoving;
  	weightJoint = valueMoving * valueFixed / dDenominatorJoint;

	weight = ( weightMoving - weightJoint );
	weight *= (*biter).MovingImageValue - (*aiter).MovingImageValue;

	totalweight += weight;
	derivative -= (*aditer) * weight;

      } // end of sample A loop

    derivative += derivB * totalweight;

    } // end of sample B loop

And you have to initialize totalweight at the beginning of the function as
double.
derivB is a constant array in the loop. Add it takes a long time. Just add
the weight during the A loop and multiply derivB at the end is faster....

   aloys