validationsamplegenerator.txx

DTMK软件开发包,此为开源软件,是一款很好的医学图像开发资源.

      for ( i = 0 ; i < m_NumberOfMeasurements ; i++ )
        {
          pixel = iter.Get() ;

          m_ImageMeans[i] += (double) pixel[i] ;

          if ( pixel[i] < m_ImageMins[i] )
            {
              m_ImageMins[i] = pixel[i] ;
            }

          if ( pixel[i] > m_ImageMaxes[i] )
            {
              m_ImageMaxes[i] = pixel[i] ;
            }
        }
      ++iter ;
    }

  for ( i = 0 ; i < m_NumberOfMeasurements ; i++ )
    {
      m_ImageMeans[i] /= (double) sampleSize ;
      std::cout << "DEBUG: image means[" << i << "] = " << m_ImageMeans[i]
                << std::endl ;
      m_MeasurementsLowerBound[i] = m_ImageMins[i] ;
      m_MeasurementsUpperBound[i] = m_ImageMaxes[i] ;
    }

  double temp ;
  m_ImageStandardDeviations.Fill(0.0) ;
  iter.GoToBegin() ;

  while (!iter.IsAtEnd())
    {
      for ( i = 0 ; i < m_NumberOfMeasurements ; i++ )
        {
          pixel = iter.Get() ;

          temp = (double) pixel[i] - m_ImageMeans[i] ;
          temp *= temp ;
          m_ImageStandardDeviations[i] += temp ;
        }
      ++iter ;
    }

  for ( i = 0 ; i < m_NumberOfMeasurements ; i++ )
    {
      m_ImageStandardDeviations[i] = 
        vcl_sqrt(m_ImageStandardDeviations[i] / (double) sampleSize) ;
      std::cout << "DEBUG: image standard deviations[" << i << "] = " 
                << m_ImageStandardDeviations[i]
                << std::endl ;
    }
}

template< class TImage, class TClassMaskImage, class TVectorImage >
void
ValidationSampleGenerator< TImage, TClassMaskImage, TVectorImage >
::Normalize()
{
  unsigned int i ;
  typename TVectorImage::PixelType pixel ;

  itk::ImageRegionIteratorWithIndex< TVectorImage >
    iter(m_VectorImage, m_VectorImage->GetLargestPossibleRegion()) ;

  while ( !iter.IsAtEnd() )
    {
      pixel = iter.Get() ;

      for ( i = 0 ; i < m_NumberOfMeasurements ; i++ )
        { 
          pixel[i] = 
            (MeasurementType) ((pixel[i] - m_ImageMeans[i]) / 
            m_ImageStandardDeviations[i] * m_NormalizationScale) ;
        }

      iter.Set(pixel) ;
      ++iter ;
    }
}

