Skip to content

Commit fa9fc54

Browse files
kyruzickyruzic
committed
Adding all files
1 parent 42f5f79 commit fa9fc54

17 files changed

+2103
-0
lines changed

‎README.txt

Lines changed: 45 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,45 @@
1+
3D Model of High Frequency Radar Systems Propagated Through the Ionosphere
2+
3+
Author:
4+
Kyle Ruzic
5+
6+
Date:
7+
August 30th 2018
8+
9+
========================================================================
10+
11+
This package generates a 3D model of SuperDARN Saskatoon's gain pattern propagated through the International Reference Ionosphere (IRI), an empirical model of the ionosphere.
12+
13+
========================================================================
14+
15+
Motivation
16+
17+
========================================================================
18+
19+
This project's motivation was to create a model that could be used to better understand high frequency (HF) radiowave propagation in the ionosphere. Specifically the goal was to create a model that could be used to compare to the data received by the radio receiver instrument (RRI), which is one of the eight scientific instruments that make up the enhanced Polar Outflow Probe (e-POP). The ability to compare the model to the data received by RRI will give us a better understanding of the data, and allow us to further explore regions of the ionosphere where we find discrepancies between the modelled data and the data received by RRI.
20+
21+
========================================================================
22+
23+
Methods
24+
25+
========================================================================
26+
27+
The first step in generating the model is to generate an ionospheric model, this is done using the IRI which is produced for a specific time, date, and region using a PHaRLAP routine. The next step in order to create the 3D model uses a technique know as raytracing, which is a computational efficient way of computing the path of light propagating through a medium with a non-uniform index of refraction, such as the ionosphere in the HF radiowave spectrum. To compute the path of the rays we use the Provision of High-frequency Raytracing LAboratory for Propagation studies (PHaRLAP), which is required to use this package but not included. To generate the model we create a 3D grid which represents latitudes, longitudes, and altitudes. We then trace the paths of roughly 100000 rays, then we implement a 3D binning algorithm which takes each point along the ray's path and places the power of the light in the bin corresponding to the points latitude, longitude, and altitude. To compute the power at each point of the ray we use the radar equation
28+
29+
P_r = P_s / (4 * pi * R^2)
30+
31+
Where P_s is the power of the radar at the source, and R is the length of the path that the ray has taken from the source. After this is completed for every ray the generation of the model is finished. In order to make comparisons of to RRI data an interpolant of the model needs to be made in order to compare the high sampling frequency of RRI, to the relatively low resolution of the model. The points along the path of RRI are queried against the interpolant and then that can be compared to RRI's data.
32+
33+
========================================================================
34+
35+
How to use?
36+
37+
========================================================================
38+
39+
The basics of generating the model involve first generating the ionospheric grid, then calling the raytracer and the binning function. An example of how this would look is
40+
41+
[iono_struct, general_struct] = gen_iono_ns(UT)
42+
rad_grid = rayCaller(dimensions, iono_struct, general_struct)
43+
44+
iono_struct is the structure containing the ionospheric grid and other necessary information, general_struct is the structure containing the general information that is necessary for ray tracing, and rad_grid is the generated 3D model. More information regarding these can be in found the documentation for these files. Examples of use are included.
45+

‎example.m

