How to align my plots correctly with the x-axis?

13 Ansichten (letzte 30 Tage)
Alex
Alex am 15 Feb. 2020
Kommentiert: Alex am 18 Feb. 2020
First, I apologize for the title of the question. Essentially, when I plot two series individually (see P7.png and P25.png) and compare them with the dashed lines at 4, 7, 9.75 and 11 degrees (where the peaks should occur) they do pretty well. As soon as I try to plot them on the same graph they seem to be shifted and no longer align (see P7and25.png). This gets worse the more dataset that I have. I have included the data as test.mat. Thank you in advance for your help.
P7.png:
P25.png:
P7and25.png:
load(test.mat)
h = cell(1,length(d));
cmap = jet(m);
figure
%xlim([-5 20])
hold all
for k=1:m
%Apply 50 degree elevation cutoff
idx2 = find(All_oelev(:,k) > 50 & All_oelev(:,k) < 140); % indices to plot of All_oelev
plot(xi2(k,idx2),VTEC2(idx2,k), 'Color', cmap(k, :))%, '--g','linewidth',2)% Fig5
sz=size(idx2);
idx2 = NaN(sz);
%h{k} =sprintf('aircraft at %3.4f degrees latitude', parnam{1,k}.data(1,1)) ;
end
plot([4 4], [0 max(max(AllTEC))],'--k','LineWidth',2);
plot([7 7], [0 max(max(AllTEC))],'--k','LineWidth',2);
plot([9.75 9.75], [0 max(max(AllTEC))],'--k','LineWidth',2);
plot([11 11], [0 max(max(AllTEC))],'--k','LineWidth',2);
%ylim([-0.4 0.8])
xlabel('Latitude (Degree)','fontsize',24, 'fontweight','bold');
ylabel('Total Electron Content (TECu = 10^{16} electrons/m^2)','fontsize',24, 'fontweight','bold')
%legend(h)
hold off
Best regards,
Alex
  2 Kommentare
