Example 05 - Generating a channel impulse response file using a channel model

In this example we use the DLR land-mobile channel model to generate a channel impulse response (CIR) file. This file will contain realistic data of a vehcile driving through an urban environment.

Please go to http://www.kn-s.dlr.de/satnav and click on "get the model" to download the latest model version.

Extract the model to snacs/examples/LandMobileMultipathChannelModel30.

You can use the file generate_DLRLMS_urban_channel_saveas_SNACS_HDF5_CIF_file.m in the example05 directory to generate the CIR file for the SNACS simulation. The file is based on LandMobileMultipathChannel_Demo_UrbanCar.m from the channel model distribution. It is extended by parts to save the resulting CIR as a SNACS CIR input file. These parts are marked by "[for SNACS]" comments.

First, define the parameters:

%%for SNACS, initialize SNACS parameters sec
% write Parameters and ChannelParams to h5 file:
snacs_cir_filename = 'snacs-cir-file_DLRLMS_example-05.h5';
% [for SNACS] map the channel model parameter names to the SNACS CIR parameter names:
simulation_parameters.c0 = Co;
simulation_parameters.cir_rate = Parameters.SampFreq;

At the beginning, we insert one second of a steady line-of-sight signal only to let the loops in the simulation settle:

% write only-LOS signal at the beginning of file.
% length 1 second. 1 s * cir_rate = number of only-LOS cirs:
only_los_offset = 1 * simulation_parameters.cir_rate; 

simulation_parameters.cir_amount = only_los_offset + Parameters.NumberOfSteps;

% [for SNACS] init the SNACS CIR HDF5 file:
init_snacs_cir_file(snacs_cir_filename, simulation_parameters);

% [for SNACS] no channel for the first second to let the loops lock.
% write one second of a A = 1 exp(j*0) component at delay = 0 s to the CIR
% file:
only_los_cir = struct('delay', {0}, 'real', {0}, 'imag', {0});
only_los_cir.delay = 0;
only_los_cir.real = 1;
only_los_cir.imag = 0;
% initialize vectors for the reference distance from transmitter 
% to receiver
ActualPR_vec = zeros(simulation_parameters.cir_amount, 1);
disp('writing only-LOS CIR snapshots...');
for k = 0:only_los_offset-1
    append_one_cir_to_cir_file(snacs_cir_filename, k, only_los_cir);
end
% the next written CIR will start with number: only_los_offset

Later the channel model parameters are written to the CIR file:

%% for SNACS, write channel scenery parameters 
% write channel scenery parameters to HDF5 file:
hdf5write(snacs_cir_filename, '/parameters/model/channel_scenery' , ChannelParams, 'WriteMode', 'append');

The main part follows: Here, the CIR as generated by the channel model is prepared to be written to the HDF5 file:

    %% [for SNACS], write cir to h5 file
    % For that we build one vector from
    % LOS and echoes first:
    % (dhv is the number of the actual CIR snapshot (starting with 1))
    clear echoes
    echoes_elements = numel(LOSDelays)+numel(DelayVec); % amount of echoes
    echo_zeros = zeros(echoes_elements, 1);
    echoes(1:echoes_elements) = struct('delay', {0}, 'real', {0}, 'imag', {0});
    for k = 1:numel(LOSDelays)
        echoes(k).delay = LOSDelays(k);
        echoes(k).real = real(LOS(k));
        echoes(k).imag = imag(LOS(k));
    end
    idx_offset = numel(LOSDelays);
    for k = 1:numel(DelayVec)
        echoes(k+idx_offset).delay = DelayVec(k);
        echoes(k+idx_offset).real = real(ComplexOutputVec(k));
        echoes(k+idx_offset).imag = imag(ComplexOutputVec(k));
    end
    cir_number = only_los_offset + dhv-1; % dhv-1: CIR number for SNACS should start at 0

Finally, the CIR is written to the file:

    % only_los_offset is added because we added 1 second of a only los
    % component to the CIR file, see above.
    append_one_cir_to_cir_file(snacs_cir_filename, cir_number, echoes);
    
    ActualPR_vec(cir_number + 1) = ActualPR; % +1 since Matlab vectors start at 1

    disp(sprintf('written step #: %d, actual range: %f\n', cir_number, ActualPR));

