% рисование часовых данных обсерваторий ИКИР

load YSS_1957-1988_hourly_WDC_20180614 Tyss Xyss Yyss Hyss Dyss Zyss
load KHB_2012-2016_hourly_WDC_20180614 Tkhb Xkhb Ykhb Hkhb Dkhb Zkhb
load PET_1969-2016_hourly_WDC_20180614 Tpet Xpet Ypet Hpet Dpet Zpet
load MGD_1966-2017_hourly_WDC_20180614 Tmgd Xmgd Ymgd Hmgd Dmgd Zmgd

load CPS_hourly_2015-2016 T H D Z F
Tcps=T; Hcps=H; Dcps=D; Zcps=Z; Fcps=F;

k=Tmgd>datenum(2012,1,1); Hmgd(k)=total_field(Xmgd(k),Ymgd(k));
k=Tpet>datenum(2003,1,1); Hpet(k)=total_field(Xpet(k),Ypet(k));
Hkhb=total_field(Xkhb,Ykhb);

% Intermagnet
load pet_hourly_2001-2023.mat
Hhpet=total_field(Xhpet,Yhpet);

figure(1), HH=20000;
plot(Tpet,Hpet-HH,'b',Tmgd,Hmgd-HH+3500,'r',...
     Tkhb,Hkhb-HH-2000,'g',Tyss,Hyss-HH-2500,'k',...
     Tcps,Hcps-HH+11000,'y',...
     Thpet,Hhpet-HH,'b')
datetick('x','yyyy'), ylabel('H-20000, nT')
legend('PET','MGD+3500','KHB-2000','YSS-2500','CPS-11000','Location','SouthWest')
title('Magnetic data (hourly values) from WDC (Edinburg + Moscow)')
% маркировка кривых
text(datenum(1989,1,1),2700,'YSS','Color','k','FontWeight','Bold')
text(datenum(1967,1,1),2250,'PET','Color','b','FontWeight','Bold')
text(datenum(1965,1,1),1700,'MGD','Color','r','FontWeight','Bold')
text(datenum(2011,1,1),2150,'KHB','Color','g','FontWeight','Bold')
text(datenum(2011,1,1),3000,'CPS','Color','y','FontWeight','Bold')
% % маркировка данных
% % PET (photo & digit)
% T1=datenum(1998,1,1); T2=datenum(2005,12,31); H1=1600; H2=1750;
% patch([T1 T1 T2 T2],[H1 H2 H2 H1],[0.5 1 1],'FaceAlpha',0.3)
% T1=datenum(2002,1,1); T2=datenum(2006,12,31); H1=1550; H2=1800;
% patch([T1 T1 T2 T2],[H1 H2 H2 H1],[1 0.2 1],'FaceAlpha',0.3)
% % MGD (photo & digit)
% T1=datenum(1998,1,1); T2=datenum(2005,12,31); H1=750; H2=1250;
% patch([T1 T1 T2 T2],[H1 H2 H2 H1],[0.5 1 1],'FaceAlpha',0.3)
% text(datenum(2000,1,1),1000,'Analogue')
% T1=datenum(2005,5,1); T2=datenum(2011,12,31); H1=600; H2=1100;
% patch([T1 T1 T2 T2],[H1 H2 H2 H1],[1 0.2 1],'FaceAlpha',0.3)
% text(datenum(2007,6,1),900,'Digital')
% % KHB (digit)
% T1=datenum(2005,1,1); T2=datenum(2011,12,31); H1=1900; H2=2100;
% patch([T1 T1 T2 T2],[H1 H2 H2 H1],[1 0.2 1],'FaceAlpha',0.3)
% % YSS (digit)
% T1=datenum(2005,11,1); T2=datenum(2007,10,31); H1=2450; H2=2650;
% patch([T1 T1 T2 T2],[H1 H2 H2 H1],[1 0.2 1],'FaceAlpha',0.3)
% 
% 
% 
% % figure(3)
% % plot(Tpet,Zpet,Tmgd,Zmgd,Tkhb,Zkhb,Tyss,Zyss)
% % datetick('x','yyyy')
% % 
% % figure(1)
