function [fragments, frames] = dicom_encode_rle(pixel_cells, bits_allocated, dims)
%DICOM_ENCODE_RLE   Run-length encode pixel cells.
%   [FRAGMENTS, LIST] = DICOM_ENCODE_RLE(PIXCELLS, BITS, DIMS) compresses
%   the pixel cells PIXCELLS using run-length encoding.  BITS is the
%   number of bits allocated for each pixel cell, and DIMS is the
%   original dimensions of the image stored in PIXCELLS, which is now a
%   vector.  FRAGMENTS is a cell array containing the encoded frames (as
%   UINT8 data) from the compressor.  LIST is a vector of indices to the
%   first fragment of each compressed frame of a multiframe image.
%
%   See also DICOM_ENCODE_JPEG_LOSSY.

%   Copyright 1993-2011 The MathWorks, Inc.


% Run-length encoded pixel cells have the following form:
%
%   RLE Header (16 UINT32 components)
%     - The first component is number of RLE segments (0 - 15)
%     - The remaining 15 components are offsets to the beginning of each
%       segment (or 0 if there are no more segments).
%   RLE Segment 1
%   RLE Segment 2
%   etc.
%
% The segments are generated by stripping off successive bytes from the
% composite pixel code, starting with the most significant byte [1].
% The stripped bytes are passed through the RLE coder to create the
% compressed segment.
%
% Individual image rows must be sent through the coder separately even if
% they are part of the same segment.  This prevents a run from spanning
% multiple rows.  The compressed rows are concatenated within each
% segment.
%
% If a composite pixel cell does not end on byte boundaries, 0 bits are
% appended to the least significant byte.  This will happen if the number
% of bits allocated for pixels and overlays (0028,0100) is not divisible
% by 8.
%
% Sequences of bytes can be compressed into an arbitrary number of
% fragments.  This allows for encoding the image into a set of "strips".
% The only restriction is that each fragment must contain the same number
% of pixel cell elements in each segment (e.g., the same number of red,
% green, and blue samples, etc.).
%
% [1] - A composite pixel code is one pixel cell in the sequence of pixel
%       cells.  The unstated conclusion of using composite pixel codes is
%       that planar configuration (0028,0006) must be 0 if more than one
%       sample is present and RLE coding is used.
%
% See Annex G in part 5 of the DICOM spec for more details.

% DICOM_ENCODE_PIXEL_CELLS explicitly creates pixel cells that are byte
% aligned.
if (rem(bits_allocated, 8) ~= 0)
    error(message('images:dicom_encode_rle:badBitsAllocated'))
end

% Reshape data to make segment creation easy.
if (length(dims) >= 3)
    num_segments = (bits_allocated / 8) * dims(3);
else
    num_segments = (bits_allocated / 8);
end

if (num_segments > 15)
    error(message('images:dicom_encode_rle:tooManySegments'))
end

if (numel(dims) == 4)
    numFrames = dims(4);
else
    numFrames = 1;
end

swap = determine_swap;

pixel_cells = dicom_typecast(pixel_cells, 'uint8', swap);
pixel_cells = reshape(pixel_cells, [num_segments, ...
                                    (numel(pixel_cells) / num_segments)]);

offset = 0;

fragments = cell(1, numFrames);
for f = 1:numFrames

    encoded_data = [];

    % Encode each segment.
    segment_lengths = zeros(1, num_segments);
    for p = 1:num_segments

        segment = [];

        % Encode each row.
        for q = 1:dims(1)

            start_idx = offset + (q - 1) * dims(2) + 1;
            end_idx = offset + q * dims(2);

            encoded_row = dicom_encode_rle_segment(pixel_cells(p, start_idx:end_idx));

            segment = [segment; encoded_row]; %#ok<AGROW>

        end

        % Store segment.
        encoded_data = [encoded_data; segment]; %#ok<AGROW>
        segment_lengths(p) = numel(segment);

    end

    offset = end_idx;

    % Create RLE header.
    header = repmat(uint32(0), [16 1]);
    header(1) = uint32(num_segments);
    header(2) = uint32(64);

    for p = 2:(num_segments)

        header(p + 1) = uint32(64 + sum(segment_lengths(1:(p - 1))));

    end

    encoded_data = [dicom_typecast(header, 'uint8', swap); encoded_data]; %#ok<AGROW>

    % Create output from fragments.
    fragments{f} = encoded_data;

