-
Notifications
You must be signed in to change notification settings - Fork 1
/
Copy pathspitspec.m
116 lines (116 loc) · 2.78 KB
/
spitspec.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
function h=spitspec(fig,ax,sacf,lag,dt,r0)
global bval d_parbl gating el site combine combhold maxsize_acf selax webtg
if nargin<6, r0=0; end
if ~isempty(selax) && selax.fig==fig && selax.axes==ax
selax.sacf=sacf; selax.lag=lag;
if nargin>4
selax.dt=dt; selax.r0=r0;
end
end
sgating=gating;
[nlags,ngates]=size(sacf);
if nlags>maxsize_acf(1)
nlags=maxsize_acf(1);
sacf=sacf(1:nlags,:); lag=lag(1:nlags);
end
if ngates/sgating>maxsize_acf(2)
sgating=ceil(ngates/maxsize_acf(2));
end
if sgating>1
ngates=floor(ngates/sgating);
sacf=reshape(sum(reshape(sacf(:,1:(ngates*sgating))',sgating,ngates,nlags)),ngates,nlags)';
end
waterlim=[4 6];
if bval(11)==1 && ngates>waterlim(2) && r0>0
[h,pars]=qfit(fig,ax,sacf,lag,dt*sgating,r0+(sgating-1)*dt/2);
if ~isempty(selax) && selax.fig==fig && selax.axes==ax
selax.pars=pars;
end
return
elseif bval(11)==1 && nargin>4 && (site==5 || site==6)
[s,w]=max(sacf(1,:));
par=qfit(fig,ax,sacf(:,w),lag,0,1e-2);
end
if ~isempty(webtg) && ~isempty(selax)
h=[];
return
end
fsc=1000;
if bval(7)~=1
[s,w]=acf2spec(sacf,lag);
s=s*1000;
if max(w)<fsc, fsc=1; end
w=w/fsc;
elseif ngates>1
s=real(sacf);
if max(lag)*fsc>1, fsc=1; end
w=lag*fsc*1e3;
else
s=[real(sacf) imag(sacf)];
if max(lag)*fsc>1, fsc=1; end
w=lag*fsc*1e3;
end
if isempty(combhold) && ngates<waterlim(1)
updateplot(fig,ax,w,s);
h=get(ax,'ylabel');
else
if bval(7)~=1, s(find(s<0))=0; end
r=(r0+(0:(ngates-1))*dt*sgating)*1e6;
if r0>0
r=r*.15;
re=6370; r=r/re; r=re*sqrt(1+r.*(r+2*sin(el/57.2957795)))-re;
end
if ~isempty(combhold)
if length(r)>1
dr=diff(r(1:2));
else
dr=diff(combhold.r(1:2));
end
%s=interp2(r,w',s,r,combhold.w');
s=interp1(w',s,combhold.w');
n=round((r(1)-combhold.r(end))/dr)-1;
r=[combhold.r (1:n)*dr+combhold.r(end) r];
s=[combhold.s NaN*ones(length(combhold.w),n) s];
ngates=length(r); combhold=[];
end
if ~isempty(combine) && any(combine==ax)
combhold.s=s; combhold.r=r; combhold.w=w;
d=find(combine==ax); combine(d(1))=1;
h=[]; return
end
if ngates<waterlim(2) || any(get(ax,'view')-[0 90])
updatewater(fig,ax,w,r,s')
else
updateimage(fig,ax,w,r,s')
end
h=get(ax,'ylabel');
if r0>0
set(h,'string','Altitude (km)')
else
set(h,'string','\mus')
end
h=get(ax,'zlabel');
end
if bval(7)~=1
set(h,'string','Power (K/kHz)')
if fsc==1000
set(get(ax,'xlabel'),'string','Frequency (kHz)')
else
set(get(ax,'xlabel'),'string','Frequency (Hz)')
end
set(ax,'xgrid','on','ygrid','on')
else
set(h,'string','Power (K)')
if fsc==1000
set(get(ax,'xlabel'),'string','Lag (\mus)')
else
set(get(ax,'xlabel'),'string','Lag (ms)')
end
set(ax,'xgrid','off','ygrid','on')
end
h=get(ax,'title');
set(ax,'visible','on')
if bval(11)==1 && nargin>4 && (site==5 || site==6)
set(h,'string',sprintf('%.0fK %.0fK %.0fms^{-1}',par))
h=get(ax,'zlabel');
end