% cements 2 or more frames together, given the CUMULATIVE cross chirp between
%
% parameters are vectors, cumulatively expressed already
%
% Examples of use: first use P=pintegrate(p,3); % make cumulative, select lref
%                : C=cement(P,l1,l2); % calls loadframe.m to get ea. frame
%                : C=cement(P,l1,l2,skip_increment);
%                : C=cement(P,1,100,10); % uses every 10th frame
%
% See also: cement_cumulative
%          (takes list of variables rather than calling loadframe)
%
% calls ``loadframe.m'' and ``pchirp2nocrop.m''

function empty=myname(P,l1,l2,skip_increment)

more off % make sure doesn't stop if screen full of text...

if nargin==3
  disp('cement: skip_increment defaulting to 1; we''ll use every frame in interval')
  skip_increment=1;
end%if

if nargin==4
  dis('cement: running from frame %d, to frame %d, in increments of %d',...
      l1,l2,skip_increment)
end%if

if max(max(size(l1))) > 1
  error('cement: starting frame number (l1) must be a scalar')
end%if
if max(max(size(l2))) > 1
  error('cement: ending frame number (l2) must be a scalar')
end%if
if isstr(l1)
  error('cement: starting frame (l1) must be a NUMBER, not a string')
end%if
if isstr(l2)
  error('cement: ending frame (l2) must be a NUMBER, not a string')
end%if

[M_P, N_P]=size(P);
if N_P ~= 8
  error('cement: input parameters must be in matrix with rows of width 8')
end%if
if isstr(P)
  error('cement: input parameters must be numbers not a string')
end%if
if ~isreal(P)
  error('cement: input parameters must be real')
end%if

[Ml,Nl,Mh,Nh]=boundingboxextrema([P(l1+1:skip_increment:l2+1,:)]);
                                                % v03 to v15 uses P(4:16,:)
% note that P includes reference parameter which may or may not be the identity

for l=l1:skip_increment:l2
  dis('doing frame %d <=  ******* %d ******* <= %d of %d frames total',l1,l,l2,l2-l1+1)
  eval(sprintf('v%02d=loadframe(%02d);',l,l))
  if l==l1 % find size from first frame, hope that all frames same size
    [M, N]=size(eval(sprintf('v%02d',l1)));
    % base on size of one of frames which should be equal so base on frame l1
    canvassize=[M*(Mh-Ml),  N*(Nh-Nl)];
    canvassize=2*ceil(canvassize/2);  % getting even with mpeg
    disp(sprintf('cement_cumulative: about to nalloc %d by %d canvas (approx %7.3g Meg)',canvassize, 8*prod(canvassize)/(1024.^2)))
    empty=NaN*ones(ceil(canvassize));
  else%if % each time, check that other frames are same size as first
    [Mnow, Nnow]=size(eval(sprintf('v%02d',l)));
    if Mnow~=M
      error(sprintf('cement: frame %d different height than first frame',l))
    end%if
    if Nnow~=N
      error(sprintf('cement: frame %d different width than first frame',l))
    end%if
  end%if l==0
  % now insert frame l:
  eval(sprintf('[v%02dd,Ml%02d,Nl%02d]=pchirp2nocrop(P(%02d+1,:),v%02d);',...
               l,l,l,l,l))  % note P(4,:) pertains to v03, etc.
  eval(sprintf('clear v%02d',l))
  eval(sprintf('[M%02d,N%02d]=size(v%02dd);',l,l,l)) % size after dechirp
  eval(sprintf('Mh%02d=Ml%02d+(M%02d-1)/M%02d; Nh%02d=Nl%02d+(N%02d-1)/N%02d;',l,l,l,l,l,l,l,l))
  eval(sprintf('Mpixell%02d=round(M*Ml%02d+1   - M*Ml);',l,l))
  eval(sprintf('Mpixelh%02d=round(M*Ml%02d+M%02d - M*Ml);',l,l,l))
  eval(sprintf('Npixell%02d=round(N*Nl%02d+1   - N*Nl);',l,l))
  eval(sprintf('Npixelh%02d=round(N*Nl%02d+N%02d - N*Nl);',l,l,l))




  % UNCOMMENT ONE OF THESE (max, mean, new, newbroken):

  %eval(sprintf('empty(Mpixell%02d:Mpixelh%02d, Npixell%02d:Npixelh%02d)=maxnan(empty(Mpixell%02d:Mpixelh%02d, Npixell%02d:Npixelh%02d),v%02dd);',l,l,l,l,l,l,l,l,l))

  %eval(sprintf('empty(Mpixell%02d:Mpixelh%02d, Npixell%02d:Npixelh%02d)=meannan(empty(Mpixell%02d:Mpixelh%02d, Npixell%02d:Npixelh%02d),v%02dd);',l,l,l,l,l,l,l,l,l))

  % meancertain seems to give parts of the image looking negated for image sequence that doesn't have agc or such changes...
  %eval(sprintf('empty(Mpixell%02d:Mpixelh%02d, Npixell%02d:Npixelh%02d)=meancertainnan(empty(Mpixell%02d:Mpixelh%02d, Npixell%02d:Npixelh%02d),v%02dd);',l,l,l,l,l,l,l,l,l))

  %eval(sprintf('empty(Mpixell%02d:Mpixelh%02d, Npixell%02d:Npixelh%02d)=newnan(empty(Mpixell%02d:Mpixelh%02d, Npixell%02d:Npixelh%02d),v%02dd);',l,l,l,l,l,l,l,l,l))
  
  eval(sprintf('empty(Mpixell%02d:Mpixelh%02d, Npixell%02d:Npixelh%02d)=feathernewnan(empty(Mpixell%02d:Mpixelh%02d, Npixell%02d:Npixelh%02d),v%02dd,21);',l,l,l,l,l,l,l,l,l)) % last arg of feathernewnan is feather radius
  
  %eval(sprintf('empty(Mpixell%02d:Mpixelh%02d, Npixell%02d:Npixelh%02d)=v%02dd;',l,l,l,l,l))





  eval(sprintf('clear v%02dd',l))
end%for


%% what used to be there before eval(sprintf(' 
%% % insert v02
%% if LENGTH_P >=02+1
%%   [v02d,Ml02,Nl02]=pchirp2nocrop(P(02+1,:),v02);
%%   clear v02
%%   [M02,N02]=size(v02d);
%%   Mh02=Ml02+(M02-1)/M02; Nh02=Nl02+(N02-1)/N02;
%%   Mpixell02=round(M*Ml02+1   - M*Ml);
%%   Mpixelh02=round(M*Ml02+M02 - M*Ml);
%%   Npixell02=round(N*Nl02+1   - N*Nl);
%%   Npixelh02=round(N*Nl02+N02 - N*Nl);
%%   empty(Mpixell02:Mpixelh02, Npixell02:Npixelh02)=maxnan(...
%%   empty(Mpixell02:Mpixelh02, Npixell02:Npixelh02),v02d);
%%   clear v02d
%% end%if

%% % insert v03
%% ...