//--- lecture des N premiers caractères d'un fichier texte (ici goriot.txt) 
//--- au format ascii, ie chaque caractère est codé 
//--- par un entier entre 0 et 255
//--- N peut-être supérieur à la taille du fichier: dans ce cas, le fichier
//--- sera lu en entier, et le nombre de caractères lus sera length(x)

N=2^12
f=mopen("goriot.txt");
x=ascii(mgetstr(N,f));
mclose(f);

//--- affichage du texte à l'écran
//--- on peut aussi utiliser char(x) tout seul, 
//--- mais les sauts de ligne sont alors hasardeux...

printf("%s\n",char(x))



//----------------------------------------------------------------------------
//                                  bin_add
//----------------------------------------------------------------------------
// calcule c + x*2^(-l), et renvoie le résultat dans c 
// c est un tableau binaire représentant le nombre 0,c1c2c3...cN en base deux
// x et l sont des entiers
// le résultat est donné sous forme binaire 0,c1c2c3...cN
// si le résultat est supérieur ou égal à 1, la partie entière est tronquée
//----------------------------------------------------------------------------

function c=bin_add(c,x,l)
  if l>length(c) then c=[c zeros(1,l-length(c))]; end
  i=l
  r=0
  while (i>0) & (x~=0 | r~=0)
    c(i)=c(i)+modulo(x,2)+r
    if (c(i)>1) then r=1; c(i)=c(i)-2; else r=0; end 
    x=int(x/2)
    i=i-1
  end
endfunction


//----------------------------------------------------------------------------
//                                  bin_leq
//----------------------------------------------------------------------------
// renvoie le résultat (booléen) du test c <= d
// c et d sont des tableaux binaires représentant des nombres 
// de l'intervalle [0,1[ sous la forme 0,b1b2b3...bN en base deux
//----------------------------------------------------------------------------

function t=bin_leq(c,d)
  l=min(length(c),length(d))
  f=find(d(1:l)-c(1:l)~=0)
  if f==[] then
    if length(c)==length(d) then t=%T  
    else t=(length(c)<length(d)) | (find(c(l+1:$)~=0)==[])
    end
  else t=(c(f(1))<d(f(1)))
  end
endfunction



//----------------------------------------------------------------------------
//                                  bruite
//----------------------------------------------------------------------------
// simulation d'un canal bruité binaire symétrique
// x est un vecteur binaire (chaque élément vaut 0 ou 1)
// y(i) vaut 1-x(i) avec probabilité p et x(i) avec probabilité 1-p
//----------------------------------------------------------------------------

function y=bruite(x,p)
  r=(rand(x)<p)
  y=x
  y(r)=1-x(r)
endfunction


// exemple d'utilisation
x=ascii("000000000000000000000000")-48;
char(bruite(x,0.1)+48)


//----------------------------------------------------------------------------
//                             code_hamming74
//----------------------------------------------------------------------------
// encodage par le code de Hamming 7,4
// en entrée s est un vecteur binaire de taille 4n 
// en sortie x est un vecteur binaire de taille 7n 
// les bits de parités sont stockés juste après chaque bloc de 4 bits
//----------------------------------------------------------------------------

function x=code_hamming74(s)
  n=length(s)/4
  I=0:n-1
  x(7*I+1)=s(4*I+1)
  x(7*I+2)=s(4*I+2)
  x(7*I+3)=s(4*I+3)
  x(7*I+4)=s(4*I+4)
  x(7*I+5)=modulo(s(4*I+2)+s(4*I+3)+s(4*I+4),2)
  x(7*I+6)=modulo(s(4*I+1)+s(4*I+3)+s(4*I+4),2)
  x(7*I+7)=modulo(s(4*I+1)+s(4*I+2)+s(4*I+4),2)
endfunction


//----------------------------------------------------------------------------
//                             decode_hamming74
//----------------------------------------------------------------------------
// decodage par le code de Hamming 7,4
// en entrée x est un vecteur binaire de taille 7n 
// en sortie y est un vecteur binaire de taille 4n 
//----------------------------------------------------------------------------

function x=decode_hamming74(y)
  n=length(y)/7
  H=[0 0 0 1 1 1 1;0 1 1 0 0 1 1;1 0 1 0 1 0 1]
  z=matrix(y,7,-1)
  E=[4 2 1]*modulo(H*z,2)
  J=find(E==1); z(1,J)=1-z(1,J)
  J=find(E==2); z(2,J)=1-z(2,J)
  J=find(E==3); z(3,J)=1-z(3,J)
  J=find(E==4); z(4,J)=1-z(4,J)
  x=matrix(z(1:4,:),1,-1)
endfunction


// exemple d'utilisation
s=ascii("011010100001010100101010")-48;
x=code_hamming74(s);
y=bruite(x,0.05);
z=decode_hamming74(y);
char(48+s) // signal avant codage
char(48+z) //signal après décodage
sum(modulo(x+y,2)) //nombre d'erreurs introduites par le canal
sum(modulo(s+z,2)) //nombre d'erreurs non corrigées