problem with BER for OFDM estimation for wimax please help

[WALA]

Hello, I am working on estimating wimax system using pilots, i applied the OFDM block and made the pilot insertion based on the block arrangement but the problem am facing is that the BER stays in its max. value it doesn't change with the SNR . i have tryed the comb type with interp, the insertion in both time and freq with using interp2, but i have the same problem.
Please help me pleas please

- - - Updated - - -

Hello, I am working on estimating wimax system using pilots, i applied the OFDM block and made the pilot insertion based on the block arrangement but the problem am facing is that the BER stays in its max. value it doesn't change with the SNR . i have tryed the comb type with interp, the insertion in both time and freq with using interp2, but i have the same problem.
Please help me pleas please

this is my code. i have another problem with the graph of mean square error please help me

N    = 128;   % fft size
Pil  = N/4; % 32
S    = N-Pil;   % number of data subcarriers  112
nIteration = 500;
L    = 30;
GI   = N/4; % 32
Ep   = 2;
SNR_V = [0:3:27];
ber = zeros(1,length(SNR_V));   % initializing bit error rate

M = 2; % modulation object

% Generating data

T_Data=randint(N,1,M); %  128  1

% modulating data

mod_data = pskmod(T_Data,M); 

% serial to parallel conversion
par_data = reshape(mod_data,8,16); %    8 16
Ip=1:4:N;
Is = setxor(1:N,Ip); 
% pilot insertion
 pilot_ins_data=[];
 for Ip=1:4:13
     pilots= Ep*par_data(:,Ip);
     xp=[pilots  par_data(:,Ip+1) par_data(:,Ip+2) par_data(:,Ip+3)];
     pilot_ins_data=[ pilot_ins_data  xp];
 
 end

 pilot_ins_data; % 8 16
 PilotData=reshape(pilot_ins_data,8*16,1); % 128 1
 LocPilot= [ PilotData(1:8,1); PilotData(33:40,1); PilotData(65:72,1); PilotData(97:104,1) ]; % Locations of pilots after reshape 32 1
 
 % fourier transform time doamain data and normalizing the data

no_of_error=[];
ratio=[];
for i = 1 : length(SNR_V)
    SNR = SNR_V(i)
    for k = 1 : nIteration
        
        IFFT_data = ifft(PilotData,N); % 128 1

        % addition cyclic prefix

        Cylic_add_Data = [IFFT_data(N- GI + 1 : N); IFFT_data]; % 96 1

        % generating random channel coefficients
        h(1:L,1)  =     random('norm',0,1,L,1) + ...
                        j * random('norm',0,1,L,1);  
        h  = h./sum(abs(h));    % normalization   30 1

        ChannEffect = filter(h,1,Cylic_add_Data);    % channel effect  160 1

        % adding awgn noise
        Signal_Noise = awgn(ChannEffect , SNR - db(std(ChannEffect))); % normalization to signal power


        Cyclic_pre_Rem =  Signal_Noise(GI+1:N+GI);  %cyclic prefix removal 
 
        FFT_recdata =fft(Cyclic_pre_Rem,N);    % freq domain transform  128 1 

       % LS Estimation & 1D interpolation
       
           F = exp(2*pi*sqrt(-1)/N .* meshgrid([0:N-1],[0:N-1])...
                           .* repmat([0:N-1]',[1,N])); % 128 128
      
       Tpilot= LocPilot;% transmitted pilots 32 1
       
       Rpilot=[ FFT_recdata(1:8,1);   FFT_recdata(33:40,1);   FFT_recdata(65:72,1);  FFT_recdata(97:104,1) ]; % received pilots 32 1
     
      % Hpls=inv(diag(LocPilot))*  Rpilot; % 32 1
      
       Ipp=[ 1:8 33:40 65:72 97:104 ]; 
       
       G = (Ep * length(Ipp))^-1 ...
                         * ctranspose(sqrt(Ep)*diag(Tpilot)*ctranspose(F(1:L,Ipp))); % 30 32
       Hpls = G*Rpilot; % 30 1
              
       Hh=interp(Hpls,4); % 120 1
 
       Demod_Data = pskdemod(  FFT_recdata./(fft(Hpls,N)),M);  %demodulating the data  128 1 
       % or  Demod_Data = pskdemod(  (fft(Hh,N)./(fft(Hpls,N)),M);

       DataEstiRec=[  Demod_Data(9:32,1); Demod_Data(41:64,1); Demod_Data(73:96,1); Demod_Data(105:128,1) ];
       DataTra=[  T_Data(9:32,1); T_Data(41:64,1); T_Data(73:96,1); T_Data(105:128,1) ];
       
       [nErr bErr(i,k)] = symerr(DataEstiRec ,  DataTra);

        
        ms_error_mat=mean(((abs(h-Hpls))/abs(h)).^2);
 for v=1:16
    if(ms_error_mat(v)~=0)
        ms_error=ms_error_mat(v);
    end
end
ls_mse(k,i)=ms_error;
ls_mse_ave=mean(ls_mse);
semilogy(SNR,ls_mse_ave,'--ro')
grid on
 end
 end



% plotting the result

xlabel('SNR in dB')
ylabel('mean squared error')
title('PLOT OF SNR V/S MSE FOR AN OFDM SYSTEM WITH LS ESTIMATOR BASEDRECEIVERS')
figure(2)
semilogy(SNR_V,mean(bErr'),'--or','linewidth',2);
hold on;
EbN0Lin = 10.^(SNR_V/10);
theoryBer = 0.5.*(1-sqrt(EbN0Lin./(EbN0Lin+1)));
semilogy(SNR_V,theoryBer,'--ob','linewidth',2);
legend('simulated','theoritical')
grid on

xlabel('SNR');
ylabel('BER')
title('Bit error probability curve for BPSK using OFDM');

 

[rkmodi1212]

help to obtain this code in time domain

please help me to obtain your code in time domain and frequency domain separately and what is the meaning of Ip+1,Ip+2 in Xp syntax means how it generate that matrix??