Lines changed: 63 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,63 @@
1+
% Purpose:
2+
% Example of how to generate model from start to finish. This
3+
% is intended to be a basic use example, more information
4+
% regarding the functions used in this example can be found in
5+
% their respective files.
6+
%
7+
% Author:
8+
% Kyle Ruzic
9+
%
10+
% Date:
11+
% August 30th 2018
12+
13+
14+
15+
16+
% general_struct - contains general information that is used
17+
% for ray tracing
18+
% .UT - 5x1 array containing UTC date and time - year, month,
19+
% day, hour, minute
20+
% .speed_of_light - speed of light
21+
% .re - radius of the Earth
22+
% .R12 - R12 index
23+
% .origin_lat - origin latitude of rays to be
24+
% traced, this should be the location
25+
% of the radar.
26+
% .origin_long - origin longitude
27+
% .origin_ht - origin height
28+
% .gain_dat - gain pattern of radar
29+
30+
general_struct.UT = [2015 5 7 8 0]; % UT - [year month day hour minute]
31+
general_struct.speed_of_light = 2.99792458e8;
32+
general_struct.re = 6376000;
33+
general_struct.R12 = 100;
34+
general_struct.origin_lat = 52.16;
35+
general_struct.origin_long = -106.52;
36+
general_struct.origin_ht = 0.0;
37+
general_struct.gain_dat = getGain('dat/SuperDARN_sas_11MHz_boresight.dat');
38+
39+
40+
% Set up dimension information for the model
41+
%
42+
% .range - [minLat, maxLat, minLon, maxLon, minAlt, maxAlt]
43+
% .spacing - [numLatBins, numLonBins, numAltBins], the dimensions of
44+
% the model
45+
46+
dimensions.range = [122, 172, 24, 124, 0, 1600];
47+
dimensions.spacing = [200, 400, 200];
48+
49+
% Generate and then save ionospheric grid, and general information
50+
% structure
51+
52+
iono_struct = gen_iono_ns(general_struct.UT);
53+
saveIonoGrid(iono_struct, general_struct)
54+
55+
% Generate and then save model and dimensions structure
56+
57+
radGrid = rayCaller_ns(dimensions, iono_struct, general_struct);
58+
saveRadGrid(radGrid, dimensions, general_struct.UT)
59+
60+
% Create plots of the generated model
61+
simplePlotter(radGrid, dimensions, UT)
62+
63+

‎gen_iono_ns.m

