finish coding and testing of online acq; add sample data

串行在线捕获程序已经完成编写,并进行了仿真和测试(对飞和814均可用);增加了在电子所使用转台测试的幅相采集数据
2023-05-08 19:13:53 +08:00
parent e07389e326
commit 2d4650a3b6
2835 changed files with 408373 additions and 181193 deletions
@@ -0,0 +1,44 @@
+close all;
+
+% 4 ch
+fileNameList = [...
+"iladata_duifei_1t2r_48msps_bw15m_trx_ext_lo_noise_1-1.csv",...
+"iladata_duifei_1t2r_48msps_bw15m_trx_ext_lo_noise_1-2.csv",...
+"iladata_duifei_1t2r_48msps_bw15m_trx_ext_lo_noise_1-3.csv",...
+"iladata_duifei_1t2r_48msps_bw15m_trx_ext_lo_noise_2-1.csv",...
+"iladata_duifei_1t2r_48msps_bw15m_trx_ext_lo_noise_2-2.csv",...
+"iladata_duifei_1t2r_48msps_bw15m_trx_ext_lo_noise_2-3.csv",...
+"iladata_duifei_1t2r_48msps_bw15m_trx_ext_lo_noise_3-1.csv",...
+"iladata_duifei_1t2r_48msps_bw15m_trx_ext_lo_noise_3-2.csv",...
+"iladata_duifei_1t2r_48msps_bw15m_trx_ext_lo_noise_3-3.csv"
+];
+
+[~,fileNum] = size(fileNameList);
+
+maxList = strings(fileNum+1,5);
+stdList = strings(fileNum+1,5);
+
+maxList(1,:) = ["Max Value", "ch1", "ch2", "ch3", "ch4"];
+maxList(2:end,1) = fileNameList;
+stdList(1,:) = ["Std Value", "ch1", "ch2", "ch3", "ch4"];
+stdList(2:end,1) = fileNameList;
+
+for i = 1:fileNum
+  rxRaw = readmatrix(fileNameList(:,i));
+  [rowNum,colNum] = size(rxRaw);
+
+  if (colNum == 11)
+    rx = [rxRaw(:, 4)+1i*rxRaw(:, 5),rxRaw(:, 6)+1i*rxRaw(:, 7),...
+          rxRaw(:, 8)+1i*rxRaw(:, 9),rxRaw(:,10)+1i*rxRaw(:,11)];
+  elseif (colNum == 7)
+    rx = [rxRaw(:, 4)+1i*rxRaw(:, 5),rxRaw(:, 6)+1i*rxRaw(:, 7),...
+          zeros(rowNum,1),zeros(rowNum,1)];
+  else
+    disp('Invalid Data Format.')
+    return
+  end
+
+  maxList(i+1,2:end) = string(max(rx));
+  stdList(i+1,2:end) = string(std(rx));
+
+end
@@ -0,0 +1,86 @@
+%% File Header/文件说明
+%
+% Beamforming with phased array derives signal characteristics from phase
+% differences between sensors. Therefore the amplitude-phase consistence is
+% vital to phased systems.
+%
+% The script analyzes the amp-phase relations of ILA captured data, the format
+% of which is square-shaped pseudorandorm bit sequence. Each channel is
+% processed seperately with different frequency shifts and different optimum
+% sampling points. The estimation of amplitude and phase are based on
+% acculumated correlation peaks.
+%
+clc; clear;
+%% initial parameters/初始参数
+interp = 8; % interpolation factor
+channelNum = 4;
+dataLen = 131072;
+fs = 48e6;
+dataFile = "iladata_duifei_1t2r_48msps_bw15m_rx_ext_lo_bpsk_2-1.csv";
+
+% standard pn sequence
+pnSequence = comm.PNSequence('Polynomial','x^7+x^6+1', ...
+  'InitialConditions',[zeros(1,6),1],'SamplesPerFrame',2^7-1);
+PRBS7 = pnSequence();
+PRBS7_polar = 1-2*[PRBS7;0];
+
+chipLen = length(PRBS7_polar);    % pn chip length
+symNum = dataLen/chipLen/interp;  % total symbol number
+symNumPerAccum = 1;               % symbol nbumber of each accumulation
+accNum = symNum/symNumPerAccum;   % total accumulation numbers, i.e. peak numbers
+
+%% data analysis
+
+% read data, and swap ch1 and ch2, because loop data locate at ch2 and ch4
+rxdata = readmatrix(dataFile);
+seqRx(:,1) = (rxdata(:, 6) + 1i*rxdata(:, 7)).';
+seqRx(:,2) = (rxdata(:, 4) + 1i*rxdata(:, 5)).';
+seqRx(:,3) = (rxdata(:, 8) + 1i*rxdata(:, 9)).';
+seqRx(:,4) = (rxdata(:,10) + 1i*rxdata(:,11)).';
+
+% compensate frequecny deviation and downsample at optimum sample position based
+% on the varience of constellation
+[seqRxDown,fd,samplePos] = suqareShapeFreqEst(seqRx,interp,fs,-50,50,5);
+
+% split [seqRxDown] into multiple groups, each group contains several symbols to
+% be correlated with PRBS7 and accumulated
+seqRxDownSplit = reshape(seqRxDown(1:dataLen/interp,:),chipLen*symNumPerAccum, ...
+  accNum,channelNum);
+
+seqRxAcc = zeros(chipLen,accNum,channelNum);
+desiredPnt = zeros(accNum,channelNum);
+ampEstSym = zeros(accNum,channelNum);
+angleEstSym = zeros(accNum,channelNum);
+
+for i = 1:accNum
+  for j = 1:channelNum
+    % correlation & multi-symbol accumulation
+    groupCorrAcc = sum(reshape(conv(flip(seqRxDownSplit(:,i,j)), ...
+      PRBS7_polar,'same'),chipLen,symNumPerAccum),2)';
+    seqRxAcc(:,i,j) = groupCorrAcc;
+    % peak position
+    desiredPnt(i,j) = groupCorrAcc(abs(groupCorrAcc)==max(abs(groupCorrAcc)));
+  end
+end
+
+for i = 1:channelNum
+  ampEstSym(:,i) = abs(desiredPnt(:,i))./abs(desiredPnt(:,1));
+  angleEstSym(:,i) = mod((angle(desiredPnt(:,i))-angle(desiredPnt(:,1)))/pi*180,360);
+end
+
+%% figures
+figure(1);
+
+subplot(121);
+plot(ampEstSym(:,[1,4]));
+legend({"ch1","ch4"});
+title("Amplitude Estimation");
+ylim([min(ampEstSym(:,[1,4]),[],'all')-0.1,max(ampEstSym(:,[1,4]),[],'all')+0.1]);
+grid on; grid minor;
+
+subplot(122);
+plot(angleEstSym(:,[1,4]));
+legend({"ch1","ch4"});
+title("Angle Estimation");
+ylim([min(angleEstSym(:,[1,4]),[],'all')-1,max(angleEstSym(:,[1,4]),[],'all')+1]);
+grid on; grid minor;
@@ -0,0 +1,50 @@
+function [dataOut, fdTarget, samplePosTarget] = suqareShapeFreqEst(dataIn, interp, fs, fmin, fmax, fd)
+%suqareShapeFreqEst - 利用星座图的集中程度在一定频率范围内使用频率步进估计方波成形信号的频偏
+%
+% data - input complex data, row vector
+% dataLen - data length
+% interp - interpolation rate
+% fs - sampling rate, in Herts
+% fmin - minimum estimation frequency, in Herts
+% fmax - maximum estimation frequency, in Herts
+% fd - estimation frequency step, in Herts
+  
+% dimension
+[dataLen,channelNum] = size(dataIn);
+
+freqList = fmin:fd:fmax; 
+freqNum = length(freqList);
+varRecord = zeros(freqNum,interp);
+
+dataOut = zeros(dataLen/interp,channelNum);
+fdTarget = zeros(channelNum,1);
+samplePosTarget = zeros(channelNum,1);
+
+for channel = 1:channelNum
+
+  for freq = 1:freqNum
+
+    % carrier
+    cxr = exp(-1i*2*pi*freqList(freq)/fs*(0:dataLen-1)');
+
+    % apply frequency shift onto input data
+    seqRx = dataIn(:,channel).*cxr;
+  
+    for samplePos = 0:interp-1
+      seqRxDown = downsample(seqRx,interp,samplePos); % 下采样
+      seqRxDownQuaterOne = abs(real(seqRxDown)) + 1i*abs(imag(seqRxDown)); % 集中到第一象限
+      varRecord(freq,samplePos+1) = var(seqRxDownQuaterOne);
+    end
+  
+  end
+
+  [fdIndex, samplePosIndex] = find(varRecord==min(min(varRecord)));
+  fdTarget(channel) = freqList(fdIndex);
+  samplePosTarget(channel) = samplePosIndex-1;
+
+  cxrTarget = exp(-1i*2*pi*fdTarget(channel)/fs*(0:dataLen-1)');
+  seqRxTarget = dataIn(:,channel).*cxrTarget;
+  dataOut(:,channel) = downsample(seqRxTarget,interp,samplePosTarget(channel));
+end
+
+end