Jakob B. Nielsen
Jakob B. Nielsen am 17 Feb. 2020
How do you generate your first two (individual) plots? Looking at your data, the peak that you want to have at 11 degrees is at 11,36 so I don't understand how you manage to make them align in your first plot.
Alex
Alex am 18 Feb. 2020
Here is the data from the first two datasets (P7_test.mat and P25_test.mat).
Here is the complete script I am using. The first part is the engine used to open the raw datasets. It should be fine, It just parses the different datasets text files into parameters such as the quantitudy 'AllTEC', elevation angles 'Alloelev' and others that are not used.
clear all
close all
clc
help VTEC_mapping_AIS
Re =6371; %Earth radius (km)
h = 350;% (250 km in Sharma example) height (km) of airglow emission (or thin shell)
% the data
proot='Plane*';% <- your FOLDER template
froot='par_o.par'; % <- your FILENAME template
hroot='Data.txt'; % <- your FILENAME template
% the engine
d=dir(proot);
d={d([d.isdir]).name}.';
fnam=cellfun(@(x) sprintf('%s%s%s',x,filesep,froot),d,'uni',false);%all data looks good (text and col headers not so much)
gnam=cellfun(@(x) sprintf('%s%s%s',x,filesep,froot),d,'uni',false);
hnam=cellfun(@(x) sprintf('%s%s%s',x,filesep,hroot),d,'uni',false);
%%%%variable for number of traces
rays=importdata(hnam{1});%rays=4; use fixed value for FFTData or rewrite similar script
rays=rays(1,1);
% - processing each file
m=length(d);
for i=1:m
cfile=fnam{i};
gfile=gnam{i};
disp(sprintf('working on file %s',cfile));
disp(sprintf('working on file %s',gfile));%should be redundant
% code to read contents of CFILE
fnam{i}=importdata(cfile, ' ', 32); %delimeter and headerlines using import wizard
% to format unknown file(s) *.par
gnam{i}=importdata(cfile, ' ', 57);%171
p_line=rays+53;
parnam{i} =importdata(cfile, ' ', p_line);%p_line should give 113 for 60 rays, 171 for 118 rays
%extract starting latitude, longitude, and altitude to corresponding
%vectors
Latrayorigin(1:rays,i)=parnam{1,i}.data(1,1);%WRONG
Lonrayorigin(1:rays,i)=parnam{1,i}.data(1,2);
Altrayorigin(1:rays,i)=parnam{1,i}.data(1,3);
Latrayend(:,i)=fnam{i,1}.data(:,3);
Lonrayend(:,i)=fnam{i,1}.data(:,2);
Altrayend(:,i)=fnam{i,1}.data(:,4);
%rawTEC{i}=fnam{i}(1:12,:);
end
Alto=0.001*Altrayorigin;%convert origin to km to match endpoint
O(1,:)=reshape(Latrayorigin.',1,numel(Latrayorigin));%x=lat
O(2,:)=reshape(Alto.',1,numel(Alto));%y=alt
O(3,:)=reshape(Lonrayorigin.',1,numel(Lonrayorigin));%z=lon
%get direction unit vectors
numx=Latrayend-Latrayorigin;
numy=Altrayend-Alto;
numz=Lonrayend-Lonrayorigin;
numex=reshape(numx.',1,numel(numx));
numey=reshape(numy.',1,numel(numy));
numez=reshape(numz.',1,numel(numz));
denominator=sqrt(numex.^2 + numey.^2 + numez.^2);%magnitude of direction vector
D(1,:)=numex./denominator;%x=lat
D(2,:)=numey./denominator;%y=alt
D(3,:)=numez./denominator;%z=lon
%direction= %direction unit vector
for i=1:m
hfile=hnam{i};
disp(sprintf('working on file %s',hfile));
% your code to read contents of CFILE
hnam{i}=importdata(hfile);
%satellite latitudes (convergence latittude)
rawsatlat{i}=hnam{i}(2:(rays+1));
ophase{i}=hnam{i}((rays+2):((2*rays)+1));
xphase{i}=hnam{i}(((2*rays)+2):((3*rays)+1));
oaspect{i}=hnam{i}(((3*rays)+2):((4*rays)+1));
xaspect{i}=hnam{i}(((4*rays)+2):((5*rays)+1));
altitude{i}=hnam{i}(((5*rays)+2):((6*rays)+1));
oelev{i}=hnam{i}(((6*rays)+2):((7*rays)+1));
xelev{i}=hnam{i}(((7*rays)+2):((8*rays)+1));
opath{i}=hnam{i}(((8*rays)+2):((9*rays)+1));
xpath{i}=hnam{i}(((9*rays)+2):((10*rays)+1));
otime{i}=hnam{i}(((10*rays)+2):((11*rays)+1));
xtime{i}=hnam{i}(((11*rays)+2):((12*rays)+1));
timediff{i}=hnam{i}(((12*rays)+2):((13*rays)+1));
FRone{i}=hnam{i}(((13*rays)+2):((14*rays)+1));
FRtwo{i}=hnam{i}(((14*rays)+2):((15*rays)+1));
rotrateone{i}=hnam{i}(((15*rays)+2):((16*rays)+1));
rotratetwo{i}=hnam{i}(((16*rays)+2):((17*rays)+1));
maxedens{i}=hnam{i}(((17*rays)+2):((18*rays)+1));
minedens{i}=hnam{i}(((18*rays)+2):((19*rays)+1));
TEC{i}=hnam{i}(((19*rays)+2):((20*rays)+1));
opower{i}=hnam{i}(((20*rays)+2):((21*rays)+1));
xpower{i}=hnam{i}(((21*rays)+2):((22*rays)+1));
orientation{i}=hnam{i}(((22*rays)+2):((23*rays)+1));
ellipticity{i}=hnam{i}(((23*rays)+2):((24*rays)+1));
%Combine TEC data from all iterations into a single array
rowrawTEC{i}=reshape(TEC{i}.',rays,1);
AllTEC(:,i)=rowrawTEC{i};
%Combine FR data from all iterations into a single array
rowrawFR{i}=reshape(FRone{i}.',rays,1);
AllFR(:,i)=rowrawFR{i};
% refine calculation by considering the elvation angle
row_oelev{i}=reshape(oelev{i}.',rays,1);
All_oelev(:,i)=abs(row_oelev{i});
% rawAlt{i}=hnam{i}(7,:);
% %Combine all Altitude data, from all iterations into a single Altitude Data
% %matrix
% rowrawAlt{i}=reshape(rawAlt{i}.',rays,1);
% AllAlt(:,i)=rowrawAlt{i};
end
TECfinal(1,:)=reshape(AllTEC.',1,numel(AllTEC));
FRfinal(1,:)=reshape(AllFR.',1,numel(AllFR));
latrayendfinal(1,:)=reshape(Latrayend.',1,numel(Latrayend));
% Find and remove any TEC ==0
for kk =1:1:m
Tidx = find(AllTEC(:,kk)==0); % indices
AllTEC(Tidx,kk)=NaN;
AllFR(Tidx,kk)=NaN;
All_oelev(Tidx,kk)=NaN;
Latrayend(Tidx,kk)=NaN;
Latrayorigin(Tidx,kk)=NaN;
end
%interpolate for NaN values
xi = linspace(min(Latrayend(:,1)), max(Latrayend(:,1)), numel(Latrayend(:,1))); % Interpolating Vector
%2nd method of calculating VTEC
zp =sind(90-All_oelev).*(Re/(Re+h));
VTEC2 =AllTEC.*sqrt(1-sind(zp).^2);
% FIGURE 3
% plot VTEC as a function of satlat with lines of different color and style
figure
hold all
for k=1:m
VTEC2(:,k) = interp1(Latrayend(~isnan(VTEC2(:,k)),k), VTEC2(~isnan(VTEC2(:,k)),k), xi, 'linear');
line(xi,VTEC2(:,k), 'Color', cmap(k, :), 'LineStyle', '--')% USE FOR LARGE NUMBER OF
end
ylabel('Total Electron Content (TECu = 10^{16} electrons/m^2)','fontsize',24, 'fontweight','bold')
xlabel('Satellite Latitude (Degree)','fontsize',24, 'fontweight','bold');
title('FIG 4- VTEC (first and second method) as a function of satellite latitude')
legend('VTEC method 1', 'VTEC method 2')
hold off
%Mapping the IPP subionopheric latitude (from VTEC_mapping_example.m)
%ground receiver located at:
gnd_lon = Lonrayorigin;%74.2;% 74.2 E longitude
L = Latrayorigin;%16.8; % gnd_lat 16.8 N latitude
% satellite transmitter located at:
sat_lon = Lonrayend;%71.2;%71.2 E longitude
% %TEST using latitude instead of longitude
for k=1:m
for ijk=1:length(VTEC)
G(k,ijk) = gnd_lon(ijk,k) - sat_lon(ijk,k); %(3 degrees) difference between satellite longitude and ground receiver longitude
Az(k,ijk) = 180 + atand(tand(G(k,ijk)./tand(L(ijk,k))));%(189.8 degrees) azimuth angle
El(k,ijk) = All_oelev(ijk,k);%= atand((cosd(G(ijk)).*cosd(L(ijk)) - 0.1512)./(sqrt(1 - cosd(G(ijk)).^2.*cosd(L(ijk)).^2))) ;% (69.97 degrees) elevation angle
Chi_pp(k,ijk) = 180/2 - El(k,ijk) - asind((Re/(Re + h)).*cosd(El(k,ijk))); % (0.79 degrees) Earth centred angle between ground and IPP
lat_pp(k,ijk) = asind(sind(L(ijk,k)).*cosd(Chi_pp(k,ijk)) + cosd(L(ijk,k)).*sind(Chi_pp(ijk)).*cosd(Az(ijk)));%(16.02 degrees) latitude of IPP
lon_pp(k,ijk) = gnd_lon(ijk,k) + asind((sind(Chi_pp(k,ijk)).*sind(Az(k,ijk)))./cosd(lat_pp(k,ijk)));% (74.06 degrees) longitude of IPP
end
%interpolate for NaN values in lat_pp
xi2(k,:) = linspace(min(lat_pp(k,:)), max(lat_pp(k,:)), numel(lat_pp(k,:))); % Interpolating Vector
end
Previous code to plot the data goes here.

Melden Sie sich an, um zu kommentieren.

Antworten (0)

Kategorien

Mehr zu Earth, Ocean, and Atmospheric Sciences finden Sie in Help Center und File Exchange

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by