Lines changed: 156 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,156 @@
1+
function [iono_struct, general_struct] = gen_iono_ns(UT)
2+
% Name:
3+
% gen_iono_ns
4+
%
5+
% Author:
6+
% Kyle Ruzic
7+
%
8+
% Date:
9+
% August 24th 2018
10+
%
11+
% Purpose:
12+
% This program is largely
13+
% based on the beginning of this PHaRLAP example,
14+
% "ray_test_3d_sp.m" The purpose of this is to generate a IRI
15+
% ionospheric grid, that is used by PHaRLAP to calculate the
16+
% path of the rays. A version of this also exists that saves
17+
% the generated ionosphere which allows for it to be easily
18+
% reloaded to prevent having to regenerate it every time you
19+
% make a new model for the same date and time.
20+
%
21+
% Inputs:
22+
% UT - Vector containing date and time information, must be
23+
% formatted like [YYYY MM DD HH MM]
24+
%
25+
% Outputs:
26+
% iono_struct
27+
% .pf_grid - 3d grid (height vs lat. vs lon.) of ionospheric plasma
28+
% frequency (MHz)
29+
% .pf_grid_5 - 3d grid (height vs lat. vs lon.) of ionospheric plasma
30+
% frequency (MHz) 5 minutes later
31+
% .collision_freq - 3d grid (height vs lat. vs lon.) of ionospheric
32+
% collision frequencies
33+
% .Bx - 3d grid of x component of geomagnetic field
34+
% .By - 3d grid of y component of geomagnetic field
35+
% .Bz - 3d grid of z component of geomagnetic field
36+
%
37+
% .iono_grid_parms - 9x1 vector containing the parameters which define the
38+
%
39+
% ionospheric grid :
40+
% (1) geodetic latitude (degrees) start
41+
% (2) latitude step (degrees)
42+
% (3) number of latitudes
43+
% (4) geodetic longitude (degrees) start
44+
% (5) lonitude step (degrees)
45+
% (6) number of longitudes
46+
% (7) geodetic height (km) start
47+
% (8) height step (km)
48+
% (9) number of heights
49+
%
50+
% .geomag_grid_parms - 9x1 vector containing the parameters which define the
51+
%
52+
% geomagnetic grid :
53+
% (1) geodetic latitude (degrees) start
54+
% (2) latitude step (degrees)
55+
% (3) number of latitudes
56+
% (4) geodetic lonitude (degrees) start
57+
% (5) lonitude step (degrees)
58+
% (6) number of longitudes
59+
% (7) geodetic height (km) start
60+
% (8) height step (km)
61+
% (9) number of heights
62+
%
63+
64+
65+
66+
67+
68+
69+
if ~exist('minute', 'var')
70+
minute = 0;
71+
end
72+
73+
UT = UT
74+
speed_of_light = 2.99792458e8;
75+
re = 6376000; % Radius of Earth in m
76+
R12 = 100;
77+
origin_lat = 52.16; % Location of SuperDARN saskatoon
78+
origin_long = -106.52;
79+
origin_ht = 0.0;
80+
doppler_flag = 1;
81+
82+
%ionospheric grid setup
83+
84+
ht_start = 50; % start height for ionospheric grid (km)
85+
ht_inc = 6; % height increment (km)
86+
num_ht = 301.0;
87+
lat_start = 30.0;
88+
lat_inc = 0.086;
89+
num_lat = 701.0;
90+
lon_start= -156.0;
91+
lon_inc = 0.15;
92+
num_lon = 701.0;
93+
iono_grid_parms = [lat_start, lat_inc, num_lat, lon_start, lon_inc, num_lon, ...
94+
ht_start, ht_inc, num_ht];
95+
96+
B_ht_start = ht_start; % start height for geomagnetic grid (km)
97+
B_num_ht = 201; %ceil(num_ht * ht_inc/B_ht_inc); % Max that can be used is 201
98+
B_ht_inc = (ht_inc*(num_ht-1))/(B_num_ht-1); % height increment (km)
99+
B_lat_start = lat_start;
100+
B_num_lat = 101; % max value that can be used is 101
101+
B_lat_inc = (lat_inc*(num_lat-1))/(B_num_lat-1);
102+
B_lon_start = lon_start;
103+
B_num_lon = 101; % max value that can be used is 101
104+
B_lon_inc = (lon_inc*(num_lon-1))/(B_num_lon-1);
105+
geomag_grid_parms = [B_lat_start, B_lat_inc, B_num_lat, B_lon_start, ...
106+
B_lon_inc, B_num_lon, B_ht_start, B_ht_inc, B_num_ht];
107+
108+
109+
% this block of code ensures that the size of both the iono_grid
110+
% and geomagnetic grid are equal
111+
112+
kill = false; % if they aren't equal terminate the program and print
113+
% an error message
114+
if (B_ht_start + (B_num_ht-1)*B_ht_inc) ~= (ht_start + (num_ht-1)*ht_inc)
115+
kill = true;
116+
fprintf("Inconsistently sized Magnetic and Ionospheric grids (height)");
117+
tot_B_height = (B_ht_start + (B_num_ht-1)*B_ht_inc)
118+
tot_iono_height = (ht_start + (num_ht-1)*ht_inc)
119+
elseif (B_lat_start + (B_num_lat-1)*B_lat_inc) ~= (lat_start + (num_lat-1)*lat_inc)
120+
kill = true;
121+
fprintf("Inconsistently sized Magnetic and Ionospheric grids (lat)");
122+
tot_B_lat = (B_lat_start + (B_num_lat-1)*B_lat_inc)
123+
tot_iono_lat = (lat_start + (num_lat-1)*lat_inc)
124+
elseif (B_lon_start + (B_num_lon-1)*B_lon_inc) ~= (lon_start + (num_lon-1)*lon_inc)
125+
kill = true;
126+
fprintf("Inconsistently sized Magnetic and Ionospheric grids (lon)");
127+
tot_B_lon = (B_lon_start + (B_num_lon-1)*B_lon_inc)
128+
tot_iono_lon = (lon_start + (num_lon-1)*lon_inc)
129+
end
130+
131+
if kill
132+
fprintf('Size mismatch of iono and geomagnetic grids')
133+
return
134+
end
135+
136+
%Generates the model ionosphere
137+
tic
138+
fprintf('Generating ionospheric and geomag grids... ')
139+
[iono_pf_grid, iono_pf_grid_5, collision_freq, Bx, By, Bz] = ...
140+
gen_iono_grid_3d(UT, R12, iono_grid_parms, ...
141+
geomag_grid_parms, doppler_flag, 'iri2016');
142+
toc
143+
fprintf('\n')
144+
145+
%Creates a structure with the outputs
146+
147+
iono_struct.pf_grid = iono_pf_grid;
148+
iono_struct.pf_grid_5 = iono_pf_grid_5;
149+
iono_struct.collision_freq = collision_freq;
150+
iono_struct.Bx = Bx;
151+
iono_struct.By = By;
152+
iono_struct.Bz = Bz;
153+
iono_struct.iono_grid_parms = iono_grid_parms;
154+
iono_struct.geomag_grid_parms = geomag_grid_parms;
155+
156+

