% Neptune's atmospheric absorption:
% =================================
clear all; close all;
pos = [1,29,1280,704]; % For figure display...

data = xlsread('NeptuneAbsorption.xls');

% We only need distance (maybe pressure?) from this sheet:
data = data(:,1:2);
distance = data(:,1)';

% Now get the attenuation contributions of each gas, for 100MHz,
% 500MHz, 1GHz and 2.3 GHz.

h2oabs = xlsread('NeptuneAbsorption.xls','H2O');
h2oabs = [h2oabs(:,1),h2oabs(:,3),h2oabs(:,5),h2oabs(:,7)];  % Cleanup empty columns

nh3abs = xlsread('NeptuneAbsorption.xls','NH3');
nh3abs = [nh3abs(:,1),nh3abs(:,3),nh3abs(:,5),nh3abs(:,7)];  % Cleanup empty columns

h2sabs = xlsread('NeptuneAbsorption.xls','H2S');
h2sabs = [h2sabs(:,1),h2sabs(:,3),h2sabs(:,5),h2sabs(:,7)];  % Cleanup empty columns

ph3abs = xlsread('NeptuneAbsorption.xls','PH3');
ph3abs = [ph3abs(:,1),ph3abs(:,3),ph3abs(:,5),ph3abs(:,7)];  % Cleanup empty columns

h2abs = xlsread('NeptuneAbsorption.xls','H2');
h2abs = [h2abs(:,1),h2abs(:,3),h2abs(:,5),h2abs(:,7)];  % Cleanup empty columns

% Plot them to get an idea:
att_h = figure(1);
set(att_h,'Position',pos);
set(att_h,'Name','Attenuations for different gases in Neptune atmosphere at different frequencies');
subplot(221);
plot(distance, h2oabs(:,1)','r',distance, nh3abs(:,1)','g',distance, h2sabs(:,1)','b',distance,ph3abs(:,1)','c',distance, h2abs(:,1)','m');
title('H2O(r), NH3(g), H2S(b), PH3(c), H2(m), attentuations at 100 MHz');
xlabel('Distance (km) away from 500 bars (towards the surface)');
ylabel('Attenuation in dB/km');

subplot(222);
plot(distance, h2oabs(:,2)','r',distance, nh3abs(:,2)','g',distance, h2sabs(:,2)','b',distance,ph3abs(:,2)','c',distance, h2abs(:,2)','m');
title('H2O(r), NH3(g), H2S(b), PH3(c), H2(m), attentuations at 500 MHz');
xlabel('Distance (km) away from 500 bars (towards the surface)');
ylabel('Attenuation in dB/km');

subplot(223);
plot(distance, h2oabs(:,3)','r',distance, nh3abs(:,3)','g',distance, h2sabs(:,3)','b',distance,ph3abs(:,3)','c',distance, h2abs(:,3)','m');
title('H2O(r), NH3(g), H2S(b), PH3(c), H2(m), attentuations at 1 GHz');
xlabel('Distance (km) away from 500 bars (towards the surface)');
ylabel('Attenuation in dB/km');

subplot(224);
plot(distance, h2oabs(:,4)','r',distance, nh3abs(:,4)','g',distance, h2sabs(:,4)','b',distance,ph3abs(:,4)','c',distance, h2abs(:,4)','m');
title('H2O(r), NH3(g), H2S(b), PH3(c), H2(m), attentuations at 2.3 GHz');
xlabel('Distance (km) away from 500 bars (towards the surface)');
ylabel('Attenuation in dB/km');

pause(4);
% Calculate total attenuation for at every distance:
% ==================================================

% Cumulative attenuation, at different frequencies:
for i=1:4
toth2o = integr(distance,h2oabs(:,i)');
totnh3 = integr(distance,nh3abs(:,i)');
toth2s = integr(distance,h2sabs(:,i)');
totph3 = integr(distance,ph3abs(:,i)');
toth2 = integr(distance,h2abs(:,i)');
if(i==1)
totatt100 = toth2o + totnh3 + toth2s + totph3 + toth2;
end
if(i==2)
totatt500 = toth2o + totnh3 + toth2s + totph3 + toth2;
end
if(i==3)
totatt1000 = toth2o + totnh3 + toth2s + totph3 + toth2;
end
if(i==4)
totatt2300 = toth2o + totnh3 + toth2s + totph3 + toth2;
end
end
tot_att_h = figure;
set(tot_att_h,'Position',pos);
set(tot_att_h, 'Name', 'Total attenuation at different frequencies');
plot(distance,totatt100,'r',distance,totatt500,'g',distance,totatt1000,'b',distance,totatt2300,'c');
title('Attenuations at 100 MHz (r), 500 MHz (g), 1 GHz (b), 2.3 GHz (c)');
xlabel('Distance (km) away from 500 bars (towards the surface)');
ylabel('Attenuation in dB');

pause(4);
% ==============================================
% Now interpolate linearly for broadband signals
% ==============================================

% For a given distance, find the attenuation at some frequency f,
% knowing the attenuation at the following frequencies:

freq = [100,500,1000,2300];
step = 5;
f_interp = 100:step:2300;    % find attenuation every 10 Mhz
att_interp = zeros(1,length(f_interp));   % Array to contain the interpolated attenuation

% find attenuation at this frequency for all distances
for j=1:length(distance)
att_f = interp1(freq,[totatt100(j),totatt500(j),totatt1000(j),totatt2300(j)],f_interp,'linear');
% Store it:
att_interp(j,:) = att_f;    % First column is for 100 MHz, second for 110 MHz etc... last column is at 2300 MHz.
end

% Let's try plotting the attenuations at 150, 750, and 1500 MHz:
att_interp_h = figure;
set(att_interp_h, 'Position', pos, 'Name','Interpolated attenuations');
f = [150,750,1500];
for i=1:length(f)
column(i) = (f(i) - freq(1))/step + 1;
end
plot(distance,totatt100,'r',distance,totatt500,'g',distance,totatt1000,'b',distance,totatt2300,'c',...
    distance,att_interp(:,column(1))','m--',distance,att_interp(:,column(2))','y--',distance,att_interp(:,column(3))','k--');
title('Interpolated attenuations at 150 (m),750 (y) and 1500 (k) MHz vs. actual data');
xlabel('Distance (km) away from 500 bars (towards the surface)');
ylabel('Attenuation in dB');
pause(4);
% Plot the attenuation at the surface (farthest distance) vs. frequency
% for broadband signal attenuation estimation:
att_broad = -att_interp(end,:); % Invert sign to display as a filter gain
att_broad_h = figure;
set(att_broad_h,'Position',pos,'Name','Filter model of Neptune atmosphere');
plot(f_interp,att_broad);
title('Neptune atmosphere filter model: Total loss of a signal originating from 500 bars in the atmosphere');
xlabel('Frequency (MHz)'); ylabel('Gain (dB)');
grid;