end

frames = 1:numFrames;



% This is essentially how dicom_rle_encode_segment() works.
% function coded_bytes = rle_coder(X)
% 
% X = double(X(:));
% 
% current_pos = 1;
% max_pos = numel(X);
% 
% % Find runs.  Where the diff is 0, there is a repeated value.  The first run
% % is not included.
% X_diff = diff(X);
% X_reps = find(X_diff == 0);
% rep_pix = [X_reps(1); X_reps(find(diff(X_reps) > 1) + 1); (max_pos + 1)];
% 
% all_runs = find(X_diff ~= 0) + 1;
% 
% rep_number = 1;
% next_rep = rep_pix(rep_number);
% 
% coded_values = [];
% more_data = 1;
% 
% while (more_data)
%     
%     if (current_pos ~= next_rep)
%         
%         % Handle literal values between runs.
%         
%         % Only 128 literal values can be encoded in one literal run.
%         while ((next_rep - current_pos) > 128)
%             
%             coded_values = [coded_values
%                             127
%                             X(current_pos:(current_pos + 127))];
%             current_pos = current_pos + 128;
%             
%         end
%         
%         run_length = next_rep - current_pos;
%         coded_values = [coded_values
%                         (run_length - 1)
%                         X(current_pos:(next_rep - 1))];
%         
%         % All literal values until the next repetitive run are encoded.
%         current_pos = next_rep;
%         
%     else
%         
%         % Handle repetitive runs.
%         
%         % Find next run of different values.
%         idx = find(all_runs == current_pos);
%         
%         if (isempty(idx))
%             
%             % This run appeared at the very beginning of the data.
%             if (isempty(all_runs))
%                 
%                 % There is only one run.
%                 next_run = max_pos + 1;
%                 
%             else
%                 
%                 next_run = all_runs(1);
%                 
%             end
%             
%         elseif (idx == numel(all_runs))
%             
%             % This is the last run.
%             next_run = max_pos + 1;
%             
%         else
%             
%             next_run = all_runs(idx + 1);
%             
%         end
%         
%         while ((next_run - current_pos) > 128)
%             
%             coded_values = [coded_values
%                             -127
%                             X(current_pos)];
%             current_pos = current_pos + 128;
%             
%         end
%         
%         run_length = next_run - current_pos;
%         coded_values = [coded_values
%                         (-run_length + 1)
%                         X(current_pos)];
%         
%         current_pos = current_pos + run_length;
%         
%         rep_number = rep_number + 1;
%         next_rep = rep_pix(rep_number);
%         
%     end
%     
%     if (current_pos > max_pos)
%         
%         more_data = 0;
%         
%     end
%     
% end
% 
% % Convert double data to bytestream.
% %
% % The next statement maps the signed values in the range [-128, 255] to 
% % [0, 255], the valid range for UINT8.  Values in the range [0, 255] are
% % unchanged; values in [-128, -1] are mapped to [128, 255].
% %
% % Converting negative values is okay because the bit pattern for a
% % converted value as a UINT8 is the same as the negative value's bit
% % pattern as an INT8.  The RLE decoder algorithm will read the next byte
% % as a signed or unsigned value based on its own logic; it suffices for
% % this function to write indistinguishable bytes (signed or unsigned).
% coded_bytes = uint8(rem((coded_values + 256), 256));



function swap = determine_swap
% Determine if endianness of data needs swapped.

[~, ~, mach] = fopen(1);

% RLE compressed data must be written to file as IEEE-LE.
if (isempty(strfind(lower(mach), 'ieee-le')))
    swap = 1;
else
    swap = 0;
end