template< class TImage, class TClassMaskImage, class TVectorImage >
void
ValidationSampleGenerator< TImage, TClassMaskImage, TVectorImage >
::GenerateListSample()
{
  unsigned int i ;
  int classIndex ;
  
  m_MeasurementsLowerBound.Fill(itk::NumericTraits< MeasurementType >::max()) ;
  m_MeasurementsUpperBound.Fill(itk::NumericTraits< MeasurementType >::min()) ;

  m_Sample->SetImage(m_VectorImage) ;

  itk::ImageRegionIteratorWithIndex< TClassMaskImage >  
    m_iter( m_ClassMaskImage, m_ClassMaskImage->GetLargestPossibleRegion()) ;
  
  ImageSampleType::Iterator s_iter = m_Sample->Begin() ;
  ImageSampleType::Iterator s_end = m_Sample->End() ;

  MeasurementVectorType temp ;
  unsigned int numberOfClasses = m_UniqueClassLabels.size() ;
  while ( s_iter != s_end )
    {
      classIndex = this->GetClassIndex(m_iter.Get()) ;
//       std::cout << "DEBUG: class index = " << classIndex 
//                 << " labels size = " << m_UniqueClassLabels.size()
//                 << std::endl ;
      if ( classIndex == -1 )
        {
          m_UniqueClassLabels.push_back(m_iter.Get()) ;
          numberOfClasses = m_UniqueClassLabels.size() ;
          classIndex = m_UniqueClassLabels.size() - 1 ;
          m_ClassMeans.resize(m_UniqueClassLabels.size()) ;
          m_ClassMeans[classIndex].Fill(0.0) ;
          m_ClassCovariances.resize(m_UniqueClassLabels.size()) ;
          m_ClassCovariances[classIndex].Fill(0.0) ;
          m_ClassSizes.resize(m_UniqueClassLabels.size()) ;
          m_ClassSizes[classIndex] = 0 ;
        }

      temp = s_iter.GetMeasurementVector() ;

      if ( std::find(m_SelectedClassLabels.begin(), m_SelectedClassLabels.end(), 
                    m_iter.Get()) != m_SelectedClassLabels.end() )
        {
          // selected class
          m_Subsample->AddInstance(s_iter.GetInstanceIdentifier()) ;

          // add the class label
          m_ClassLabels[s_iter.GetInstanceIdentifier()] = m_iter.Get() ;

          // update the sample space bound
          for ( i = 0 ; i < m_NumberOfMeasurements ; i++ )
            {
              if ( temp[i] < m_MeasurementsLowerBound[i] )
                {
                  m_MeasurementsLowerBound[i] = temp[i] ;
                }
              
              if ( temp[i] > m_MeasurementsUpperBound[i] )
                {
                  m_MeasurementsUpperBound[i] = temp[i] ;
                }
            }
        }

      for ( i = 0 ; i < m_NumberOfMeasurements ; i++ )
        {
          // used as a temporary sum storage
          m_ClassMeans[classIndex][i] += temp[i] ;
        }

      m_ClassSizes[classIndex] += 1 ;

      ++m_iter ;
      ++s_iter ;
    }

  // calculate the means
  for ( i = 0 ; i < numberOfClasses ; i++ )
    {
      m_ClassMeans[i] /= (double) m_ClassSizes[i] ;
      std::cout << " class labels = " << m_UniqueClassLabels[i] 
                << " mean = " << m_ClassMeans[i] 
                << " size = " << m_ClassSizes[i] << std::endl ;
    }

  m_iter.GoToBegin() ;
  s_iter = m_Sample->Begin() ;

  unsigned int row, col ;
  MeanType diff ;
  // calculate the covariances
  while (s_iter != s_end)
    {
      classIndex = this->GetClassIndex(m_iter.Get()) ;
      temp = s_iter.GetMeasurementVector() ;

      for (i = 0 ; i < m_NumberOfMeasurements ; i++)
        {
          diff[i] = temp[i] - m_ClassMeans[classIndex][i] ;
        }

      for ( row = 0; row < m_NumberOfMeasurements ; row++)
        {
          for ( col = 0; col < row + 1 ; col++)
            {
              m_ClassCovariances[classIndex].GetVnlMatrix()(row,col) += 
                diff[row] * diff[col] ;
            }
        }
      ++m_iter ;
      ++s_iter ;
    }

  // fills the upper triangle using the lower triangle  

  for ( i = 0 ; i < numberOfClasses ; i++ )
    {
      for (row = 1 ; row < m_NumberOfMeasurements ; row++)
        {
          for (col = 0 ; col < row ; col++)
            {
              m_ClassCovariances[i].GetVnlMatrix()(col, row) = 
                m_ClassCovariances[i].GetVnlMatrix()(row, col) ;
            } 
        }
      m_ClassCovariances[i].GetVnlMatrix() /= m_ClassSizes[i] ;
    }
}



template< class TImage, class TClassMaskImage, class TVectorImage >
void
ValidationSampleGenerator< TImage, TClassMaskImage, TVectorImage >
::GenerateHistogram()
{
  unsigned int i ;
  HistogramType::SizeType histogramSize ;
  HistogramType::MeasurementVectorType lowerBound ;
  HistogramType::MeasurementVectorType upperBound ;
  
  for ( i = 0 ; i < m_NumberOfMeasurements ; i++ )
    {
      histogramSize[i] = 
        (unsigned long) (m_MeasurementsUpperBound[i] - 
                         m_MeasurementsLowerBound[i] + 1) ;
      lowerBound[i] = m_MeasurementsLowerBound[i] - 0.5 ;
      upperBound[i] = m_MeasurementsUpperBound[i] + 0.5 ;
    }
  
  // creats equal size bins ;
  m_Histogram->Initialize(histogramSize, lowerBound, upperBound) ;

  // imports the subsample to histogram
  m_Importer->SetListSample(m_Subsample.GetPointer()) ;
  m_Importer->SetHistogram(m_Histogram.GetPointer()) ;
  time_t begin = clock() ;
  m_Importer->Run() ;
  time_t end = clock() ;
  std::cout << "DEBUG: histogram genererated in " 
            << double(end - begin) /CLOCKS_PER_SEC 
            << " seconds" << std::endl ;
}

template< class TImage, class TClassMaskImage, class TVectorImage >
void
ValidationSampleGenerator< TImage, TClassMaskImage, TVectorImage >
::GenerateWeightedCentroidKdTree()
{
  m_TreeGenerator->SetSample(m_Subsample) ;
  m_TreeGenerator->SetBucketSize(m_KdTreeBucketSize) ;
  time_t begin = clock() ;
  m_TreeGenerator->GenerateData() ;
  time_t end = clock() ;
  m_KdTree = m_TreeGenerator->GetOutput() ;
  std::cout << "DEBUG: k-d tree genererated in " 
            << double(end - begin) /CLOCKS_PER_SEC 
            << " seconds" << std::endl ;
}

#endif