# OK to post

# PART 1: WtEnumerator(q,M,t) PROCEDURE:
# -----------------------

WtEnumerator:=proc(q,M,t) local v,n,C:
    n:=nops(M[1]);
    C:=LtoC(q,M);

    return add([seq(t^HD(v, [0$n]), v in C)]);
end;

# PART 2: RELATIONSHIP BETWEEN WtEnumerator(M) and WtEnumerator(H):
# -----------------------

# could not find a pattern (maybe my WtEnumerator function was wrong?)

############################################################## OLD CODE ##############################################################

# LtoC(q,M): inputs a list of basis vectors for our linear code over GF(q)
# outputs all the codewords (the actual subset of GF(q)^n with q^k elements)
LtoC:=proc(q,M) local n,k,C,i,M1,c:
    option remember;
    
    k:=nops(M);
    n:=nops(M[1]);

    if k=1 then
        return {seq(i*M[1] mod q, i=0..q-1)};
    fi;

    M1:=M[1..k-1];
    C:=LtoC(q, M1);
    return {seq(seq(c + i*M[k] mod q, i=0..q-1), c in C)};
end;

# hamming distance
HD:=proc(u,v) local i, c:
	# TODO: type checking
	
	if nops(u) <> nops(v) then
		return FAIL;
	fi;

	c:=0;
	for i from 1 to nops(u) do
		if u[i] <> v[i] then
			c:=c+1;
		fi;
	od;

	return c;
end;

# MinW(q,M): The minimal weight of the linear code generated by M over GF(q)
MinW:=proc(q,M) local n,c,C:
    n:=nops(M[1]);
    C:=LtoC(q,M);

    return min(seq(HD(c, [0$n]), c in C minus {[0$n]}));
end;

Fqn:=proc(q,n) local S,a,v:
	option remember:
	if n=0 then
		return {[]};
	fi;

	S:=Fqn(q,n-1);
	return { seq(seq([op(v),a],a=0..q-1), v in S) };
end;

# GLC1: tries to add a new number to the current basis M out of the other vectors
# that still has min weight d
GLC1:=proc(q,M,d) local n,A,c,M1:
    n:=nops(M[1]);
    A:=Fqn(q,n) minus LtoC(q,M);
    
    for c in A do
        M1:=[op(M),c];
        if MinW(q,M1)=d then
            return M1;
        fi;
    od;

    return FAIL;
end;

# GLC(q,n,d) inputs q and n and d (for min. distance), greedily and randomly tries to
# get as large a linear code give by a basis as possible
GLC:=proc(q,n,d) local M,M1:
    M:=[[1$d, 0$(n-d)]];
    M1:=GLC1(q,M,d):

    while M1<>FAIL do
        M:=M1;
        M1:=GLC1(q,M,d);
    od;

    return M;
end;

# PCM(q,M) inputs a non-empty basis describing some linear code (assumed in standard form)
# over GF(q)^n outputs the parity check matrix
PCM:=proc(q,M) local k,n,i,j,H:
    option remember:
    k:=nops(M):
    n:=nops(M[1]):

    if [seq([op(1..k,M[i])],i=1..k)]<>[seq([0$(i-1),1,0$(k-i)],i=1..k)] then
        print(`Not standard form , please use SFde(q,M) first `):
        return FAIL;
    fi:

    for i from 1 to n-k do:
        for j from 1 to k do
            H[i,j]:=-M[j][i+k]  mod q:
        od:

        for j from k+1 to n do
            if j-k=i then
                H[i,j]:=1:
            else
                H[i,j]:=0:
            fi:
        od:
    od:

    return [seq([seq(H[i,j],j=1..n)],i=1..n-k)]:
end:

# gets a row from a matrix
GetRow:=proc(M,r) local v,i,n:
    v:=[];
    n:=nops(M[1]);

    for i from 1 to n do:
        v:=[op(v),M[r,i]];
    od;

    return v;
end;

# sets a row in a matrix
SetRow:=proc(M,r,v) local i,M1:
    M1:=copy(M);
    
    for i from 1 to nops(v) do:
        M1[r,i]:=v[i];
    od;

    return M1;
end;

# gets a column from a matrix
GetCol:=proc(M,c) local v,i,n:
    v:=[];
    n:=nops(M);
    
    for i from 1 to n do:
        v:=[op(v),M[i,c]];
    od;

    return v;
end;

# sets a column in a matrix
SetCol:=proc(M,c,v) local i,M1:
    M1:=copy(M);
    
    for i from 1 to nops(v) do:
        M1[i,c]:=v[i];
    od;

    return M1;
end;

# returns matrix M in standard form
SF:=proc(q,M) local k,n,i,j,S,rj,cj,ri:
    k:=nops(M);
    n:=nops(M[1]);
    S:=copy(M);

    for i from 1 to k do: # algorithm is iterated from 1 to k
        # if S_ii = 0, then we need to perform a swap:
        if S[i,i] = 0 then 
            for j from i+1 to k while S[j,i] = 0 do od; # look for available row

            if j<=k then # swap rows
                rj:=GetRow(S,j);
                S:=SetRow(S,j, GetRow(S,i));
                S:=SetRow(S,i,rj);

            else # look to swap columns

                for j from i+1 to n while S[i,j] = 0 do od; # look for available col

                # swap cols
                cj:=GetCol(S,j);
                S:=SetCol(S,j, GetCol(S,i));
                S:=SetCol(S,i,cj);
            fi;
        fi;

        # scale row to have leading entry 1
        ri:=GetRow(S,i);
        ri:=(ri*(ri[i]&^(-1) mod q)) mod q;
        S:=SetRow(S,i,ri);

        for j from 1 to k do:
            if j <> i then
                rj:=GetRow(S,j);
                rj:=(rj - (rj[i] * ri mod q)) mod q;
                S:=SetRow(S,j,rj);
            fi;
        od;
    od;

    return S;
end;