close all;
clear all;
clc;
%*************************************************************************
%*  Peter Miščík
%*  Program Genbit - kompresia sekvencie
%*  vstup - DNA sekvencia (FASTA, EMBL, GENBANK)
%*  výstup - *.genbit
%*************************************************************************
reply1 = input('Zvoľte typ formátu DNA sekvencie (FASTA,EMBL,GENBANK):', 's');
switch reply1
    case 'FASTA'
        sequence = input('\nZadajte názov sekvencie s priponou fasta: ','s');
        [header seq_in] = fastaread (sequence);   % Načítanie DNA sekvencie fasta formatu
        N = length(seq_in); % Dĺžka sekvencie
        alpha = 'ACGT'; % abeceda sekvencie
    case 'EMBL'
        sequence = input('\nZadajte názov sekvencie s priponou: ','s');
        seq_in=emblread (sequence, 'SequenceOnly',true); 
        N = length(seq_in); % načítanie z EMBL
        alpha = 'acgt';
    case 'GENBANK'
        sequence = input('\nZadajte názov sekvencie s priponou: ','s');
        seq=genbankread(sequence); % načítanie z genbank
        seq_in=seq.Sequence;
        N = length(seq_in);
        alpha = 'acgt';
end

reply2 = input('\nChcete zadať úsek sekvencie?(A/N): ', 's');
if reply2=='A'
    str1=sprintf('\nZvoľte interval sekvencie od 1 do %d. Od: ',N);
    a=input(str1);
    str2=sprintf('\nDo: ');
    b=input(str2);
    seq_in=seq_in(a:b);
    N=length(seq_in);
end

str=sprintf('\nPrebieha kompresia...');
disp(str)

code = ({'00','01','10','11'}); % Kódy daným bázam

encseq = '';

tic;

for i=1:N
    idx = find(seq_in(i)==alpha);  % for cyklus nahradenie báz 2 bitmy
    encseq = strcat(encseq,cell2mat(code(idx))); 
end

Z = length(encseq);
ninebin='';
osembin=[];

j = 1;
k = 1;
n = 2;

for i = 1:8:(length(encseq)-mod(Z,8))
    osembin(j) = bin2dec(encseq(i:i+7)); % rozdelenie encseq na 8 bit. čísla
    j=j+1;
end

while n <= length(osembin)          % while cyklus, vyhľadávanie 
    if (osembin(n)==osembin(n-1))   % po sebe zhodných 8 bitových úsekov
       ninebin{k}=strcat('1',dec2bin(osembin(n-1),8));
       n = n+2;
    else
       ninebin{k}=strcat('0',dec2bin(osembin(n-1),8));
       n = n+1;
       if (n==length(osembin)+1)
        ninebin{k+1}=strcat('0',dec2bin(osembin(n-1),8));
       end
    end
k = k+1;
end
if (mod(Z,8)>0)     % určenie zvyšku 
    z = mod(Z,8);
    ninebin{k+1} = (encseq((Z-z+1):Z));
end

ninebinmat=cell2mat(ninebin); % úprava výstupného reťazca
M=length(ninebinmat);
LengthBin=M;
Z=mod(M,8);
if Z>0
    D='00000000';
    ninebinmat=strcat(ninebinmat(1:end-Z),D(1:8-Z),ninebinmat(end-Z+1:end));
    M=length(ninebinmat);
    LengthBin=M;
    ninebin=(reshape(ninebinmat(1:end),8,M/8))';
    ninedec=bin2dec(ninebin);
    ninedec(end+1)=Z;
else
    ninebin=(reshape(ninebinmat(1:end),8,M/8))';
    ninedec=bin2dec(ninebin);
end

CompressionRatio=LengthBin/N

uloz=input('Zadajte názov súboru pre uloženie: ','s');
uloz=strcat(uloz,'.genbit');

fid = fopen(uloz,'wb'); % uloženie súboru
fwrite(fid,ninedec,'uint8');
fclose(fid);

str=sprintf('\nUložený súbor je: %s',uloz);
disp(str)