order = 5;

poly1 = [1 0 0 0 0 1]; % X6+X5+X2+X+1

poly2 = [1 1 0 1 1 1]; %X6 + X + 1

poly3 = [1 1 1 1 0 1 ]; %X6+X5+X4+X+1

poly4 = [1 0 0 1 0 1];

%Generating Sequences

polyF = fliplr(poly1); % Flipped version of polynomial

register1 = [1 zeros(1,order-1)];

for i=1:(2^order -1)

    Seq1(i) = register1(1);

    register1 = [register1(2:order) polyF(order+1)* rem(sum(polyF(1:order).*register1(1:order)),2)];

end

register2 = [1 zeros(1,order-1)];

polyF = fliplr(poly2);

for i=1:(2^order -1)

    Seq2(i) = register2(1);

    register2 = [register2(2:order) polyF(order+1)* rem(sum(polyF(1:order).*register2(1:order)),2)];

end

register3 = [1 zeros(1,order-1)];

polyF = fliplr(poly3);

for i=1:(2^order -1)

    Seq3(i) = register3(1);

    register3 = [register3(2:order) polyF(order+1)* rem(sum(polyF(1:order).*register3(1:order)),2)];

end

register4 = [1 zeros(1,order-1)];

polyF = fliplr(poly4);

for i=1:(2^order -1)

    Seq4(i) = register4(1);

    register4 = [register4(2:order) polyF(order+1)* rem(sum(polyF(1:order).*register4(1:order)),2)];

end

%mapping sequences to -1 and 1

x1 = transpose(2*Seq1 -1);

x2 = transpose(2*Seq2 -1);

x3 = transpose(2*Seq3 -1);

x4 = transpose(2*Seq4-2);

GPN =[x1 x2 x3 x4]

N = 256000; % Number of bits
N_Users = 4;
ip_bits = randi([0 1],N,1) ; % Bits Generation
ip_bits = reshape(ip_bits,[(N/N_Users),4]);
BPSK = (ip_bits.*2)-1; %Reshaping to have 4 separate in BPSK
Tx_Data = kron(BPSK(:,1),GPN(:,1)) + kron(BPSK(:,2),GPN(:,2)) + kron(BPSK(:,3),GPN(:,3)) + kron(BPSK(:,4),GPN(:,4));

%% Adding AWGN Noise %%
Eb_BPSK = 1;
BER_BPSK = zeros(4,9);
RAND_VAR  = (randn(size(Tx_Data,1),size(Tx_Data,2)));
BER1= [ ] ;
BER2 = [ ] ;
BER3 = [ ] ;
BER4 = [ ] ;
BER_Th = [ ];
for SNR= 1:9
    No_BPSK = size(GPN,1) / (10^(SNR/10));
    Variance_BPSK = sqrt(No_BPSK/2);
    Noise_BPSK = Variance_BPSK * (randn(size(Tx_Data,1),size(Tx_Data,2)));
    Tx_Data_Noisy = Tx_Data + Noise_BPSK;
    Rx_Data_RS = reshape(Tx_Data_Noisy,31,(N/N_Users));

    Rx_Data_RS = Rx_Data_RS' * GPN ;

    Demapping = sign(Rx_Data_RS);
    BER1 = [BER1 sum((BPSK(:,1)) ~= (Demapping(:,1)))/(N/N_Users)];
    BER2 = [BER2 sum((BPSK(:,2)) ~= (Demapping(:,2)))/(N/N_Users)];
    BER3 = [BER3 sum((BPSK(:,3)) ~= (Demapping(:,3)))/(N/N_Users)];
    BER4 = [BER4 sum((BPSK(:,4)) ~= (Demapping(:,4)))/(N/N_Users)];
    BER_Th = [BER_Th 0.5*erfc(sqrt((10^(SNR/10))))];
end
%% P1otting %%
%figure (1)
semilogy(1:9,BER1,'-rx',1:9,BER2,'-bx',1:9,BER3,'-kx',1:9,BER4,'-gx','LineWidth',1.5);
title ('Users BER');
legend ('User 1','User 2','User 3', 'User 4', 'Classical BPSK');
xlabel ('SNR (dB)');
ylabel ('BER');
grid on;