1 |
function [fig1]=confplot(EMEP,rangeE,radialE,spec,sysinfo) |
---|
2 |
|
---|
3 |
%set up the wavesmon data in a structure |
---|
4 |
|
---|
5 |
%what is the magnetic variation to nearest degree |
---|
6 |
magvar=10; |
---|
7 |
|
---|
8 |
%first change to m^2/Hz/deg |
---|
9 |
wmon.S=spec/(360*1000*1000); |
---|
10 |
|
---|
11 |
%what is the start angle |
---|
12 |
heading=sysinfo(18,:); |
---|
13 |
heading=heading/100; |
---|
14 |
sangle=heading+magvar; |
---|
15 |
|
---|
16 |
%what is the frequency and dir resolution for those generated in DIWASP |
---|
17 |
freqres=0.01; |
---|
18 |
freqs=[0.01:0.01:.4]; |
---|
19 |
dirres=2; |
---|
20 |
dirs=[-180:2:180]; |
---|
21 |
|
---|
22 |
%set up the directions |
---|
23 |
adj_angle=90-sangle+360+180; |
---|
24 |
wmon.dirs=[adj_angle:-4:-(356-adj_angle)]; |
---|
25 |
wmondirres=4; |
---|
26 |
wmon.xaxisdir=90; |
---|
27 |
wmon.freqs=[0.00781250:0.00781250:1]; |
---|
28 |
wmonfreqres=0.00781250; |
---|
29 |
|
---|
30 |
%calculate just the frequency energy spectrum |
---|
31 |
EMEPfreq=sum(EMEP.S')*dirres; |
---|
32 |
EMEPrangefreq=sum(rangeE.S')*dirres; |
---|
33 |
EMEPradialfreq=sum(real(radialE.S)')*dirres; |
---|
34 |
wmonfreq=sum(wmon.S')*wmondirres; |
---|
35 |
|
---|
36 |
%Compute the coefficient for the upper and lower error bounds for the power |
---|
37 |
%spectrum assuming 95% confidence. Added on 9/17/08 |
---|
38 |
degF=radialE.degF; |
---|
39 |
chiUp=chi2inv(.975,degF); |
---|
40 |
chiLow=chi2inv(.025,degF); |
---|
41 |
coeffUp=degF/chiLow; |
---|
42 |
coeffLow=degF/chiUp; |
---|
43 |
logCoeffUp=log10(coeffUp); |
---|
44 |
logCoeffLow=log10(coeffLow); |
---|
45 |
|
---|
46 |
%calculate the conf limits throughout the frequency spectrum |
---|
47 |
EMEPradialfreqUP=EMEPradialfreq*coeffUp; |
---|
48 |
EMEPradialfreqLOW=EMEPradialfreq*coeffLow; |
---|
49 |
|
---|
50 |
%same thing for the Puvw and range |
---|
51 |
degF2=25; |
---|
52 |
chiUp2=chi2inv(.975,degF2); |
---|
53 |
chiLow2=chi2inv(.025,degF2); |
---|
54 |
coeffUp2=degF2/chiLow2; |
---|
55 |
coeffLow2=degF2/chiUp2; |
---|
56 |
logCoeffUp2=log10(coeffUp2) |
---|
57 |
logCoeffLow2=log10(coeffLow2) |
---|
58 |
|
---|
59 |
%calculate the conf limits throughout the frequency spectrum |
---|
60 |
EMEPfreqUP=EMEPfreq*coeffUp2; |
---|
61 |
EMEPfreqLOW=EMEPfreq*coeffLow2; |
---|
62 |
EMEPrangefreqUP=EMEPrangefreq*coeffUp2; |
---|
63 |
EMEPrangefreqLOW=EMEPrangefreq*coeffLow2; |
---|
64 |
|
---|
65 |
scrsz = get(0,'ScreenSize'); |
---|
66 |
fig1=figure('Position',[scrsz]); |
---|
67 |
|
---|
68 |
%plot the frequency energy spectrum |
---|
69 |
h1=semilogy(freqs,EMEPfreq,'b'); |
---|
70 |
hold on |
---|
71 |
plot([0.1,0.1],[10^logCoeffLow2,10^logCoeffUp2],'b.-'); |
---|
72 |
%plot(freqs,EMEPfreqUP,'b--'); |
---|
73 |
%plot(freqs,EMEPfreqLOW,'b--'); |
---|
74 |
|
---|
75 |
h2=semilogy(freqs,EMEPrangefreq,'r'); |
---|
76 |
plot([0.11,0.11],[10^logCoeffLow2,10^logCoeffUp2],'r.-'); |
---|
77 |
%plot(freqs,EMEPrangefreqUP,'r--'); |
---|
78 |
%plot(freqs,EMEPrangefreqLOW,'r--'); |
---|
79 |
|
---|
80 |
h3=semilogy(freqs,EMEPradialfreq,'g'); |
---|
81 |
plot([0.12,0.12],[10^logCoeffLow,10^logCoeffUp],'g.-'); |
---|
82 |
%plot(freqs,EMEPradialfreqUP,'g--'); |
---|
83 |
%plot(freqs,EMEPradialfreqLOW,'g--'); |
---|
84 |
ylim([10^-2,10^1]); |
---|
85 |
%set(gca,'ytick',[-4 -3 -2 -1 0 1]); |
---|
86 |
axis(axis); |
---|
87 |
h4=semilogy(wmon.freqs,wmonfreq,'k'); |
---|
88 |
legend([h1,h2,h3,h4],'EMEP uvw','EMEP range','EMEP radial','wavesmon','location','best'); |
---|
89 |
title('directional wave spectrum integrated over direction'); |
---|
90 |
xlabel('frequency in Hz'); |
---|
91 |
ylabel('m^2 / hz'); |
---|
92 |
|
---|
93 |
|
---|
94 |
|
---|
95 |
|
---|
96 |
|
---|
97 |
|
---|
98 |
|
---|
99 |
|
---|