20 #include "TPaletteAxis.h" 24 isLogz(false),
title(
"") {
31 if(ztype.CompareTo(
"amplitude")==0) iztype=1;
32 if(ztype.CompareTo(
"energy")==0) iztype=2;
33 if(iztype==0) {cout <<
"CWB::STFT::STFT: Error wrong ztype (amplitude/energy)" << endl;
exit(1);}
35 if (nfft<=0) {cout <<
"CWB::STFT::STFT: Error nfft must be positive" << endl;
exit(1);}
36 if (noverlap<=0||noverlap>=nfft) {cout <<
"CWB::STFT::STFT: Error noverlap must be > 0 & < nfft" << endl;
exit(1);}
37 if (x.
rate()<=0) {cout <<
"CWB::STFT::STFT: Error sample rate must be > 0" << endl;
exit(1);}
38 if (x.
start()<0) {cout <<
"CWB::STFT::STFT: Error start time must be positive" << endl;
exit(1);}
39 if (
int(x.
size())<nfft) {cout <<
"CWB::STFT::STFT: Error size must be > nfft " << endl;
exit(1);}
46 CWB::Window wnd(const_cast<char*>(fwindow.Data()),nfft,fparam);
50 double df=(double)x.
rate()/(double)(nfft);
54 int N=
int(nfft/nshift);
66 double norm=sqrt(x.
rate());
67 for (
int n=0;
n<N*loops;
n++) {
69 if(shift+nfft>
int(x.
size()))
break;
73 double time = (shift+nfft/2)*dt;
75 double psd=sqrt(pow(y.data[
i],2)+pow(y.data[
i+1],2))*sqrt(1/df);
78 if(iztype==2) psd*=
psd;
79 h2->SetBinContent(
int(nT*time/(Tmax-x.
start())),
int(nF*frequency/Fmax),
psd);
88 if(
h2!=NULL)
delete h2;
95 int dpaletteId, Option_t* goption) {
101 canvas =
new TCanvas(
name,
"LVC experiment", 300,40, 800, 600);
103 canvas->ToggleEventStatus();
107 canvas->SetFillColor(kWhite);
108 canvas->SetRightMargin(0.10);
109 canvas->SetLeftMargin(0.10);
110 canvas->SetBottomMargin(0.13);
116 gStyle->SetFrameBorderMode(0);
119 gStyle->SetTitleH(0.050);
120 gStyle->SetTitleW(0.95);
121 gStyle->SetTitleY(0.98);
122 gStyle->SetTitleFont(12,
"D");
123 gStyle->SetTitleColor(kBlue,
"D");
124 gStyle->SetTextFont(12);
125 gStyle->SetTitleFillColor(kWhite);
126 gStyle->SetLineColor(kWhite);
127 gStyle->SetNumberContours(256);
130 gStyle->SetCanvasColor(kWhite);
131 gStyle->SetStatBorderSize(1);
137 gStyle->SetPalette(1,0);
143 gStyle->SetPalette(1,0);
150 h2->SetStats(kFALSE);
151 h2->SetTitleFont(12);
153 h2->SetFillColor(kWhite);
155 h2->GetXaxis()->SetNdivisions(506);
156 h2->GetXaxis()->SetLabelFont(42);
157 h2->GetXaxis()->SetLabelOffset(0.014);
158 h2->GetXaxis()->SetTitleOffset(1.4);
159 h2->GetYaxis()->SetTitleOffset(1.2);
160 h2->GetYaxis()->SetNdivisions(506);
161 h2->GetYaxis()->SetLabelFont(42);
162 h2->GetYaxis()->SetLabelOffset(0.01);
163 h2->GetZaxis()->SetLabelFont(42);
164 h2->GetZaxis()->SetNoExponent(
false);
165 h2->GetZaxis()->SetNdivisions(506);
167 h2->GetXaxis()->SetTitleFont(42);
168 h2->GetXaxis()->SetTitle(
"Time (sec)");
169 h2->GetXaxis()->CenterTitle(
true);
171 h2->GetYaxis()->SetTitleFont(42);
172 h2->GetYaxis()->SetTitle(
"Frequency (Hz)");
173 h2->GetYaxis()->CenterTitle(
true);
178 h2->GetZaxis()->SetTitleOffset(0.6);
179 h2->GetZaxis()->SetTitleFont(42);
181 h2->GetZaxis()->CenterTitle(
true);
183 h2->GetXaxis()->SetLabelSize(0.03);
184 h2->GetYaxis()->SetLabelSize(0.03);
185 h2->GetZaxis()->SetLabelSize(0.03);
187 if(
title.Sizeof()==1) {
189 sprintf(stitle,
"Spectrogram (Normalized tile %s)",
ztype.Data());
190 h2->SetTitle(stitle);
192 double dt=
h2->GetXaxis()->GetBinWidth(0);
193 double df=
h2->GetYaxis()->GetBinWidth(0);
196 int nt2=
h2->GetNbinsX();
198 int nf2=
h2->GetNbinsY();
200 if(t2>t1) nt1=
int((t1-
h2->GetXaxis()->GetXmin())/dt);
201 if(t2>t1) nt2=
int((t2-
h2->GetXaxis()->GetXmin())/dt);
202 if(f2>f1) nf1=
int(f1/df);
203 if(f2>f1) nf2=
int(f2/df);
206 for (
int i=nt1;
i<nt2;
i++) {
207 for (
int j=nf1;
j<nf2;
j++) {
208 double binc=
h2->GetBinContent(
i,
j);
209 if(binc>h2max) h2max=binc;
213 h2->GetZaxis()->SetRangeUser(0,h2max);
215 if(t2>t1)
h2->GetXaxis()->SetRangeUser(t1,t2);
216 if(f2>f1)
h2->GetYaxis()->SetRangeUser(f1,f2);
217 if(z2>z1)
h2->GetZaxis()->SetRangeUser(z1,z2);
223 TPaletteAxis *
palette = (TPaletteAxis*)
h2->GetListOfFunctions()->FindObject(
"palette");
224 palette->SetX1NDC(0.91);
225 palette->SetX2NDC(0.933);
226 palette->SetTitleOffset(0.92);
227 palette->GetAxis()->SetTickSize(0.01);
238 const Int_t
NRGBs = 5;
239 const Int_t
NCont = 255;
241 if (
fabs(paletteId)==1) {
242 Double_t
stops[
NRGBs] = { 0.00, 0.34, 0.61, 0.84, 1.00 };
243 Double_t
red[
NRGBs] = { 0.00, 0.00, 0.87, 1.00, 0.51 };
244 Double_t
green[
NRGBs] = { 0.00, 0.81, 1.00, 0.20, 0.00 };
245 Double_t
blue[
NRGBs] = { 0.51, 1.00, 0.12, 0.00, 0.00 };
247 TColor::CreateGradientColorTable(NRGBs, stops, blue, green, red, NCont);
249 TColor::CreateGradientColorTable(NRGBs, stops, red, green, blue, NCont);
252 if (
fabs(paletteId)==2) {
253 Double_t
stops[
NRGBs] = { 0.00, 0.34, 0.61, 0.84, 1.00 };
254 Double_t
red[
NRGBs] = { 0.00, 0.00, 0.87, 1.00, 0.51 };
256 Double_t
green[
NRGBs] = { 0.00, 0.81, 1.00, 0.20, 0.00 };
258 Double_t
blue[
NRGBs] = { 1.00, 1.00, 0.00, 0.00, 0.00 };
260 TColor::CreateGradientColorTable(NRGBs, stops, blue, green, red, NCont);
262 TColor::CreateGradientColorTable(NRGBs, stops, red, green, blue, NCont);
265 if (
fabs(paletteId)==3) {
266 Double_t
stops[
NRGBs] = { 0.00, 0.34, 0.61, 0.84, 1.00 };
267 Double_t
red[
NRGBs] = { 0.00, 0.09, 0.18, 0.09, 0.00 };
268 Double_t
green[
NRGBs] = { 0.01, 0.02, 0.39, 0.68, 0.97 };
269 Double_t
blue[
NRGBs] = { 0.17, 0.39, 0.62, 0.79, 0.97 };
271 TColor::CreateGradientColorTable(NRGBs, stops, blue, green, red, NCont);
273 TColor::CreateGradientColorTable(NRGBs, stops, red, green, blue, NCont);
276 if (
fabs(paletteId)==4) {
277 Double_t
stops[
NRGBs] = { 0.00, 0.50, 0.75, 0.875, 1.00 };
278 Double_t
red[
NRGBs] = { 1.00, 1.00, 1.00, 1.00, 1.00 };
279 Double_t
green[
NRGBs] = { 1.00, 0.75, 0.50, 0.25, 0.00 };
280 Double_t
blue[
NRGBs] = { 0.00, 0.00, 0.00, 0.00, 0.00 };
282 TColor::CreateGradientColorTable(NRGBs, stops, blue, green, red, NCont);
284 TColor::CreateGradientColorTable(NRGBs, stops, red, green, blue, NCont);
287 if (
fabs(paletteId)==5) {
288 Double_t
stops[
NRGBs] = { 0.00, 0.34, 0.61, 0.84, 1.00 };
289 Double_t
red[
NRGBs] = { 1.00, 0.84, 0.61, 0.34, 0.00 };
290 Double_t
green[
NRGBs] = { 1.00, 0.84, 0.61, 0.34, 0.00 };
291 Double_t
blue[
NRGBs] = { 1.00, 0.84, 0.61, 0.34, 0.00 };
293 TColor::CreateGradientColorTable(NRGBs, stops, blue, green, red, NCont);
295 TColor::CreateGradientColorTable(NRGBs, stops, red, green, blue, NCont);
298 if (
fabs(paletteId)==57) {
299 Double_t
stops[9] = { 0., 0.125, 0.25, 0.375, 0.5, 0.625, 0.75, 0.875, 1.0};
300 Double_t
red[9] = { 0.2082, 0.0592, 0.0780, 0.0232, 0.1802, 0.5301, 0.8186, 0.9956, 0.9764};
301 Double_t
green[9] = { 0.1664, 0.3599, 0.5041, 0.6419, 0.7178, 0.7492, 0.7328, 0.7862, 0.9832};
302 Double_t
blue[9] = { 0.5293, 0.8684, 0.8385, 0.7914, 0.6425, 0.4662, 0.3499, 0.1968, 0.0539};
304 TColor::CreateGradientColorTable(9, stops, blue, green, red, NCont);
306 TColor::CreateGradientColorTable(9, stops, red, green, blue, NCont);
309 gStyle->SetNumberContours(NCont);
virtual void rate(double r)
void Print(TString pname)
cout<< endl;cout<< "ts size = "<< ts.size()<< " ts rate = "<< ts.rate()<< endl;tf.Forward(ts, wdm);int levels=tf.getLevel();cout<< "tf size = "<< tf.size()<< endl;double dF=tf.resolution();double dT=1./(2 *dF);cout<< "rate(hz) : "<< RATE<< "\ layers : "<< nLAYERS<< "\ dF(hz) : "<< dF<< "\ dT(ms) : "<< dT *1000.<< endl;int itime=TIME_PIXEL_INDEX;int ifreq=FREQ_PIXEL_INDEX;int index=(levels+1) *itime+ifreq;double time=itime *dT;double freq=(ifreq >0) ? ifreq *dF :dF/4;cout<< endl;cout<< "PIXEL TIME = "<< time<< " sec "<< endl;cout<< "PIXEL FREQ = "<< freq<< " Hz "<< endl;cout<< endl;wavearray< double > x
wavearray< double > psd(33)
virtual void start(double s)
double GetValue(unsigned i)
void Draw(double t1=0.0, double t2=0.0, double f1=0.0, double f2=0.0, double z1=0.0, double z2=0.0, int dpaletteId=DUMMY_PALETTE_ID, Option_t *option="colfz")
virtual size_t size() const
void SetPlotStyle(int paletteId=1)
double fabs(const Complex &x)
sprintf(tfres,"(1/%g)x(%g) (sec)x(Hz)", 2 *df, df)
STFT(wavearray< double > x, int nfft, int noverlap, TString ztype="amplitude", TString fwindow="hann", double fparam=0.0, TString name="stft")