// Author:        Bill Christmas
// Date:          Jan 1999
// Descripton:   Calculates motion vectors for a region from a pair of images using LMS fit
//
// Synopsis:
//   exLMSRegion <image1> <image2> <region_image> [-r region_label]
//
// In the region image, all pixels of each region are set to the label for that
// region, starting from 1.  Label 0 is reserved for "not a region" areas.


#include "amma/Option.hh"
#include "amma/DP/FileFormatIO.hh"
#include "amma/Pair.hh"
#include "amma/Image/IPConvSym.hh"
#include "amma/Motion/LMSRegionMatch.hh"

int main (int argc, char **argv) 
{
  // process command-line parameters
  OptionC opt(argc, argv);
  FilenameC ImageFileName1 (opt.String("", NULL, "Input image 1"));
  FilenameC ImageFileName2 (opt.String("", NULL, "Input image 2"));
  FilenameC RegionName (opt.String("", NULL, "Image of regions"));
  IntT region_id (opt.Int("r", 0, "Region I.D. (default: all regions)")); 
  opt.CompulsoryArgs(3);
  opt.Check();

  // load pair of images
  PairC<NumImageC<RealT> > image, filtered;
  StdIO::Load(ImageFileName1, image[0]);
  StdIO::Load(ImageFileName2, image[1]);
  // filter images with antialias filter
  SArray1dC<RealT> mask(3); mask[0]=1; mask[1]=0.5; mask[2]=0.25;
  IPConvSymC<RealT, RealT> filter(mask, TwoD);
  for (UIntT i=0; i<=1; ++i) {
    filtered[i] = filter.Apply(image[i]);
  }
  // load image of regions & centre it
  NumImageC<IntT> regions;
  StdIO::Load(RegionName, regions);
  Index2dC offset(image[0].Rectangle().Center()-regions.Rectangle().Center());
  regions.ShiftRowIndexes(offset.Row());
  regions.ShiftColIndexes(offset.Col());
  // compute motion
  LMSRegionMatchC motion(filtered);
  if (opt.IsOnCommandLine("r"))
    cout << "Motion vector: " << motion.Estimate(regions, region_id) << '\n';
  else {
    IntT region_id(1);
    do {
      Vector2dC vec (motion.Estimate(regions, region_id));
      if (vec.SqrNorm() > -1)  // (i.e. not NaN)
      cout << "Region #" << region_id << "; motion vector: " << vec << '\n';
      else break;
      ++region_id;
    } while (true);
  }
}