%$ (c) Jan Modersitzki 2009/03/25, see FAIR.2 and FAIRcopyright.m.
%$ note: '%$' is used so hide comments in the book                              
%$ ==  CALL MLPIR with different distance measures  ========================

% set-up some data
setupPETCTdata; %close all;
level = 6; omega = MLdata{level}.omega; m = MLdata{level}.m; 

DM = {'SSD','NCC','MIcc','NGFdot'};

FAIRprint('reset','folder',fullfile(fairPath,'temp','PETCT'),...
  'prefix',[mfilename,'-'],...
  'pause',1,'obj','gcf','format','jpg','draft','off');
FAIRprint('disp');

str = @(w) sprintf('T(y(%s^o))',num2str(w*180/pi));
variable = @(k)['DM',DM{k}];

theta = [0,10];


for q=1:length(theta),
    fprintf('============== %s ====================\n\n',...
      sprintf('theta=%s',num2str(theta(q))));
  
    % initialize interpolation [spline] and troansformation [rotation]
    inter('reset','inter','splineInter2D','regularizer','moments','theta',theta(q));
    [T,R] = inter('coefficients',MLdata{level}.T,MLdata{level}.R,omega);
    X  = getCenteredGrid(omega,m);
    Rc = inter(R,omega,X);
    center = (omega(2:2:end)-omega(1:2:end))'/2;
    trafo('reset','trafo','rotation2D','c',center);
    trafo('w0');

    % run loo over the following distance measures
    edge = 10;
    filename = sprintf('%s-%s-theta=%s.mat',mfilename,'rot',num2str(theta(q)));
    if ~exist(filename,'file'),
      w = pi/2*linspace(-1,1,101)';
      save(filename,'w','DM');
    else
      clear DM*
      load(filename)
    end;


    for k=1:length(DM),
      var = whos('-file',filename);
      j = find(strcmp({var(:).name},variable(k))==1);

      if isempty(j),
        fprintf('============== %s ====================\n\n',variable(k))
        disp([variable(k),'=dm;']);
        dm = zeros(size(w));
        for j=1:length(w),
          Y = trafo(w(j),X(:));
          Tc = inter(T,omega,Y);
          dm(j) = feval(DM{k},Tc,Rc,omega,m,'edge',edge);
          if j == 1,
            figure(k); clf;
            subplot(1,3,1); viewImage(Rc,omega,m); title('reference');
            subplot(1,3,2); ph = viewImage(Tc,omega,m); th = title(str(w(j)));
            subplot(1,3,3); rh = plot(w(j),dm(j),'r.','markersize',20);
            axis([w(1),w(end),-inf,inf]); hold on;      title(DM{k})
          else
            set(ph,'cdata',reshape(Tc,m)'); set(th,'string',str(w(j)));
            subplot(1,3,3); set(rh,'visible','off');
            plot(w(1:j),dm(1:j),'k-','linewidth',2);
            rh = plot(w(j),dm(j),'r.','markersize',20);    pause(1/100)
          end;
          fprintf('.'); if ~rem(j,50) || j==length(w), fprintf('\n'); end;
        end;
        eval([variable(k),'=dm;']);
        save(filename,'-append',variable(k));
        fprintf('============== %s done ===============\n\n',variable(k))
      end;
    end;

    Name = @(k) sprintf('%s-theta=%s',DM{k},num2str(theta(q)));
    for k=1:4,
      eval(['dm=',variable(k),';']);
      dm = dm(1:length(w));
      [wOpt,j] = min(dm);
      FAIRfigure(k);
      ph = plot(w,dm,'-',w(j),dm(j),'*');
      set(ph,'linewidth',2,'color','k','markersize',20);
      set(gca,'fontsize',30);
      a = max(dm)-0.2*(max(dm)-min(dm));
      th = text(w(j),a,['$w^*=',num2str(w(j)),'$']);
      set(th,'fontsize',30,'interpreter','latex','horizontalalignment','center');
      FAIRprint(Name(k));
    end;
    fprintf('=====================================================================\n');
    
end;
return
% =========================================================================
