disp('Running MATLAB script file nmft.f') ;
%
%    Both makes plot of Normal mixture densities,
%      from Marron and Wand (1992),
%      plus number 16 from Janssen, et. al.
%      and also tests and develops the MATLAB proc nmf.m
%      which evaluates normal mixture densities
%
%    Reads in parameters from file nmpar.dat
% 
%    Makes plot of all on one 


itest = 5 ;     %  1 - Makes original Marron and Wand 15,  on 3 x 5 plot
                %  2 - Adds in JMVS # 16,                  on 4 x 4 plot
                %  3 - tests a single entry
                %  4 - tests another single entry
                %  5 - tests another single entry

if itest <= 2 ;    % then load in usual parameter matrices

  %  Load Normal mixtures parameters from nnpar.mat
  load \matlab\steve\nmstuff\nmpar.mat

  if itest == 1 ;    %  then work with original 15
    idiste = 15 ;
    plotr = 3 ;
    plotc = 5 ;
  elseif itest == 2 ;    %  then work with original 16
    idiste = 16 ;
    plotr = 4 ;
    plotc = 4 ;
  end ;

  xgrid = linspace(-3,3,401)' ;

  %  Loop through distributions
  for idist = 1:idiste ;
    eval(['parmat = nmpar' num2str(idist) ' ;']) ;
    eval(['titstr = nmtit' num2str(idist) ' ;']) ;
    disp(['    Working on ' titstr]) ;

    %  do a quick check of the wt vector
    if abs(sum(parmat(:,3)) - 1) > 10^(-14) ;
      disp('!!!    Warning: weights don''t sum to 1 !!!') ;
    end ;

    subplot(plotr,plotc,idist) ;
    f = nmf(xgrid,parmat) ;
    plot(xgrid,f) ;
      title(['#' num2str(idist) ' ' titstr]) ;

  end ;

else ;     %  then work with single entry versions 

  if itest == 3 ;
    vmu = [-2; 2] ;
    vsig2 = [.25; .25] ;
    vw = [.3; .7] ;

    xgrid = linspace(-4,4,801)' ;
    f = nmf(xgrid,vmu,vsig2,vw) ;
    disp(['  Check that this is 1:   ' num2str(sum(f * .01))]) ;
    disp(['  Size of vector is:  ' num2str(size(f))]) ;

  elseif itest == 4 ;
    vmu = [-2; 0; 2] ;
    vsig2 = [.25; 1; .25] ;
    vw = [.3; .4; .3] ;

    xgrid = linspace(-4,4,801)' ;
    f = nmf(xgrid,vmu,vsig2,vw) ;
    disp(['  Check that this is 1:   ' num2str(sum(f * .01))]) ;
    disp(['  Size of vector is:  ' num2str(size(f))]) ;

  elseif itest == 5 ;
    vmu = [-2] ;
    vsig2 = [.25] ;
    vw = [1] ;

    xgrid = linspace(-4,4,801)' ;
    f = nmf(xgrid,vmu,vsig2,vw) ;
    disp(['  Check that this is 1:   ' num2str(sum(f * .01))]) ;
    disp(['  Size of vector is:  ' num2str(size(f))]) ;

  end ;

end ;
