function [] = iterativeRegistering_memory(vidname, NCparamStruct, regname)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Consecutive registrations with increasing grid number (2x2 -> 3x3 -> 6x6)
% using memory with peaks identified from the watershed algorithm
% Inputs:
%     - vidname: The name of the video. Assumes that the videos are h5
%     files
%     - NCparamStruct: A struct of the desired parameters for NormCorre
%     - regname: The desired name for the registered video
% Output: (none)
%     - 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%% load in the data
M2 = h5read(vidname, '/video');
data_type = class(M2);
[d1 d2 d3] = size(M2);
%%% 'nMAD', NCparamStruct.nMAD(1),  removed
options_nr = NoRMCorreSetParms_memory('d1',d1(1),'d2',d2(1),'bin_width', 50, ...
                    'grid_size',NCparamStruct.gridsize(1, :),'mot_uf',NCparamStruct.motuf(1),'correct_bidir',false, ...
                    'overlap_pre',NCparamStruct.overlap(1),'overlap_post', NCparamStruct.overlap(1),'max_shift',....
                    NCparamStruct.maxshift(1), 'upd_template', false, 'shifts_method', NCparamStruct.interp,...
                    'us_fac', NCparamStruct.usfac(1), 'Quant', NCparamStruct.Quant(1),'nPeaks', NCparamStruct.nPeaks(1),...
                    'init_batch', 50);

nframe = min([ceil(d3*0.05) 50]);                

svid = normcorre_batch_memory(M2(:, :, 1:nframe),options_nr);                
templateim = mean(svid, 3);

%% Calculate initial metrics
[cM2f(:, 1),mM2f,vM2f(1)] = motion_metrics(M2, 0);
% Check to see if the h5 file exists. If it doesn't create it. 
clear video
iteration = 0; 

if ~isfile(regname)
    h5create(regname, '/iteration', size(iteration));
    h5write(regname, '/iteration', iteration);
    h5create(regname, '/originalname', size(string(vidname)), 'Datatype', 'string');
    h5write(regname, '/originalname', string(vidname));
end

%% Go through the different stages of registration
maxmove = [0 0 0 0]; 
assocdat = strcat(regname(1:end-3), '_associated.mat');

for fn = 1:size(NCparamStruct.gridsize, 1)
    if fn == 1
    else
        
        templateim = mean(M2(:, :, 1:nframe), 3);
    end
    disp(regname)
    completed = h5read(regname, '/iteration');
    if completed < fn
        % Check to see if the registration completed successfully because
        % sometimes it can fail in the upsampling factor is too high
        count = 0;
        display(num2str(count))
        err_count = 0;
        M2 = M2(maxmove(fn, 3)+1:end-maxmove(fn, 4), maxmove(fn, 1)+1:end-maxmove(fn, 2), :);
        templateim = templateim(maxmove(fn, 3)+1:end-maxmove(fn, 4), maxmove(fn, 1)+1:end-maxmove(fn, 2), :);
        [d1(fn+1) d2(fn+1) d3(fn+1)] = size(M2);
        
        options_nr = NoRMCorreSetParms_memory('d1',d1(fn+1),'d2',d2(fn+1),'bin_width', 50, ...
            'grid_size',NCparamStruct.gridsize(fn, :),'mot_uf',NCparamStruct.motuf(fn),'correct_bidir',false, ...
            'overlap_pre',NCparamStruct.overlap(fn),'overlap_post', NCparamStruct.overlap(fn),'max_shift',....
            NCparamStruct.maxshift(fn), 'upd_template', false, 'shifts_method', NCparamStruct.interp,...
            'us_fac', NCparamStruct.usfac(fn), 'Quant', NCparamStruct.Quant(fn), 'nPeaks', NCparamStruct.nPeaks(fn),...
            'init_batch', 10, 'max_dev', NCparamStruct.max_dev(fn));
        
        tic; [M_n,shifts2{fn+1},template2, opt, peaks{fn+1}, transKey{fn+1}] = normcorre_batch_memory(M2,options_nr, templateim); toc
        save(assocdat, 'peaks', 'transKey', 'maxmove', 'd1', 'd2', 'shifts2', 'NCparamStruct');
        itname = strcat('/video', num2str(fn));

        peakshift{fn+1} = DijkstraPeakSearch(peaks{fn+1}, shifts2{fn+1}, transKey{fn+1}, NCparamStruct.distthresh(fn), NCparamStruct.fillthresh(fn), NCparamStruct.buffer(fn), NCparamStruct.maxshift(fn), 0);

        save(assocdat, 'peakshift', 'peaks', 'transKey', 'maxmove', 'd1', 'd2', 'shifts2', 'NCparamStruct');
        cshift = convertShifts(peakshift{fn+1});
        M2 = apply_shifts(M2,peakshift{fn+1},options_nr,0, 0);
        maxmove(fn+1, :) = ceil([abs(max([0 max(cshift.x(:))])) abs(min([0 min(cshift.x(:))])) abs(max([0 max(cshift.y(:))])) abs(min([0 min(cshift.y(:))]))]);

        infcheck = h5info(regname);
        if ~contains(strcat([infcheck.Datasets.Name]), itname(2:end))
            h5create(regname, itname, size(M2), 'Datatype', data_type);
        end
        h5write(regname, itname, M2);
        iteration = fn;
        h5write(regname, '/iteration', iteration);
        
    end
end

%% Save the associated data into a .mat file for easier use

h5create(regname, '/datname', size(string(assocdat)), 'Datatype', 'string');
h5write(regname, '/datname', string(assocdat));
save(assocdat, 'peakshift', 'peaks', 'transKey', 'maxmove', 'd1', 'd2', 'shifts2', 'NCparamStruct');

end