At the end, the reference distances are written to the HDF5 CIR file:

%% [for SNACS], write all reference range values
hdf5write(snacs_cir_filename, '/reference_range/range_absolut', ActualPR_vec, 'WriteMode', 'append');

Also the file LandMobileMultipathChannelModel30//generate.m of the model distribution has to be changed to return the absolute reference range. Append the file by these two lines:

this.LastLOSPR = this.LastLOSPR + (this.SatVector*WayVec');
varargout{9}=this.LastLOSPR;

Use the MATLAB script plot_DLRLMS_CIR_result_example05.m to display the output of the generated CIR:

DLRLMS_cir-result.png
Resulting CIR of the DLR LMS channel model simulation

This is the full MATLAB script:

1 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2 %
3 % This file saves a channel impulse response series as HDF5 file
4 %
5 % based on
6 % LandMobileMultipathChannel_Demo_UrbanCar.m
7 % Version 3.0
8 % distributed by http://www.kn-s.dlr.de/satnav
9 %
10 % adapted by Frank Schubert, 2010
11 %
12 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
13 
14 close all
15 clear all
16 clear classes
17 clc
18 
19 addpath('../helper_functions');
20 
21 % channel model start:
22 pth = '../LandMobileMultipathChannelModel30';
23 addpath(pth);
24 addpath ([pth,'/Knife_Edge']);
25 addpath ([pth,'/Parameter_Files']);
26 
27 Parameters.SampFreq=500; % Hz
28 Parameters.MaximumSpeed=50; % km/h
29 Parameters.SatElevation=30; % Deg
30 Parameters.SatAzimut=-135; % Deg (North == 0, East == 90, South == 180, West == 270)
31 Parameters.NumberOfSteps=10000; % 5000 / 1000 = 1 second
32 Co=2.99e8; % Speed of Light
33 
34 %%for SNACS, initialize SNACS parameters sec
35 % write Parameters and ChannelParams to h5 file:
36 snacs_cir_filename = 'snacs-cir-file_DLRLMS_example-05-2.h5';
37 % [for SNACS] map the channel model parameter names to the SNACS CIR parameter names:
38 simulation_parameters.c0 = Co;
39 simulation_parameters.cir_rate = Parameters.SampFreq;
40 
41 % write only-LOS signal at the beginning of file.
42 % length 1 second. 1 s * cir_rate = number of only-LOS cirs:
43 only_los_offset = 1 * simulation_parameters.cir_rate;
44 
45 simulation_parameters.cir_amount = only_los_offset + Parameters.NumberOfSteps;
46 
47 % [for SNACS] init the SNACS CIR HDF5 file:
48 init_snacs_cir_file(snacs_cir_filename, simulation_parameters);
49 
50 % [for SNACS] no channel for the first second to let the loops lock.
51 % write one second of a A = 1 exp(j*0) component at delay = 0 s to the CIR
52 % file:
53 only_los_cir = struct('delay', {0}, 'real', {0}, 'imag', {0});
54 only_los_cir.delay = 0;
55 only_los_cir.real = 1;
56 only_los_cir.imag = 0;
57 % initialize vectors for the reference distance from transmitter
58 % to receiver
59 ActualPR_vec = zeros(simulation_parameters.cir_amount, 1);
60 disp('writing only-LOS CIR snapshots...');
61 for k = 0:only_los_offset-1
62  append_one_cir_to_cir_file(snacs_cir_filename, k, only_los_cir);
63 end
64 % the next written CIR will start with number: only_los_offset
65 
66 %% channel model continue...
67 FontSize=12;
68 
69 % ---- General Parameters ----
70 
71 ChannelParams.CarrierFreq=1.57542e9; % Hz
72 ChannelParams.SampFreq=Parameters.SampFreq;
73 ChannelParams.EnableDisplay=0; % 3D visualization is not available in the free version
74 ChannelParams.EnableCIRDisplay=0; % enables CIR display
75 
76 % ---- Mode Parameters ----
77 
78 ChannelParams.UserType = 'Car';
79 ChannelParams.Surrounding = 'Urban';
80 ChannelParams.AntennaHeight = 2; % m Height of the Antenna
81 ChannelParams.MinimalPowerdB=-40; % Echos below this Limit are not initialised
82 
83 % ---- UserParameters ---
84 
85 ChannelParams.DistanceFromRoadMiddle=-5; % negative: continental (right), positive: England (left)
86 
87 % ---- Graphics Parameters ---
88 
89 ChannelParams.GraphicalPlotArea=50; %
90 ChannelParams.ViewVector = [-60,20]; % 3D visualization is not available in the free version
91 ChannelParams.RoadWidth = 15; %
92 
93 % --- Building Params ---
94 
95 ChannelParams.BuildingRow1=1; % logigal to switch Building Row right(heading 0 deg) on
96 ChannelParams.BuildingRow2=1; % logigal to switch Building Row left (heading 0 deg) on
97 ChannelParams.BuildingRow1YPosition=-12; % m
98 ChannelParams.BuildingRow2YPosition=12; % m
99 
100 ChannelParams.HouseWidthMean=22; % m
101 ChannelParams.HouseWidthSigma=25; % m
102 ChannelParams.HouseWidthMin=10; % m
103 ChannelParams.HouseHeightMin=4; % m
104 ChannelParams.HouseHeightMax=50; % m
105 ChannelParams.HouseHeightMean=16; % m
106 ChannelParams.HouseHeightSigma=6.4; % m
107 ChannelParams.GapWidthMean=27; % m
108 ChannelParams.GapWidthSigma=25; % m
109 ChannelParams.GapWidthMin=10; % m
110 ChannelParams.BuildingGapLikelihood=0.18;% lin Value
111 
112 % --- Tree Params ---
113 
114 ChannelParams.TreeHeight = 8; % m
115 ChannelParams.TreeDiameter = 5; % m
116 ChannelParams.TreeTrunkLength=2; % m
117 ChannelParams.TreeTrunkDiameter=.2; % m
118 
119 ChannelParams.TreeAttenuation = 1.1; % dB/m
120 
121 ChannelParams.TreeRow1Use=1; % logical switches tree row 1 on
122 ChannelParams.TreeRow2Use=1; % logical switches tree row 2 on
123 
124 ChannelParams.TreeRow1YPosition=-8; % m
125 ChannelParams.TreeRow2YPosition=8; % m
126 
127 ChannelParams.TreeRow1YSigma=2; % m
128 ChannelParams.TreeRow2YSigma=2; % m
129 
130 ChannelParams.TreeRow1MeanDistance=60; % m
131 ChannelParams.TreeRow2MeanDistance=40; % m
132 
133 ChannelParams.TreeRow1DistanceSigma=20; % m
134 ChannelParams.TreeRow2DistanceSigma=20; % m
135 
136 % --- Pole Params ---
137 
138 ChannelParams.PoleHeight = 10; % m
139 ChannelParams.PoleDiameter = .2; % m
140 
141 ChannelParams.PoleRow1Use=1; % logical switches Pole row 1 on
142 ChannelParams.PoleRow2Use=0; % logical switches Pole row 2 on
143 
144 ChannelParams.PoleRow1YPosition=0; % m
145 ChannelParams.PoleRow2YPosition=10; % m
146 
147 ChannelParams.PoleRow1YSigma=1; % m
148 ChannelParams.PoleRow2YSigma=1; % m
149 
150 ChannelParams.PoleRow1MeanDistance=25; % m
151 ChannelParams.PoleRow2MeanDistance=10; % m
152 
153 ChannelParams.PoleRow1DistanceSigma=10; % m
154 ChannelParams.PoleRow2DistanceSigma=10; % m
155 
156 %% for SNACS, write channel scenery parameters
157 % write channel scenery parameters to HDF5 file:
158 hdf5write(snacs_cir_filename, '/parameters/model/channel_scenery' , ChannelParams, 'WriteMode', 'append');
159 
160 %% channel model continue...
161 
162 % ------------ Initial Settings -------------
163 % - Anything Below here must not be changed -
164 % -------------------------------------------
165 
166 MaximumPossibleSpeed=Co*Parameters.SampFreq/ChannelParams.CarrierFreq/2; % To fulfill the sampling Theorem
167 SamplingTime=1/Parameters.SampFreq;
168 
169 % --- Initialising the channel object ---
170 
171 pause(1)
172 disp('Initialising the channel ...')
173 TheChannelObject=LandMobileMultipathChannel(ChannelParams);
174 
175 TimeVec=0;
176 ComplexOutputVec=[];
177 
178 % --- Specify power and delay bins for output statistics ---
179 
180 pwrvec = [0:-1:-30]; % power bins in dB
181 dlyvec = [0:10e-9:500e-9]; % delay bins in s
182 
183 PowerDelayProfile(1:length(pwrvec),1:length(dlyvec)) = 0; % allocate memory
184 pwrstp = (pwrvec(end)-pwrvec(1))/(length(pwrvec)-1); % get step size
185 dlystp = (dlyvec(end)-dlyvec(1))/(length(dlyvec)-1); % get step size
186 
187 % --- start simulation ---
188 
189 h = waitbar(0,'Simulation running ...');
190 
191 if ChannelParams.EnableCIRDisplay
192 
193  % --- init CIR figure ---
194 
195  hh = figure;
196  subplot(211)
197  xlabel('Delay in s')
198  ylabel('Power in dB')
199  axis([-2e-8,40e-8,0,50])
200  set(get(hh,'Children'),'YTickLabel',[-40 -30 -20 -10 0 10]);
201  hold on
202  grid on
203  plot (0,0,'r')
204  plot (0,0)
205  legend ('Direct paths','Echo paths')
206 
207  subplot(212)
208  xlabel('Delay in s')
209  ylabel('Phase in rad')
210  axis([-2e-8,40e-8,-pi,pi])
211  hold on
212  grid on
213 
214 end
215 
216 ActualPR = 0;
217 
218 for dhv=1:Parameters.NumberOfSteps
219 
220  TimeVec(end)=dhv/Parameters.SampFreq;
221 
222  % --- "drunken" driver movement example ---
223 
224  ActualSpeed=Parameters.MaximumSpeed/2/3.6*(1+sin(TimeVec(end)/3));
225  SpeedVec(dhv)=ActualSpeed; % m/s
226  ActualHeading=20*sin(TimeVec(end)/3); % Deg (North == 0, East == 90, South == 180, West == 270)
227 
228  % --- generate CIR ---
229 
230  [TheChannelObject,LOS,LOSDelays,ComplexOutputVec,DelayVec,EchoNumberVec,WayVec(dhv),TimeVec(dhv),ActualPR_change]=generate(TheChannelObject,ActualSpeed,ActualHeading,Parameters.SatElevation,Parameters.SatAzimut);
231 
232  waitbar(dhv/Parameters.NumberOfSteps,h)
233 
234  % --- binning LOS ---
235 
236  for sfg = 1:length(LOSDelays)
237 
238  dlybin = round(LOSDelays(sfg)/dlystp) + 1;
239  pwrbin = round(20*log10(abs(LOS(sfg)))/pwrstp) + 1;
240 
241  if pwrbin<=length(pwrvec) & pwrbin>0 & dlybin<=length(dlyvec)
242  PowerDelayProfile(pwrbin,dlybin) = PowerDelayProfile(pwrbin,dlybin) + 1;
243  end
244 
245  end
246 
247  % --- binning echoes ---
248 
249  for sfg = 1:length(DelayVec)
250 
251  dlybin = round(DelayVec(sfg)/dlystp) + 1;
252  pwrbin = round(20*log10(abs(ComplexOutputVec(sfg)))/pwrstp) + 1;
253 
254  if pwrbin<=length(pwrvec) & pwrbin>0 & dlybin<=length(dlyvec)
255  PowerDelayProfile(pwrbin,dlybin) = PowerDelayProfile(pwrbin,dlybin) + 1;
256  end
257 
258  end
259 
260  if ChannelParams.EnableCIRDisplay
261 
262  % --- display CIR ---
263 
264  figure(hh);
265  subplot(211)
266  cla
267  Time = dhv/Parameters.SampFreq;
268  title(['Channel Impulse Response (CIR), T = ',num2str(Time,'%5.2f'),' s, v = ',num2str(ActualSpeed*3.6,'%4.1f'),' km/h'])
269  stem(LOSDelays,40 + 20*log10(abs(LOS)),'r');
270  stem(DelayVec,40 + 20*log10(abs(ComplexOutputVec)));
271 
272  subplot(212)
273  cla
274  stem(LOSDelays,angle(LOS),'r');
275  stem(DelayVec,angle(ComplexOutputVec));
276 
277  end
278 
279  %% [for SNACS], write cir to h5 file
280  % For that we build one vector from
281  % LOS and echoes first:
282  % (dhv is the number of the actual CIR snapshot (starting with 1))
283  clear echoes
284  ActualPR = ActualPR + ActualPR_change;
285  ActualLOSDelay = ActualPR / Co;
286  echoes_elements = numel(LOSDelays)+numel(DelayVec); % amount of echoes
287  echo_zeros = zeros(echoes_elements, 1);
288  echoes(1:echoes_elements) = struct('delay', {0}, 'real', {0}, 'imag', {0});
289  for k = 1:numel(LOSDelays)
290  echoes(k).delay = LOSDelays(k) + ActualLOSDelay;
291  echoes(k).real = real(LOS(k));
292  echoes(k).imag = imag(LOS(k));
293  end
294  idx_offset = numel(LOSDelays);
295  for k = 1:numel(DelayVec)
296  echoes(k+idx_offset).delay = DelayVec(k) + ActualLOSDelay;
297  echoes(k+idx_offset).real = real(ComplexOutputVec(k));
298  echoes(k+idx_offset).imag = imag(ComplexOutputVec(k));
299  end
300  cir_number = only_los_offset + dhv-1; % dhv-1: CIR number for SNACS should start at 0
301 
302  % only_los_offset is added because we added 1 second of a only los
303  % component to the CIR file, see above.
304  append_one_cir_to_cir_file(snacs_cir_filename, cir_number, echoes);
305 
306  ActualPR_vec(cir_number + 1) = ActualPR; % +1 since Matlab vectors start at 1
307 
308  disp(sprintf('written step #: %d, actual range: %f\n', cir_number, ActualPR));
309 
310  %% channel model continue...
311 
312 end%
313 
314 %% [for SNACS], write all PR values
315 hdf5write(snacs_cir_filename, '/reference_range/range_absolut', ActualPR_vec, 'WriteMode', 'append');
316 
317 %% channel model continue...
318 
319 close(h);
320 
321 % --- calculate probability density function ---
322 
323 PowerDelayProfile = PowerDelayProfile/sum(sum(PowerDelayProfile));
324 
325 % --- display PowerDelayProfile ---
326 
327 figure
328 surf(dlyvec,pwrvec,10*log10(PowerDelayProfile+eps),'LineStyle','none','FaceColor','interp','EdgeLighting','phong');
329 caxis([-70,-10]);
330 view(2)
331 xlabel('delay in s')
332 ylabel('power in dB')
333 title('Power delay profile - probability density function')
334 hc = colorbar;
335 set(hc,'YLim',[-70,-10])
336 clear newYTic YTic
337 YTic = get(hc,'YTickLabel');
338 axes(hc);
339 ax = axis;
340 dta = (ax(3)-ax(4))/(size(YTic,1)-1);
341 for kk = 1:size(YTic,1)
342  oldYTic(kk,:) = [' '];
343  newYTic(kk,:) = ['10^',num2str(str2num(YTic(kk,:))/10),''];
344  text(1.2,ax(3)-dta*(kk-1),['10^{',num2str(str2num(YTic(kk,:))/10),'}'],'interpreter','tex','horizontalalignment','left','verticalalignment','middle');
345 end
346 set(hc,'YTickLabel',oldYTic);
347 
348 % ---------------------------------
349 
350 disp(' ');
351 disp('Simulation finished');
352