load BasicMonostaticRadarExampleData.mat; antenna = phased.IsotropicAntennaElement... ('FrequencyRange',[8e9 12e9],'BackBaffled',true); fc = radiator.OperatingFrequency; c = radiator.PropagationSpeed; lambda = c/fc; ula = phased.ULA('Element',antenna,'NumElements',6,... 'ElementSpacing', lambda/2); pattern(ula,fc,'PropagationSpeed',c,'Type','powerdb') title('6-element Baffled ULA Response Pattern') view(60,50) radiator.Sensor = ula; collector.Sensor = ula; sensormotion = phased.Platform('InitialPosition',[0; 0; 1000]); arrayAxis = [0; 1; 0]; prf = waveform.PRF; vr = ula.ElementSpacing*prf; sensormotion.Velocity = vr/2*arrayAxis; target = phased.RadarTarget('Model','Nonfluctuating','MeanRCS',1, ... 'OperatingFrequency', fc); tgtmotion = phased.Platform('InitialPosition',[1000; 1000; 0],... 'Velocity',[30; 30; 0]); % *Jammer* jammer = phased.BarrageJammer('ERP',100); Fs = waveform.SampleRate; rngbin = c/2*(0:1/Fs:1/prf-1/Fs).'; jammer.SamplesPerFrame = numel(rngbin); jammermotion = phased.Platform('InitialPosition',[1000; 1732; 1000]); %% % *Clutter* clutter = phased.ConstantGammaClutter('Sensor',ula,'SampleRate',Fs,... 'Gamma',-15,'PlatformHeight',1000,... 'OperatingFrequency',fc,... 'PropagationSpeed',c,... 'PRF',prf,... 'TransmitERP',transmitter.PeakPower*db2pow(transmitter.Gain),... 'PlatformSpeed',norm(sensormotion.Velocity),... 'PlatformDirection',[90;0],... 'BroadsideDepressionAngle',0,... 'MaximumRange',37500,... 'AzimuthCoverage',180,... 'PatchAzimuthWidth',10,... 'OutputFormat','Pulses'); % *Propagation tgtchannel = phased.FreeSpace('TwoWayPropagation',true,'SampleRate',Fs,... 'OperatingFrequency', fc); jammerchannel = phased.FreeSpace('TwoWayPropagation',false,... 'SampleRate',Fs,'OperatingFrequency', fc); jammer.SeedSource = 'Property'; jammer.Seed = 5; clutter.SeedSource = 'Property'; clutter.Seed = 5; numpulse = 10; % Number of pulses tsig = zeros(size(rngbin,1), ula.NumElements, numpulse); jsig = tsig; tjcsig = tsig; tcsig = tsig; csig = tsig; for m = 1:numpulse % Update sensor, target and jammer positions [sensorpos,sensorvel] = sensormotion(1/prf); [tgtpos,tgtvel] = tgtmotion(1/prf); [jampos,jamvel] = jammermotion(1/prf); % Calculate the target and jammer angles as seen by the sensor [~,tgtang] = rangeangle(tgtpos,sensorpos); [~,jamang] = rangeangle(jampos,sensorpos); % Simulate propagation of pulse in direction of targets pulse = waveform(); [pulse,txstatus] = transmitter(pulse); pulse = radiator(pulse,tgtang); pulse = tgtchannel(pulse,sensorpos,tgtpos,sensorvel,tgtvel); % Collect target returns at sensor pulse = target(pulse); tsig(:,:,m) = collector(pulse,tgtang); % Collect jammer signal at sensor jamsig = jammer(); jamsig = jammerchannel(jamsig,jampos,sensorpos,jamvel,sensorvel); jsig(:,:,m) = collector(jamsig,jamang); % Collect clutter signal at sensor csig(:,:,m) = clutter(); % Receive collected signals tjcsig(:,:,m) = receiver(tsig(:,:,m)+jsig(:,:,m)+csig(:,:,m),... ~(txstatus>0)); % Target + jammer + clutter tcsig(:,:,m) = receiver(tsig(:,:,m)+csig(:,:,m),... ~(txstatus>0)); % Target + clutter tsig(:,:,m) = receiver(tsig(:,:,m),... ~(txstatus>0)); % Target echo only end %% % *True Target Range, Angle and Doppler* % % The target azimuth angle is 45 degrees, and the elevation angle is about % -35.27 degrees. tgtLocation = global2localcoord(tgtpos,'rs',sensorpos); tgtAzAngle = tgtLocation(1) %% tgtElAngle = tgtLocation(2) %% % The target range is 37000 m. tgtRng = tgtLocation(3) %% % The target Doppler normalized frequency is about 0.21. sp = radialspeed(tgtpos, tgtmotion.Velocity, ... sensorpos, sensormotion.Velocity); tgtDp = 2*speed2dop(sp,lambda); % Round trip Doppler tgtDp/prf ReceivePulse = tjcsig; plot([tgtRng tgtRng],[0 0.01],rngbin,abs(ReceivePulse(:,:,1))); xlabel('Range (m)'), ylabel('Magnitude'); title('Signals collected by the ULA within the first pulse interval') tgtCellIdx = val2ind(tgtRng,c/(2*Fs)); snapshot = shiftdim(ReceivePulse(tgtCellIdx,:,:)); % Remove singleton dim angdopresp = phased.AngleDopplerResponse('SensorArray',ula,... 'OperatingFrequency',fc, 'PropagationSpeed',c,... 'PRF',prf, 'ElevationAngle',tgtElAngle); plotResponse(angdopresp,snapshot,'NormalizeDoppler',true); text(tgtAzAngle,tgtDp/prf,'+ Target') rxmainlobedir = [0; 0]; stapdpca = phased.DPCACanceller('SensorArray',ula,'PRF',prf,... 'PropagationSpeed',c,'OperatingFrequency',fc,... 'Direction',rxmainlobedir,'Doppler',tgtDp,... 'WeightsOutputPort',true) ReceivePulse = tcsig; [y,w] = stapdpca(ReceivePulse,tgtCellIdx); plot([tgtRng tgtRng],[0 1.2e-5],rngbin,abs(y)); xlabel('Range (m)'), ylabel('Magnitude'); title('DPCA Canceller Output (no Jammer)') angdopresp.ElevationAngle = 0; plotResponse(angdopresp,w,'NormalizeDoppler',true); title('DPCA Weights Angle Doppler Response at 0 degrees Elevation') ReceivePulse = tjcsig; [y,w] = stapdpca(ReceivePulse,tgtCellIdx); plot([tgtRng tgtRng],[0 8e-4],rngbin,abs(y)); xlabel('Range (m)'), ylabel('Magnitude'); title('DPCA Canceller Output (with Jammer)') plotResponse(angdopresp,w,'NormalizeDoppler',true); title('DPCA Weights Angle Doppler Response at 0 degrees Elevation') tgtAngle = [tgtAzAngle; tgtElAngle]; stapsmi = phased.STAPSMIBeamformer('SensorArray', ula, 'PRF', prf, ... 'PropagationSpeed', c, 'OperatingFrequency', fc, ... 'Direction', tgtAngle, 'Doppler', tgtDp, ... 'WeightsOutputPort', true,... 'NumGuardCells', 4, 'NumTrainingCells', 100) [y,w] = stapsmi(ReceivePulse,tgtCellIdx); plot([tgtRng tgtRng],[0 2e-6],rngbin,abs(y)); xlabel('Range (m)'), ylabel('Magnitude'); title('SMI Beamformer Output (with Jammer)') plotResponse(angdopresp,w,'NormalizeDoppler',true); title('SMI Weights Angle Doppler Response at 0 degrees Elevation') stapadpca = phased.ADPCACanceller('SensorArray', ula, 'PRF', prf, ... 'PropagationSpeed', c, 'OperatingFrequency', fc, ... 'Direction', rxmainlobedir, 'Doppler', tgtDp, ... 'WeightsOutputPort', true,... 'NumGuardCells', 4, 'NumTrainingCells', 100) [y,w] = stapadpca(ReceivePulse,tgtCellIdx); plot([tgtRng tgtRng],[0 2e-6],rngbin,abs(y)); xlabel('Range (m)'), ylabel('Magnitude'); title(' Output (with Jammer)') plotResponse(angdopresp,w,'NormalizeDoppler',true); title(' Doppler Response at 0 degrees Elevation') displayEndOfDemoMessage(mfilename) updateRate = 0.5; N = 150; airplane1range = 37000; airplane1Azimuth = 60.0; airplane1alt = 6.0e3; airplane1Pos0 = [cosd(airplane1Azimuth)*airplane1range;... sind(airplane1Azimuth)*airplane1range;airplane1alt]; airplane1Vel0 = [400.0;-100.0;-20]; airplane1Accel = [-10.0;0.0;0.0]; airplane1platform = phased.Platform('MotionModel','Acceleration',... 'AccelerationSource','Input port','InitialPosition',airplane1Pos0,... 'InitialVelocity',airplane1Vel0,'OrientationAxesOutputPort',true,... 'InitialOrientationAxes',eye(3)); groundRadarPos = [0,0,0]'; groundRadarVel = [0,0,0]'; groundradarplatform = phased.Platform('MotionModel','Velocity',... 'InitialPosition',groundRadarPos,'Velocity',groundRadarVel,... 'InitialOrientationAxes',eye(3)); groundVehiclePos = [5e3,2e3,0]'; groundVehicleVel = [50,50,0]'; groundvehicleplatform = phased.Platform('MotionModel','Velocity',... 'InitialPosition',groundVehiclePos,'Velocity',groundVehicleVel,... 'InitialOrientationAxes',eye(3)); airplane2Pos = [8.5e3,1e3,6000]'; airplane2Vel = [-300,100,20]'; airplane2platform = phased.Platform('MotionModel','Velocity',... 'InitialPosition',airplane2Pos,'Velocity',airplane2Vel,... 'InitialOrientationAxes',eye(3)); BeamSteering = [0;50]; viewer = phased.ScenarioViewer('BeamRange',8.0e3,'BeamWidth',[2;30],'UpdateRate',updateRate,... 'PlatformNames',{' Radar',' jammer','airplan3','Airplane 2'},'ShowPosition',true,... 'ShowSpeed',true,'ShowAltitude',true,'ShowLegend',true,'ShowRange',true,... 'Title','with jammer Scenario','BeamSteering',BeamSteering); for n = 1:N [groundRadarPos,groundRadarVel] = groundradarplatform(updateRate); [airplane1Pos,airplane1Vel,airplane1Axes] = airplane1platform(updateRate,airplane1Accel); [vehiclePos,vehicleVel] = groundvehicleplatform(updateRate); [airplane2Pos,airplane2Vel] = airplane2platform(updateRate); viewer(groundRadarPos,groundRadarVel,[airplane1Pos,vehiclePos,airplane2Pos],... [airplane1Vel,vehicleVel,airplane2Vel]); BeamSteering = viewer.BeamSteering(1); BeamSteering = mod(BeamSteering + 4,360.0); if BeamSteering > 180.0 BeamSteering = BeamSteering - 360.0; end viewer.BeamSteering(1) = BeamSteering; pause(0.2); end antenna = phased.ULA(4); Fs = 1e6; fc = 1e9; rng('default') waveform = phased.RectangularWaveform('PulseWidth',100e-6,... 'PRF',1e3,'NumPulses',5,'SampleRate',Fs); transmitter = phased.Transmitter('PeakPower',1e4,'Gain',20,... 'InUseOutputPort',true); radiator = phased.Radiator('Sensor',antenna,'OperatingFrequency',fc); jammer = phased.BarrageJammer('ERP',1000,... 'SamplesPerFrame',waveform.NumPulses*waveform.SampleRate/waveform.PRF); target = phased.RadarTarget('Model','Nonfluctuating',... 'MeanRCS',1,'OperatingFrequency',fc); targetchannel = phased.FreeSpace('TwoWayPropagation',true,... 'SampleRate',Fs,'OperatingFrequency', fc); jammerchannel = phased.FreeSpace('TwoWayPropagation',false,... 'SampleRate',Fs,'OperatingFrequency', fc); collector = phased.Collector('Sensor',antenna,... 'OperatingFrequency',fc); amplifier = phased.ReceiverPreamp('EnableInputPort',true); targetloc = [1000 ; 500; 0]; jammerloc = [2000; 2000; 100]; [~,tgtang] = rangeangle(targetloc); [~,jamang] = rangeangle(jammerloc); wav = waveform(); [wav,txstatus] = transmitter(wav); wav = radiator(wav,tgtang); wav = targetchannel(wav,[0;0;0],targetloc,[0;0;0],[0;0;0]); wav = target(wav); wav = collector(wav,tgtang); %% jamsig = jammer(); jamsig = jammerchannel(jamsig,jammerloc,[0;0;0],[0;0;0],[0;0;0]); jamsig = collector(jamsig,jamang); pulsewave = amplifier(wav,~txstatus); pulsewave_jamsig = amplifier(wav + jamsig,~txstatus); %% subplot(2,1,1) t = unigrid(0,1/Fs,size(pulsewave,1)*1/Fs,'[)'); plot(t*1000,abs(pulsewave(:,1))) title('f and t Without Jamming') ylabel('ferquency') subplot(2,1,2) plot(t*1000,abs(pulsewave_jamsig(:,1))) title('f and t with jamming') xlabel('millisec') ylabel('ferquency') prmQPSKTxRx = commqpsktxrx_init % QPSK system parameters useScopes = true; % true if scopes are to be used printReceivedData = false; %true if the received data is to be printed compileIt = false; % true if code is to be compiled useCodegen = false; % true to run the generated mex file if compileIt codegen -report runQPSKSystemUnderTest.m -args {coder.Constant(prmQPSKTxRx),coder.Constant(useScopes),coder.Constant(printReceivedData)} %#ok end if useCodegen BER = runQPSKSystemUnderTest_mex(prmQPSKTxRx, useScopes, printReceivedData); else BER = runQPSKSystemUnderTest(prmQPSKTxRx, useScopes, printReceivedData); end fprintf('Error rate = %f.\n',BER(1)); fprintf('Number of detected errors = %d.\n',BER(2)); fprintf('Total number of compared samples = %d.\n',BER(3)); displayEndOfDemoMessage(mfilename)