Electric and Magnetic Fields

%% --- ECE 390 - Final Project (Fall 2020) --- %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% close all clear all clc %--- Refer to Dobson 1985 for details about the model ---- %% --- 1. Enter Frequency, Moisutre Content and Temperature ---- % Model is valid between 1.4GHz and 18 GHz freq = STARTFREQ:FREQSTEP:STOPFREQ; % Frequency range of interest in Hz. mv = 0.0:STEPSIZE:0.8;% Volumetric moisture content (between 0 and 1). T = TEMPERATURE; %-- Temperature of the soil in degree Celsius %% --- 2. Other Soil Parameters and Constants ---- %%%%%%%%%%% DO NOT CHANGE %%%%%%%%%%%%% S = 50.73 ; %-- Percentage of sand in the soil C = 9.66;% percentage of clay in the soil rho_s = 2.63;% specific density of the soil Porosity = 0.4; rho_b = (1 - Porosity)*rho_s; % bulk density alpha = 0.65; eps0 = 8.854e-12; %permittivity of free space mu0 = pi*4e-7; %permeability of free space %% -- 3. Calculations -- %%%% 3a Calculate the path loss from Friis equation %%%%%%% %%%% 3b Now calculate the Dobson permittivity of soil %%%%% for ffi = 1:length(freq) %Looping frequency f = freq(ffi); eps_s = (1.01+0.44*rho_s).^2 - 0.062; beta_eps_prime = (127.48 - 0.519*S - 0.152*C)/100; beta_eps_double_prime =(133.797 - 0.603*S - 0.166*C)/100; beta = beta_eps_prime+1j*beta_eps_double_prime; sigma_eff = -1.645 + 1.939*rho_b - 0.02013*S + 0.01594*C; Ave_Sal = 0.64 *sigma_eff*10; %-- Average salinity of water from eq 19 in Dobson in parts per thousand N = Ave_Sal.*(1.707e-2 + 1.205e-5*Ave_Sal + 4.058e-9*Ave_Sal.^2); %-- Normality of water calculated from salinity Eq(8) in Stogryn a_N = 1 - 0.2551*N + 5.151e-2 *N.^2 - 6.889e-3*N.^3; %-- Eq (4) in Stogryn ew_inf = 4.9; for ii = 1:length(mv) %Looping volumetric moisture b_N = 0.1463e-2*N*T + 1 - 0.04896*N - 0.02967*N.^2 + 5.644e-3*N.^3; %-- Eq (5) in Stogryn e0_purewater = 87.74 - 0.40008*T + 9.398e-4*T^2 + 1.410e-6*T^3; %-- Eq 6 in Stogryn ew0 = e0_purewater*a_N; tau0_purewater = 1.1109e-10 - 3.824e-12*T + 6.938e-14*T^2 - 5.096e-16*T^3; %-- Eq 7 in Stogryn tau_w = (tau0_purewater*b_N)/(2*pi); efw = ew_inf + (ew0 - ew_inf)/(1+1j*2*pi*f*tau_w) - 1j*sigma_eff/(2*pi*f*eps0)*(rho_s - rho_b)/(rho_s*mv(ii)); efw_real = ew_inf + (ew0 - ew_inf)/(1+(2*pi*f*tau_w.^2)); efw_imag_1 = 2*pi*f*tau_w*(ew0 - ew_inf)/(1+(2*pi*f*tau_w.^2)); efw_imag_2 = sigma_eff/(2*pi*eps0*f)*(rho_s-rho_b)/(rho_s*mv(ii)); efw_imag = efw_imag_1+efw_imag_2; em_real_Dob(ffi,ii) = (1 + rho_b/rho_s *(eps_s^alpha - 1) + mv(ii)^beta_eps_prime*(efw_real)^alpha - mv(ii))^(1/alpha); em_imag_Dob(ffi,ii) = (mv(ii)^beta_eps_double_prime*(efw_imag)^alpha)^(1/alpha); end end %% --- 4. Plot Pathloss versus Distance ----- %%%%%%%%%%% Plot the pathloss db calculated in 3a %% --- 5. Plot Permittivity vs Frequency ----- % The data is compiled as epsilon_r(frequency, mv) figure %--- for tt = 1:length(mv) subplot(2,1,1) hold on plot(freq/1e9,em_real_Dob(:,tt)) subplot(2,1,2) hold on plot(freq/1e9,em_imag_Dob(:,tt)) end % Adding legends and axis to the plots subplot(2,1,1) legend(strcat('m_v=',num2str(mv'))) ylabel('Real(\epsilon_r)') xlabel('Frequency (GHz)') title('Permittivity vs Frequency') subplot(2,1,2) legend(strcat('m_v=',num2str(mv'))) xlabel('Frequency (GHz)') ylabel('Imag(\epsilon_r)') %% --- 6. Plot Permittivity vs Moisture Content ----- % The data is compiled as epsilon_r(frequency, mv) figure %--- %for tt = 1:length(freq) % can remove use this for loop (and first plot %statements below) to plot for all frequencies in the freq vector. subplot(2,1,1) hold on %plot(mv,em_real_Dob(tt,:)) plot(mv,em_real_Dob(length(freq),:)) % plots at a single frequency (last value in freq vector) subplot(2,1,2) hold on %plot(mv,em_imag_Dob(tt,:)) plot(mv,em_imag_Dob(length(freq),:)) % plots at a single frequency (last value in freq vector) %end % Adding legends and axis to the plots subplot(2,1,1) %legend(strcat('f(GHz)=',num2str(freq'/1e9))) % plots legend if all %frequencies are used for the plot. ylabel('Real(\epsilon_r)') xlabel('m_v') title(strcat('Permittivity vs Moisture Content-',num2str(freq(length(freq))/1e9),' GHz')) subplot(2,1,2) %legend(strcat('f(GHz)=',num2str(freq'/1e9))) xlabel('m_v') ylabel('Imag(\epsilon_r)') %% --- 7. Calculate and Plot the Reflection Coefficient, Return Loss and Absorbed Power %%%%%%%%%%%%% Calculate the reflection coefficient as a function of moisture content here %%%%%%%%%%%%% Then plot the reflection coefficient, return loss and absorbed power as a function of moisture content. %% --- 8. Calculate and Plot Received Power from Metal Plate with sigma=1 %%%%%%%%%%%%% Calculate the power and then display its value. %% --- 9. Read in permittivity data and create image %%%%%%%%%%%%% Read the permittivity matrices in the text file. Calculate the recieved power in dB %%%%%%%%% %%%%%%%%%%%%% Normalize the recieved power (in db units) with respect to the maximum and minimum received power and then multiply it with 255 (The reason to multiply the result by 255 is to show it in the figure) %%%%%%%%%%%%% Display the figure.