‎getGain.m

Lines changed: 22 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,22 @@
1+
function gainDat = getGain(pathToGain)
2+
% Name:
3+
% getGain
4+
%
5+
% Author:
6+
% Kyle Ruzic
7+
%
8+
% Date:
9+
% August 22nd 2018
10+
%
11+
% Purpose:
12+
% load gain from data file this is a seperate function to make
13+
% it easier to change how this is handled in the future, so
14+
% that other radars gain patterns can be used
15+
%
16+
% Inputs:
17+
% pathToGain - path to where gain file is saved
18+
%
19+
% Outputs:
20+
% gainDat - array of gain pattern
21+
22+
gainDat = load(pathToGain);

‎gridInterp.m

Lines changed: 60 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,60 @@
1+
function [pathData, pathVec] = gridInterp(path, radGrid, dimensions)
2+
% Name:
3+
% gridInterp
4+
%
5+
% Author:
6+
% Kyle Ruzic
7+
%
8+
% Date:
9+
% August 24th 2018
10+
%
11+
% Purpose:
12+
% Creates an interpolant grid based upon the modelled grid
13+
% passed to it, and queries the points in the interpolant
14+
% structure defined by the path. The path is transformed from
15+
% units of lat/lon/alt to 'grid' units.
16+
%
17+
% Inputs:
18+
% path - path of cassiope loaded from data file
19+
% radGrid - generated model
20+
% dimensions - A structure containing the dimensions for which the model
21+
% will be produced, it should be in this form:
22+
%
23+
% dimensions.range = [minLat, maxLat, minLon, maxLon, minAlt, maxAlt];
24+
% dimensions.spacing = [numLat, numLon, numAlt];
25+
% .range - Contains the ranges for which the model will
26+
% be generated
27+
% .spacing - more aptly named size, contains sizes of
28+
% each dimension
29+
30+
31+
radGrid = radGrid(:,:,:,1); % removing the second dimension as
32+
% it isn't used currently, when
33+
% pointing direction of ray is
34+
% included this line will need to
35+
% be removed
36+
37+
38+
% puts the size of the model grid into a grid format, this is
39+
% needed to be done in order to create the interpolant structure
40+
latVec = [1:size(radGrid, 1)];
41+
lonVec = [1:size(radGrid, 2)];
42+
altVec = [1:size(radGrid, 3)];
43+
44+
% creates grid of size data for generating interpolant
45+
[latGrid, lonGrid, altGrid] = ndgrid(latVec, lonVec, altVec);
46+
47+
% a gridded interpolant structure that can be queried at positions
48+
F = griddedInterpolant(latGrid, lonGrid, altGrid, radGrid);
49+
50+
PathVec = transformPath(path, Dimensions);
51+
52+
% retrives the interpolated values at the points
53+
pathData = F(PathVec.lat, PathVec.lon, PathVec.height);
54+
55+
56+
57+
58+
59+
60+

‎pathGen.m

Lines changed: 28 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,28 @@
1+
function path = pathGen(cassiopeDatPath)
2+
% Loads Cassiope path data from supplied file, then interpolates the
3+
% supplied data. The path is returned as a path structure, that includes
4+
% path.lat, path.lon, and path.alt
5+
%
6+
% Inputs:
7+
% cassiopeDatPath - Path to Cassiope data
8+
9+
pathData = importdata(cassiopeDatPath, ' ',2);
10+
11+
path.lat = pathData.data(:, 12);
12+
path.lon = pathData.data(:, 13);
13+
path.alt = pathData.data(:, 14);
14+
15+
templat = [];
16+
templon = [];
17+
tempalt = [];
18+
19+
for i = 1:(length(path.lat)-2) %this needs to be improved with interpolation
20+
21+
templat = [templat, linspace(path.lat(i), path.lat(i+1), 10)];
22+
templon = [templon, linspace(path.lon(i), path.lon(i+1), 10)];
23+
tempalt = [tempalt, linspace(path.alt(i), path.alt(i+1), 10)];
24+
end
25+
26+
path.lat = templat;
27+
path.lon = templon;
28+
path.alt = tempalt;

0 commit comments

Comments
 (0)