?  search:
summary | log | graphiclog1 | graphiclog2 | commit | commitdiff | tree | refs | edit | fork
blame | history | raw | HEAD
Recognizes if input is ogg or not.
[xiph.git] / spectrum / process.c
blobf265eb0203470681ba15fed0309181849c41e10f
1 /*
2 *
3 * gtk2 spectrum analyzer
4 *
5 * Copyright (C) 2004 Monty
6 *
7 * This analyzer is free software; you can redistribute it and/or modify
8 * it under the terms of the GNU General Public License as published by
9 * the Free Software Foundation; either version 2, or (at your option)
10 * any later version.
11 *
12 * The analyzer is distributed in the hope that it will be useful,
13 * but WITHOUT ANY WARRANTY; without even the implied warranty of
14 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 * GNU General Public License for more details.
16 *
17 * You should have received a copy of the GNU General Public License
18 * along with Postfish; see the file COPYING. If not, write to the
19 * Free Software Foundation, 675 Mass Ave, Cambridge, MA 02139, USA.
20 *
21 *
22 */
23
24 #include "analyzer.h"
25
26 static int blockslice[MAX_FILES]= {-1,-1,-1,-1,-1,-1,-1,-1, -1,-1,-1,-1,-1,-1,-1,-1};
27
28 static float **blockbuffer=0;
29 static int blockbufferfill[MAX_FILES]={0,0,0,0,0,0,0,0, 0,0,0,0,0,0,0,0};
30 static float *window;
31 static float *freqbuffer=0;
32 static fftwf_plan plan;
33
34 static unsigned char readbuffer[MAX_FILES][readbuffersize];
35 static int readbufferfill[MAX_FILES]={0,0,0,0,0,0,0,0, 0,0,0,0,0,0,0,0};
36 static int readbufferptr[MAX_FILES]={0,0,0,0,0,0,0,0, 0,0,0,0,0,0,0,0};
37
38 static FILE *f[MAX_FILES]={0,0,0,0,0,0,0,0, 0,0,0,0,0,0,0,0};
39 static off_t offset[MAX_FILES]={0,0,0,0,0,0,0,0, 0,0,0,0,0,0,0,0};
40 static off_t length[MAX_FILES]= {-1,-1,-1,-1,-1,-1,-1,-1, -1,-1,-1,-1,-1,-1,-1,-1};
41 static off_t bytesleft[MAX_FILES]= {-1,-1,-1,-1,-1,-1,-1,-1, -1,-1,-1,-1,-1,-1,-1,-1};
42 int seekable[MAX_FILES]={0,0,0,0,0,0,0,0, 0,0,0,0,0,0,0,0};
43 int global_seekable=0;
44
45 pthread_mutex_t feedback_mutex=PTHREAD_RECURSIVE_MUTEX_INITIALIZER_NP;
46 int feedback_increment=0;
47
48 float *feedback_count;
49 float **plot_data;
50 float *plot_floor=NULL;
51 float **work_floor=NULL;
52 float *process_work;
53
54 float **feedback_acc;
55 float **feedback_max;
56 float **feedback_instant;
57
58 /* Gentlemen, power up the Variance hammer */
59 float **floor_y;
60 float **floor_yy;
61 int floor_count;
62
63 float **ph_acc;
64 float **ph_max;
65 float **ph_instant;
66
67 float **xmappingL;
68 float **xmappingH;
69 int metascale = -1;
70 int metawidth = -1;
71 int metares = -1;
72 int metanoise = 0;
73
74 sig_atomic_t acc_clear=0;
75 sig_atomic_t acc_rewind=0;
76 sig_atomic_t acc_loop=0;
77
78 sig_atomic_t process_active=0;
79 sig_atomic_t process_exit=0;
80
81 static int host_is_big_endian() {
82 int32_t pattern = 0xfeedface; /* deadbeef */
83 unsigned char *bytewise = (unsigned char *)&pattern;
84 if (bytewise[0] == 0xfe) return 1;
85 return 0;
86 }
87
88 /* Macros to read header data */
89 #define READ_U32_LE(buf) \
90 (((buf)[3]<<24)|((buf)[2]<<16)|((buf)[1]<<8)|((buf)[0]&0xff))
91
92 #define READ_U16_LE(buf) \
93 (((buf)[1]<<8)|((buf)[0]&0xff))
94
95 #define READ_U32_BE(buf) \
96 (((buf)[0]<<24)|((buf)[1]<<16)|((buf)[2]<<8)|((buf)[3]&0xff))
97
98 #define READ_U16_BE(buf) \
99 (((buf)[0]<<8)|((buf)[1]&0xff))
100
101 double read_IEEE80(unsigned char *buf){
102 int s=buf[0]&0xff;
103 int e=((buf[0]&0x7f)<<8)|(buf[1]&0xff);
104 double f=((unsigned long)(buf[2]&0xff)<<24)|
105 ((buf[3]&0xff)<<16)|
106 ((buf[4]&0xff)<<8) |
107 (buf[5]&0xff);
108
109 if(e==32767){
110 if(buf[2]&0x80)
111 return HUGE_VAL; /* Really NaN, but this won't happen in reality */
112 else{
113 if(s)
114 return -HUGE_VAL;
115 else
116 return HUGE_VAL;
117 }
118 }
119
120 f=ldexp(f,32);
121 f+= ((buf[6]&0xff)<<24)|
122 ((buf[7]&0xff)<<16)|
123 ((buf[8]&0xff)<<8) |
124 (buf[9]&0xff);
125
126 return ldexp(f, e-16446);
127 }
128
129 static int find_chunk(FILE *in, char *type, unsigned int *len, int endian){
130 unsigned int i;
131 unsigned char buf[8];
132
133 while(1){
134 if(fread(buf,1,8,in) <8)return 0;
135
136 if(endian)
137 *len = READ_U32_BE(buf+4);
138 else
139 *len = READ_U32_LE(buf+4);
140
141 if(memcmp(buf,type,4)){
142
143 if((*len) & 0x1)(*len)++;
144
145 for(i=0;i<*len;i++)
146 if(fgetc(in)==EOF)return 0;
147
148 }else return 1;
149 }
150 }
151
152 int input_load(void){
153
154 int stdinp=0,i,fi;
155 if(inputs==0){
156 /* look at stdin... is it a file, pipe, tty...? */
157 if(isatty(STDIN_FILENO)){
158 fprintf(stderr,
159 "Spectrum requires either an input file on the command line\n"
160 "or stream data piped|redirected to stdin. spectrum -h will\n"
161 "give more details.\n");
162 return 1;
163 }
164 stdinp=1; /* file coming in via stdin */
165 inputname[0]=strdup("stdin");
166 inputs++;
167 }
168
169 for(fi=0;fi<inputs;fi++){
170
171 if(stdinp && fi==0){
172 int newfd=dup(STDIN_FILENO);
173 f[fi]=fdopen(newfd,"rb");
174 }else{
175 f[fi]=fopen(inputname[fi],"rb");
176 }
177 seekable[fi]=0;
178
179 /* Crappy! Use a lib to do this for pete's sake! */
180 if(f[fi]){
181 char headerid[12];
182 off_t filelength;
183
184 /* parse header (well, sort of) and get file size */
185 seekable[fi]=(fseek(f[fi],0,SEEK_CUR)?0:1);
186
187 if(!seekable[fi]){
188 filelength=-1;
189 }else{
190 fseek(f[fi],0,SEEK_END);
191 filelength=ftello(f[fi]);
192 fseek(f[fi],0,SEEK_SET);
193 global_seekable=1;
194 }
195
196 fread(headerid,1,12,f[fi]);
197 if(!strncmp(headerid,"RIFF",4) && !strncmp(headerid+8,"WAVE",4)){
198 unsigned int chunklen;
199
200 if(find_chunk(f[fi],"fmt ",&chunklen,0)){
201 int ltype;
202 int lch;
203 int lrate;
204 int lbits;
205 unsigned char *buf=alloca(chunklen);
206
207 fread(buf,1,chunklen,f[fi]);
208
209 ltype = READ_U16_LE(buf);
210 lch = READ_U16_LE(buf+2);
211 lrate = READ_U32_LE(buf+4);
212 lbits = READ_U16_LE(buf+14);
213
214 if(ltype!=1){
215 fprintf(stderr,"%s:\n\tWAVE file not PCM.\n",inputname[fi]);
216 return 1;
217 }
218
219 if(bits[fi]==-1)bits[fi]=lbits;
220 if(channels[fi]==-1)channels[fi]=lch;
221 if(signedp[fi]==-1){
222 signedp[fi]=0;
223 if(bits[fi]>8)signedp[fi]=1;
224 }
225 if(bigendian[fi]==-1)bigendian[fi]=0;
226 if(rate[fi]==-1){
227 if(lrate<4000 || lrate>192000){
228 fprintf(stderr,"%s:\n\tSampling rate out of bounds\n",inputname[fi]);
229 return 1;
230 }
231 rate[fi]=lrate;
232 }
233
234 if(find_chunk(f[fi],"data",&chunklen,0)){
235 off_t pos=ftello(f[fi]);
236 int bytes=(bits[fi]+7)/8;
237 if(seekable[fi])
238 filelength=
239 (filelength-pos)/
240 (channels[fi]*bytes)*
241 (channels[fi]*bytes)+pos;
242
243 if(chunklen==0UL ||
244 chunklen==0x7fffffffUL ||
245 chunklen==0xffffffffUL){
246 if(filelength==-1){
247 length[fi]=-1;
248 fprintf(stderr,"%s: Incomplete header; assuming stream.\n",inputname[fi]);
249 }else{
250 length[fi]=(filelength-pos)/(channels[fi]*bytes);
251 fprintf(stderr,"%s: Incomplete header; using actual file size.\n",inputname[fi]);
252 }
253 }else if(filelength==-1 || chunklen+pos<=filelength){
254 length[fi]=(chunklen/(channels[fi]*bytes));
255 fprintf(stderr,"%s: Using declared file size (%ld).\n",
256 inputname[fi],(long)length[fi]*channels[fi]*bytes);
257
258 }else{
259
260 length[fi]=(filelength-pos)/(channels[fi]*bytes);
261 fprintf(stderr,"%s: File truncated; Using actual file size.\n",inputname[fi]);
262 }
263 offset[fi]=ftello(f[fi]);
264 } else {
265 fprintf(stderr,"%s: WAVE file has no \"data\" chunk following \"fmt \".\n",inputname[fi]);
266 return 1;
267 }
268 }else{
269 fprintf(stderr,"%s: WAVE file has no \"fmt \" chunk.\n",inputname[fi]);
270 return 1;
271 }
272
273 }else if(!strncmp(headerid,"FORM",4) && !strncmp(headerid+8,"AIF",3)){
274 unsigned int len;
275 int aifc=0;
276 if(headerid[11]=='C')aifc=1;
277 unsigned char *buffer;
278 char buf2[8];
279
280 int lch;
281 int lbits;
282 int lrate;
283 int bytes;
284
285 /* look for COMM */
286 if(!find_chunk(f[fi], "COMM", &len,1)){
287 fprintf(stderr,"%s: AIFF file has no \"COMM\" chunk.\n",inputname[fi]);
288 return 1;
289 }
290
291 if(len < 18 || (aifc && len<22)) {
292 fprintf(stderr,"%s: AIFF COMM chunk is truncated.\n",inputname[fi]);
293 return 1;
294 }
295
296 buffer = alloca(len);
297
298 if(fread(buffer,1,len,f[fi]) < len){
299 fprintf(stderr, "%s: Unexpected EOF in reading AIFF header\n",inputname[fi]);
300 return 1;
301 }
302
303 lch = READ_U16_BE(buffer);
304 lbits = READ_U16_BE(buffer+6);
305 lrate = (int)read_IEEE80(buffer+8);
306
307 if(bits[fi]==-1)bits[fi]=lbits;
308 bytes=(bits[fi]+7)/8;
309 if(signedp[fi]==-1)signedp[fi]=1;
310 if(rate[fi]==-1){
311 if(lrate<4000 || lrate>192000){
312 fprintf(stderr,"%s:\n\tSampling rate out of bounds\n",inputname[fi]);
313 return 1;
314 }
315 rate[fi]=lrate;
316 }
317 if(channels[fi]==-1)channels[fi]=lch;
318
319 if(bigendian[fi]==-1){
320 if(aifc){
321 if(!memcmp(buffer+18, "NONE", 4)) {
322 bigendian[fi] = 1;
323 }else if(!memcmp(buffer+18, "sowt", 4)) {
324 bigendian[fi] = 0;
325 }else{
326 fprintf(stderr, "%s: Spectrum supports only linear PCM AIFF-C files.\n",inputname[fi]);
327 return 1;
328 }
329 }else
330 bigendian[fi] = 1;
331 }
332 if(!find_chunk(f[fi], "SSND", &len, 1)){
333 fprintf(stderr,"%s: AIFF file has no \"SSND\" chunk.\n",inputname[fi]);
334 return 1;
335 }
336
337 if(fread(buf2,1,8,f[fi]) < 8){
338 fprintf(stderr,"%s: Unexpected EOF reading AIFF header\n",inputname[fi]);
339 return 1;
340 }
341
342 {
343 int loffset = READ_U32_BE(buf2);
344 int lblocksize = READ_U32_BE(buf2+4);
345
346 /* swallow some data */
347 for(i=0;i<loffset;i++)
348 if(fgetc(f[fi])==EOF)break;
349
350 if( lblocksize == 0 && (bits[fi] == 24 || bits[fi] == 16 || bits[fi] == 8)){
351
352 off_t pos=ftello(f[fi]);
353
354 if(seekable[fi])
355 filelength=
356 (filelength-pos)/
357 (channels[fi]*bytes)*
358 (channels[fi]*bytes)+pos;
359
360 if(len==0UL ||
361 len==0x7fffffffUL ||
362 len==0xffffffffUL){
363 if(filelength==-1){
364 length[fi]=-1;
365 fprintf(stderr,"%s: Incomplete header; assuming stream.\n",inputname[fi]);
366 }else{
367 length[fi]=(filelength-pos)/(channels[fi]*bytes);
368 fprintf(stderr,"%s: Incomplete header; using actual file size.\n",inputname[fi]);
369 }
370 }else if(filelength==-1 || (len+pos-loffset-8)<=filelength){
371 length[fi]=((len-loffset-8)/(channels[fi]*bytes));
372 fprintf(stderr,"%s: Using declared file size.\n",inputname[fi]);
373
374 }else{
375 length[fi]=(filelength-pos)/(channels[fi]*bytes);
376 fprintf(stderr,"%s: File truncated; Using actual file size.\n",inputname[fi]);
377 }
378 offset[fi]=pos;
379 }else{
380 fprintf(stderr, "%s: Spectrum supports only linear PCM AIFF-C files.\n",inputname[fi]);
381 return 1;
382 }
383 }
384 } else {
385 /* must be raw input */
386 fprintf(stderr,"Input has no header; assuming raw stream/file.\n");
387
388 if(channels[fi]==-1)channels[fi]=1;
389 if(rate[fi]==-1)rate[fi]=44100;
390 if(bits[fi]==-1)bits[fi]=16;
391 if(signedp[fi]==-1)signedp[fi]=1;
392 if(bigendian[fi]==-1)bigendian[fi]=host_is_big_endian();
393
394 offset[fi]=0;
395 length[fi]=-1;
396 if(seekable[fi])length[fi]=filelength/(channels[fi]*((bits[fi]+7)/8));
397
398 memcpy(readbuffer[fi],headerid,12);
399 readbufferfill[fi]=12;
400
401 }
402
403 /* select the full-block slice size: ~10fps */
404 blockslice[fi]=rate[fi]/10;
405 while(blockslice[fi]>blocksize/2)blockslice[fi]/=2;
406 total_ch += channels[fi];
407
408 if(length[fi]!=-1)
409 bytesleft[fi]=length[fi]*channels[fi]*((bits[fi]+7)/8);
410
411 }else{
412 fprintf(stderr,"Unable to open %s: %s\n",inputname[fi],strerror(errno));
413 exit(1);
414 }
415 }
416
417 blockbuffer=malloc(total_ch*sizeof(*blockbuffer));
418 process_work=calloc(blocksize+2,sizeof(*process_work));
419 feedback_count=calloc(total_ch,sizeof(*feedback_count));
420 plot_data=calloc(total_ch,sizeof(*plot_data));
421
422 feedback_acc=malloc(total_ch*sizeof(*feedback_acc));
423 feedback_max=malloc(total_ch*sizeof(*feedback_max));
424 feedback_instant=malloc(total_ch*sizeof(*feedback_instant));
425 floor_y=malloc(total_ch*sizeof(*floor_y));
426 floor_yy=malloc(total_ch*sizeof(*floor_yy));
427
428 ph_acc=malloc(total_ch*sizeof(*ph_acc));
429 ph_max=malloc(total_ch*sizeof(*ph_max));
430 ph_instant=malloc(total_ch*sizeof(*ph_instant));
431
432 freqbuffer=fftwf_malloc((blocksize+2)*sizeof(*freqbuffer));
433 for(i=0;i<total_ch;i++){
434 blockbuffer[i]=calloc(blocksize,sizeof(**blockbuffer));
435
436 floor_y[i]=calloc(blocksize/2+1,sizeof(**floor_y));
437 floor_yy[i]=calloc(blocksize/2+1,sizeof(**floor_yy));
438 feedback_acc[i]=calloc(blocksize/2+1,sizeof(**feedback_acc));
439 feedback_max[i]=calloc(blocksize/2+1,sizeof(**feedback_max));
440 feedback_instant[i]=calloc(blocksize/2+1,sizeof(**feedback_instant));
441
442 ph_acc[i]=calloc(blocksize+2,sizeof(**ph_acc));
443 ph_max[i]=calloc(blocksize+2,sizeof(**ph_max));
444 ph_instant[i]=calloc(blocksize+2,sizeof(**ph_instant));
445 }
446
447 plan=fftwf_plan_dft_r2c_1d(blocksize,freqbuffer,
448 (fftwf_complex *)freqbuffer,
449 FFTW_ESTIMATE);
450
451 /* construct proper window (sin^4 I'd think) */
452 window = calloc(blocksize,sizeof(*window));
453 for(i=0;i<blocksize;i++)window[i]=sin(M_PIl*i/blocksize);
454 for(i=0;i<blocksize;i++)window[i]*=window[i];
455 for(i=0;i<blocksize;i++)window[i]=sin(window[i]*M_PIl*.5);
456 for(i=0;i<blocksize;i++)window[i]*=window[i]/(blocksize/4)*.778;
457
458 return 0;
459
460 }
461
462 /* Convert new data from readbuffer into the blockbuffers until the
463 blockbuffer is full */
464 static void LBEconvert(void){
465 float scale=1./2147483648.;
466 int ch=0,fi;
467
468 for(fi=0;fi<inputs;fi++){
469 int bytes=(bits[fi]+7)/8;
470 int j;
471 int32_t xor=(signedp[fi]?0:0x80000000UL);
472
473 int readlimit=(readbufferfill[fi]-readbufferptr[fi])/
474 channels[fi]/bytes*channels[fi]*bytes+readbufferptr[fi];
475
476 int bfill = blockbufferfill[fi];
477 int rptr = readbufferptr[fi];
478 unsigned char *rbuf = readbuffer[fi];
479
480 if(readlimit){
481
482 switch(bytes){
483 case 1:
484
485 while(bfill<blocksize && rptr<readlimit){
486 for(j=ch;j<channels[fi]+ch;j++)
487 blockbuffer[j][bfill]=((rbuf[rptr++]<<24)^xor)*scale;
488 bfill++;
489 }
490 break;
491
492 case 2:
493
494 if(bigendian[fi]){
495 while(bfill<blocksize && rptr<readlimit){
496 for(j=ch;j<channels[fi]+ch;j++){
497 blockbuffer[j][bfill]=
498 (((rbuf[rptr+1]<<16)| (rbuf[rptr]<<24))^xor)*scale;
499 rptr+=2;
500 }
501 bfill++;
502 }
503 }else{
504 while(bfill<blocksize && rptr<readlimit){
505 for(j=ch;j<channels[fi]+ch;j++){
506 blockbuffer[j][bfill]=
507 (((rbuf[rptr]<<16)| (rbuf[rptr+1]<<24))^xor)*scale;
508 rptr+=2;
509 }
510 bfill++;
511 }
512 }
513 break;
514
515 case 3:
516
517 if(bigendian[fi]){
518 while(bfill<blocksize && rptr<readlimit){
519 for(j=ch;j<channels[fi]+ch;j++){
520 blockbuffer[j][bfill]=
521 (((rbuf[rptr+2]<<8)|(rbuf[rptr+1]<<16)|(rbuf[rptr]<<24))^xor)*scale;
522 rptr+=3;
523 }
524 bfill++;
525 }
526 }else{
527 while(bfill<blocksize && rptr<readlimit){
528 for(j=ch;j<channels[fi]+ch;j++){
529 blockbuffer[j][bfill]=
530 (((rbuf[rptr]<<8)|(rbuf[rptr+1]<<16)|(rbuf[rptr+2]<<24))^xor)*scale;
531 rptr+=3;
532 }
533 bfill++;
534 }
535 }
536 break;
537 case 4:
538
539 if(bigendian[fi]){
540 while(bfill<blocksize && rptr<readlimit){
541 for(j=ch;j<channels[fi]+ch;j++){
542 blockbuffer[j][bfill]=
543 (((rbuf[rptr+3])|(rbuf[rptr+2]<<8)|(rbuf[rptr+1]<<16)|(rbuf[rptr+3]<<24))^xor)*scale;
544 rptr+=4;
545 }
546 bfill++;
547 }
548 }else{
549 while(bfill<blocksize && rptr<readlimit){
550 for(j=ch;j<channels[fi]+ch;j++){
551 blockbuffer[j][bfill]=
552 (((rbuf[rptr])|(rbuf[rptr+1]<<8)|(rbuf[rptr+2]<<16)|(rbuf[rptr+3]<<24))^xor)*scale;
553 rptr+=4;
554 }
555 bfill++;
556 }
557 }
558 break;
559 }
560 }
561 ch+=channels[fi];
562 blockbufferfill[fi]=bfill;
563 readbufferptr[fi]=rptr;
564 }
565 }
566
567 /* when we return, the blockbuffer is full or we're at EOF */
568 /* EOF cases:
569 loop set: return EOF if all seekable streams have hit EOF
570 loop unset: return EOF if all streams have hit EOF
571 pad individual EOF streams out with zeroes until global EOF is hit */
572
573 static int input_read(void){
574 int i,fi,ch=0;
575 int eof=1;
576 int notdone=1;
577
578 for(fi=0;fi<inputs;fi++){
579
580 /* shift according to slice */
581 if(blockbufferfill[fi]==blocksize){
582 if(blockslice[fi]<blocksize){
583 for(i=0;i<channels[fi];i++)
584 memmove(blockbuffer[i+ch],blockbuffer[i+ch]+blockslice[fi],
585 (blocksize-blockslice[fi])*sizeof(**blockbuffer));
586 blockbufferfill[fi]-=blockslice[fi];
587 }else
588 blockbufferfill[fi]=0;
589 }
590 ch+=channels[fi];
591 }
592
593 while(notdone){
594 notdone=0;
595
596 /* if there's data left to be pulled out of a readbuffer, do that */
597 LBEconvert();
598
599 ch=0;
600 for(fi=0;fi<inputs;fi++){
601 if(blockbufferfill[fi]!=blocksize){
602
603 /* shift the read buffer before fill if there's a fractional
604 frame in it */
605 if(readbufferptr[fi]!=readbufferfill[fi] && readbufferptr[fi]>0){
606 memmove(readbuffer[fi],readbuffer[fi]+readbufferptr[fi],
607 (readbufferfill[fi]-readbufferptr[fi])*sizeof(**readbuffer));
608 readbufferfill[fi]-=readbufferptr[fi];
609 readbufferptr[fi]=0;
610 }else{
611 readbufferfill[fi]=0;
612 readbufferptr[fi]=0;
613 }
614
615 /* attempt to top off the readbuffer */
616 {
617 long actually_readbytes=0,readbytes=readbuffersize-readbufferfill[fi];
618
619 if(readbytes>0)
620 actually_readbytes=fread(readbuffer[fi]+readbufferfill[fi],1,readbytes,f[fi]);
621
622 if(actually_readbytes<0){
623 fprintf(stderr,"Input read error from %s: %s\n",inputname[fi],strerror(errno));
624 }else if (actually_readbytes==0){
625 /* don't process any partially-filled blocks; the
626 stairstep at the end could pollute results badly */
627
628 memset(readbuffer[fi],0,readbuffersize);
629 bytesleft[fi]=0;
630 readbufferfill[fi]=0;
631 readbufferptr[fi]=0;
632 blockbufferfill[fi]=0;
633
634 }else{
635 bytesleft[fi]-=actually_readbytes;
636 readbufferfill[fi]+=actually_readbytes;
637
638 /* conditionally clear global EOF */
639 if(acc_loop){
640 if(seekable[fi])eof=0;
641 }else{
642 eof=0;
643 }
644 notdone=1;
645 }
646 }
647 }
648 ch += channels[fi];
649 }
650 }
651 return eof;
652 }
653
654 void rundata_clear(){
655 int i,j;
656 for(i=0;i<total_ch;i++){
657 feedback_count[i]=0;
658 memset(feedback_acc[i],0,(blocksize/2+1)*sizeof(**feedback_acc));
659 memset(feedback_max[i],0,(blocksize/2+1)*sizeof(**feedback_max));
660 memset(feedback_instant[i],0,(blocksize/2+1)*sizeof(**feedback_instant));
661
662 for(j=0;j<blocksize+2;j++){
663 ph_acc[i][j]=0;
664 ph_max[i][j]=0;
665 ph_instant[i][j]=0;
666 }
667 }
668 acc_clear=0;
669 }
670
671 extern int plot_noise;
672
673 /* return 0 on EOF, 1 otherwise */
674 static int process(){
675 int fi,i,j,ch;
676 int eof_all;
677 int noise=plot_noise;
678
679 /* for each file, FOR SCIENCE! */
680 for(fi=0;fi<inputs;fi++){
681 if(acc_rewind && seekable[fi]){
682
683 blockbufferfill[fi]=0;
684 readbufferptr[fi]=0;
685 readbufferfill[fi]=0;
686 fseek(f[fi],offset[fi],SEEK_SET);
687 if(length[fi]!=-1)bytesleft[fi]=length[fi]*channels[fi]*((bits[fi]+7)/8);
688 }
689 }
690
691 eof_all=input_read();
692
693 if(eof_all){
694 if(acc_loop && !acc_rewind){
695 acc_rewind=1;
696 return process();
697 } else {
698 acc_rewind=0;
699 return 0;
700 }
701 }
702 acc_rewind=0;
703
704 if(acc_clear)
705 rundata_clear();
706
707 /* by channel */
708 ch=0;
709 for(fi=0;fi<inputs;fi++){
710 if(blockbufferfill[fi]){
711 for(i=ch;i<ch+channels[fi];i++){
712
713 float *data=blockbuffer[i];
714
715 /* window the blockbuffer into the FFT buffer */
716 for(j=0;j<blocksize;j++){
717 freqbuffer[j]=data[j]*window[j];
718 }
719
720 /* transform */
721 fftwf_execute(plan);
722
723 pthread_mutex_lock(&feedback_mutex);
724
725 /* perform desired accumulations */
726 for(j=0;j<blocksize+2;j+=2){
727 float R = freqbuffer[j];
728 float I = freqbuffer[j+1];
729 float sqR = R*R;
730 float sqI = I*I;
731 float sqM = sqR+sqI;
732
733 if(noise==1){
734 floor_yy[i][j>>1]+=sqM*sqM;
735 floor_y[i][j>>1]+=sqM;
736 }
737
738 /* deal with phase accumulate/rotate */
739 if(i==ch){
740 /* normalize/store ref for later rotation */
741 process_work[j] = R;
742 process_work[j+1] = -I;
743
744 }else{
745 /* rotate signed square phase according to ref for phase calculation */
746 float pR;
747 float pI;
748 float rR = process_work[j];
749 float rI = process_work[j+1];
750 pR = (rR*R - rI*I);
751 pI = (rR*I + rI*R);
752
753 ph_instant[i][j]=pR;
754 ph_instant[i][j+1]=pI;
755
756 ph_acc[i][j]+=pR;
757 ph_acc[i][j+1]+=pI;
758
759 if(feedback_max[i][j>>1]<sqM){
760 ph_max[i][j]=pR;
761 ph_max[i][j+1]=pI;
762 }
763 }
764
765 feedback_instant[i][j>>1]=sqM;
766 feedback_acc[i][j>>1]+=sqM;
767
768 if(feedback_max[i][j>>1]<sqM)
769 feedback_max[i][j>>1]=sqM;
770
771 }
772 feedback_count[i]++;
773
774 pthread_mutex_unlock(&feedback_mutex);
775 }
776 }
777 ch+=channels[fi];
778 }
779 if(noise==1)
780 floor_count++;
781 feedback_increment++;
782 write(eventpipe[1],"",1);
783 return 1;
784 }
785
786 void *process_thread(void *dummy){
787 while(!process_exit && process());
788 process_active=0;
789 write(eventpipe[1],"",1);
790 return NULL;
791 }
792
793 void process_dump(int mode){
794 int fi,i,j,ch;
795 FILE *out;
796
797 {
798 out=fopen("accumulate.m","w");
799 ch = 0;
800 for(fi=0;fi<inputs;fi++){
801 for(i=0;i<blocksize/2+1;i++){
802 fprintf(out,"%f ",(double)i*rate[fi]/blocksize);
803
804 for(j=ch;j<ch+channels[fi];j++)
805 fprintf(out,"%f ",todB(feedback_acc[j][i])*.5);
806 fprintf(out,"\n");
807 }
808 fprintf(out,"\n");
809 ch+=channels[fi];
810 }
811 fclose(out);
812 }
813
814 {
815 out=fopen("max.m","w");
816 ch = 0;
817 for(fi=0;fi<inputs;fi++){
818 for(i=0;i<blocksize/2+1;i++){
819 fprintf(out,"%f ",(double)i*rate[fi]/blocksize);
820
821 for(j=ch;j<ch+channels[fi];j++)
822 fprintf(out,"%f ",todB(feedback_max[j][i])*.5);
823 fprintf(out,"\n");
824 }
825 fprintf(out,"\n");
826 ch+=channels[fi];
827 }
828 fclose(out);
829 }
830
831 {
832 out=fopen("instant.m","w");
833 ch = 0;
834 for(fi=0;fi<inputs;fi++){
835 for(i=0;i<blocksize/2+1;i++){
836 fprintf(out,"%f ",(double)i*rate[fi]/blocksize);
837
838 for(j=ch;j<ch+channels[fi];j++)
839 fprintf(out,"%f ",todB(feedback_instant[j][i])*.5);
840 fprintf(out,"\n");
841 }
842 fprintf(out,"\n");
843 ch+=channels[fi];
844 }
845 fclose(out);
846 }
847
848 {
849 out=fopen("accphase.m","w");
850 ch = 0;
851 for(fi=0;fi<inputs;fi++){
852
853 /* phase */
854 for(i=0;i<blocksize+2;i+=2){
855 fprintf(out,"%f ",(double)i*.5*rate[fi]/blocksize);
856 fprintf(out,"%f ",atan2(ph_acc[ch+1][i+1],ph_acc[ch+1][i])*57.29);
857 fprintf(out,"\n");
858 }
859 fprintf(out,"\n");
860 ch+=channels[fi];
861 }
862 fclose(out);
863 }
864
865 }
866
867 void clear_noise_floor(){
868 int i;
869 for(i=0;i<total_ch;i++){
870 memset(floor_y[i],0,(blocksize/2+1)*sizeof(**floor_y));
871 memset(floor_yy[i],0,(blocksize/2+1)*sizeof(**floor_yy));
872 }
873 floor_count=0;
874 }
875
876 /* how many bins to 'trim' off the edge of calculated data when we
877 know we've hit a boundary of marginal measurement */
878 #define binspan 5
879
880 float **process_fetch(int res, int scale, int mode, int link,
881 int *active, int width,
882 float *ymax, float *pmax, float *pmin,
883 float **yfloor,int noise){
884 int ch,ci,i,j,fi;
885 float **data;
886 float **ph;
887
888 *yfloor=NULL;
889
890 /* are our scale mappings up to date? */
891 if(res != metares || scale != metascale || width != metawidth){
892 if(!xmappingL) xmappingL = calloc(inputs, sizeof(*xmappingL));
893 if(!xmappingH) xmappingH = calloc(inputs, sizeof(*xmappingH));
894 metanoise=-1;
895
896 if(!work_floor)
897 work_floor = calloc(total_ch,sizeof(*work_floor));
898 for(i=0;i<total_ch;i++){
899 if(work_floor[i])
900 work_floor[i] = realloc(work_floor[i],(width+1)*sizeof(**work_floor));
901 else
902 work_floor[i] = calloc((width+1),sizeof(**work_floor));
903 }
904
905 for(fi=0;fi<inputs;fi++){
906
907 /* if mapping preexists, resize it */
908 if(xmappingL[fi]){
909 xmappingL[fi] = realloc(xmappingL[fi],(width+1)*sizeof(**xmappingL));
910 }else{
911 xmappingL[fi] = malloc((width+1)*sizeof(**xmappingL));
912 }
913 if(xmappingH[fi]){
914 xmappingH[fi] = realloc(xmappingH[fi],(width+1)*sizeof(**xmappingH));
915 }else{
916 xmappingH[fi] = malloc((width+1)*sizeof(**xmappingH));
917 }
918
919 metascale = scale;
920 metawidth = width;
921 metares = res;
922
923
924 /* generate new numbers */
925 for(i=0;i<width;i++){
926 float off=0;
927 float loff=1.;
928 float hoff=1.;
929 float lfreq,hfreq;
930
931 switch(res){
932 case 0: /* screen-resolution */
933 off=1.;
934 break;
935 case 1: /* 1/24th octave */
936 loff = .95918945710913818816;
937 hoff = 1.04254690518999138632;
938 break;
939 case 2: /* 1/12th octave */
940 loff = .94387431268169349664;
941 hoff = 1.05946309435929526455;
942 break;
943 case 3: /* 1/3th octave */
944 loff = .79370052598409973738;
945 hoff = 1.25992104989487316475;
946 break;
947 }
948
949 switch(scale){
950 case 0: /* log */
951 lfreq= pow(10.,(i-off)/(width-1)
952 * (log10(100000.)-log10(5.))
953 + log10(5.)) * loff;
954 hfreq= pow(10.,(i+off)/(width-1)
955 * (log10(100000.)-log10(5.))
956 + log10(5.)) * hoff;
957 break;
958 case 1: /* ISO */
959 lfreq= pow(2.,(i-off)/(width-1)
960 * (log2(20000.)-log2(25.))
961 + log2(25.)) * loff;
962 hfreq= pow(2.,(i+off)/(width-1)
963 * (log2(20000.)-log2(25.))
964 + log2(25.)) *hoff;
965 break;
966 case 2: /* screen-resolution linear */
967 lfreq=(i-off)*20000./(width-1)*loff;
968 hfreq=(i+off)*20000./(width-1)*hoff;
969 break;
970 }
971
972 xmappingL[fi][i]=lfreq/(rate[fi]*.5)*(blocksize/2);
973 xmappingH[fi][i]=hfreq/(rate[fi]*.5)*(blocksize/2);
974
975 }
976
977 for(i=0;i<width;i++){
978 if(xmappingL[fi][i]<0.)xmappingL[fi][i]=0.;
979 if(xmappingL[fi][i]>blocksize/2.)xmappingL[fi][i]=blocksize/2.;
980 if(xmappingH[fi][i]<0.)xmappingH[fi][i]=0.;
981 if(xmappingH[fi][i]>blocksize/2.)xmappingH[fi][i]=blocksize/2.;
982 }
983 }
984
985 for(i=0;i<total_ch;i++)
986 if(plot_data[i]){
987 plot_data[i] = realloc(plot_data[i],(width+1)*sizeof(**plot_data));
988 }else{
989 plot_data[i] = malloc((width+1)*sizeof(**plot_data));
990 }
991 }
992
993 /* 'illustrate' the noise floor */
994 if(noise){
995 if(plot_floor)
996 plot_floor=realloc(plot_floor,(width+1)*sizeof(*plot_floor));
997 else
998 plot_floor=calloc((width+1),sizeof(*plot_floor));
999
1000 if(metanoise!=link){
1001 float *y = plot_floor;
1002 int ch=0;
1003 metanoise=link;
1004 for(i=0;i<width;i++)
1005 y[i]=-300;
1006
1007 for(fi=0;fi<inputs;fi++){
1008 float *L = xmappingL[fi];
1009 float *H = xmappingH[fi];
1010 float d = 1./floor_count;
1011
1012 for(ci=0;ci<channels[fi];ci++){
1013 float *fy = floor_y[ci+ch];
1014 float *fyy = floor_yy[ci+ch];
1015 float *w = work_floor[ci+ch];
1016
1017 for(i=0;i<width;i++){
1018 int first=floor(L[i]);
1019 int last=floor(H[i]);
1020 float esum;
1021 float vsum;
1022 float v = fyy[first]*floor_count - fy[first]*fy[first];
1023
1024 if(first==last){
1025 float del=H[i]-L[i];
1026 esum=fy[first]*del;
1027 vsum=v*del;
1028 }else{
1029 float del=1.-(L[i]-first);
1030 esum=fy[first]*del;
1031 vsum=v*del;
1032
1033 for(j=first+1;j<last;j++){
1034 v = fyy[j]*floor_count - fy[j]*fy[j];
1035 esum+=fy[j];
1036 vsum+=v;
1037 }
1038
1039 v = fyy[last]*floor_count - fy[last]*fy[last];
1040 del=(H[i]-last);
1041 esum+=fy[last]*del;
1042 vsum+=v*del;
1043 }
1044 vsum = 10*sqrt(vsum)*d;
1045 esum*=d;
1046 w[i] = esum+vsum*10;
1047 esum = todB_a(w+i)*.5;
1048
1049 if(esum>y[i])y[i]=esum;
1050 }
1051 }
1052 ch+=channels[fi];
1053 }
1054 }
1055 if(link == LINK_INDEPENDENT && mode==0)
1056 *yfloor=plot_floor;
1057 }else{
1058 for(i=0;i<total_ch;i++)
1059 memset(work_floor[i],0,width*sizeof(**work_floor));
1060 metanoise=-1;
1061 }
1062
1063 /* mode selects the base data set */
1064 switch(mode){
1065 case 0: /* independent / instant */
1066 data=feedback_instant;
1067 ph=ph_instant;
1068 break;
1069 case 1: /* independent / max */
1070 data=feedback_max;
1071 ph=ph_max;
1072 break;
1073 case 2:
1074 data=feedback_acc;
1075 ph=ph_acc;
1076 break;
1077 }
1078
1079 ch=0;
1080 *ymax = -150.;
1081 *pmax = -180.;
1082 *pmin = 180.;
1083 for(fi=0;fi<inputs;fi++){
1084 float *L = xmappingL[fi];
1085 float *H = xmappingH[fi];
1086
1087 switch(link){
1088 case LINK_INDEPENDENT:
1089
1090 for(ci=0;ci<channels[fi];ci++){
1091 float *y = plot_data[ci+ch];
1092 float *m = data[ci+ch];
1093 if(active[ch+ci]){
1094 for(i=0;i<width;i++){
1095 int first=floor(L[i]);
1096 int last=floor(H[i]);
1097 float sum;
1098
1099 if(first==last){
1100 float del=H[i]-L[i];
1101 sum=m[first]*del;
1102 }else{
1103 float del=1.-(L[i]-first);
1104 sum=m[first]*del;
1105
1106 for(j=first+1;j<last;j++)
1107 sum+=m[j];
1108
1109 del=(H[i]-last);
1110 sum+=m[last]*del;
1111 }
1112
1113 sum=todB_a(&sum)*.5;
1114 if(sum>*ymax)*ymax=sum;
1115 y[i]=sum;
1116 }
1117 }
1118 }
1119 break;
1120
1121 case LINK_SUMMED:
1122 {
1123 float *y = plot_data[ch];
1124 memset(y,0,(width+1)*sizeof(*y));
1125
1126 for(ci=0;ci<channels[fi];ci++){
1127 float *m = data[ci+ch];
1128 if(active[ch+ci]){
1129 for(i=0;i<width;i++){
1130 int first=floor(L[i]);
1131 int last=floor(H[i]);
1132
1133 if(first==last){
1134 float del=H[i]-L[i];
1135 y[i]+=m[first]*del;
1136 }else{
1137 float del=1.-(L[i]-first);
1138 y[i]+=m[first]*del;
1139
1140 for(j=first+1;j<last;j++)
1141 y[i]+=m[j];
1142
1143 del=(H[i]-last);
1144 y[i]+=m[last]*del;
1145 }
1146 }
1147 }
1148 }
1149
1150 for(i=0;i<width;i++){
1151 float sum=todB_a(y+i)*.5;
1152 if(sum>*ymax)*ymax=sum;
1153 y[i]=sum;
1154 }
1155 }
1156 break;
1157
1158 case LINK_SUB_FROM:
1159 {
1160 float *y = plot_data[ch];
1161 if(active[ch]==0){
1162 for(i=0;i<width;i++)
1163 y[i]=-300;
1164 }else{
1165 for(ci=0;ci<channels[fi];ci++){
1166 float *m = data[ci+ch];
1167 if(ci==0 || active[ch+ci]){
1168 for(i=0;i<width;i++){
1169 int first=floor(L[i]);
1170 int last=floor(H[i]);
1171 float sum;
1172
1173 if(first==last){
1174 float del=H[i]-L[i];
1175 sum=m[first]*del;
1176 }else{
1177 float del=1.-(L[i]-first);
1178 sum=m[first]*del;
1179
1180 for(j=first+1;j<last;j++)
1181 sum+=m[j];
1182
1183 del=(H[i]-last);
1184 sum+=m[last]*del;
1185 }
1186
1187 if(ci==0){
1188 y[i]=sum;
1189 }else{
1190 y[i]-=sum;
1191 }
1192 }
1193 }
1194 }
1195
1196 for(i=0;i<width;i++){
1197 float v = (y[i]>0?y[i]:0);
1198 float sum=todB_a(&v)*.5;
1199 if(sum>*ymax)*ymax=sum;
1200 y[i]=sum;
1201 }
1202 }
1203 }
1204 break;
1205 case LINK_SUB_REF:
1206 {
1207 float *r = plot_data[ch];
1208 for(ci=0;ci<channels[fi];ci++){
1209 float *y = plot_data[ch+ci];
1210 float *m = data[ci+ch];
1211 if(ci==0 || active[ch+ci]){
1212 for(i=0;i<width;i++){
1213 int first=floor(L[i]);
1214 int last=floor(H[i]);
1215 float sum;
1216
1217 if(first==last){
1218 float del=H[i]-L[i];
1219 sum=m[first]*del;
1220 }else{
1221 float del=1.-(L[i]-first);
1222 sum=m[first]*del;
1223
1224 for(j=first+1;j<last;j++)
1225 sum+=m[j];
1226
1227 del=(H[i]-last);
1228 sum+=m[last]*del;
1229 }
1230
1231 if(ci==0){
1232 r[i]=sum;
1233 }else{
1234 sum=(r[i]>sum?0.f:sum-r[i]);
1235 y[i]=todB_a(&sum)*.5;
1236 if(y[i]>*ymax)*ymax=y[i];
1237 }
1238 }
1239 }
1240 }
1241 }
1242 break;
1243
1244 case LINK_IMPEDENCE_p1:
1245 case LINK_IMPEDENCE_1:
1246 case LINK_IMPEDENCE_10:
1247 {
1248 float shunt = (link == LINK_IMPEDENCE_p1?.1:(link == LINK_IMPEDENCE_1?1:10));
1249 float *r = plot_data[ch];
1250
1251 for(ci=0;ci<channels[fi];ci++){
1252 float *y = plot_data[ci+ch];
1253 float *m = data[ch+ci];
1254
1255 if(ci==0 || active[ch+ci]){
1256 for(i=0;i<width;i++){
1257 int first=floor(L[i]);
1258 int last=floor(H[i]);
1259 float sum;
1260
1261 if(first==last){
1262 float del=H[i]-L[i];
1263 sum=m[first]*del;
1264 }else{
1265 float del=1.-(L[i]-first);
1266 sum=m[first]*del;
1267
1268 for(j=first+1;j<last;j++)
1269 sum+=m[j];
1270
1271 del=(H[i]-last);
1272 sum+=m[last]*del;
1273 }
1274
1275 if(ci==0){
1276 /* stash the reference in the work vector */
1277 r[i]=sum;
1278 }else{
1279 /* the shunt */
1280 /* 'r' collected at source, 'sum' across the shunt */
1281 float V=sqrt(r[i]);
1282 float S=sqrt(sum);
1283
1284 if(S>(1e-5) && V>S){
1285 y[i] = shunt*(V-S)/S;
1286 }else{
1287 y[i] = NAN;
1288 }
1289 }
1290 }
1291 }
1292 }
1293
1294 /* scan the resulting buffers for marginal data that would
1295 produce spurious output. Specifically we look for sharp
1296 falloffs of > 40dB or an original test magnitude under
1297 -70dB. */
1298 {
1299 float max = -140;
1300 for(i=0;i<width;i++){
1301 float v = r[i] = todB_a(r+i)*.5;
1302 if(v>max)max=v;
1303 }
1304
1305 for(ci=1;ci<channels[fi];ci++){
1306 if(active[ch+ci]){
1307 float *y = plot_data[ci+ch];
1308 for(i=0;i<width;i++){
1309 if(r[i]<max-40 || r[i]<-70){
1310 int j=i-binspan;
1311 if(j<0)j=0;
1312 for(;j<i;j++)
1313 y[j]=NAN;
1314 for(;j<width;j++){
1315 if(r[j]>max-40 && r[j]>-70)break;
1316 y[j]=NAN;
1317 }
1318 i=j+3;
1319 for(;j<i && j<width;j++){
1320 y[j]=NAN;
1321 }
1322 }
1323 if(!isnan(y[i]) && y[i]>*ymax)*ymax = y[i];
1324 }
1325 }
1326 }
1327 }
1328 }
1329 break;
1330
1331 case LINK_PHASE: /* response/phase */
1332
1333 if(channels[fi]>=2){
1334 float *om = plot_data[ch];
1335 float *op = plot_data[ch+1];
1336
1337 float *r = data[ch];
1338 float *rn = work_floor[ch];
1339 float *m = data[ch+1];
1340 float *mn = work_floor[ch+1];
1341 float *p = ph[ch+1];
1342 float mag[width];
1343
1344 if(feedback_count[ch]==0){
1345 memset(om,0,width*sizeof(*om));
1346 memset(op,0,width*sizeof(*op));
1347 }else{
1348 /* two vectors only; response and phase */
1349 /* response */
1350 if(active[ch] || active[ch+1]){
1351 for(i=0;i<width;i++){
1352 int first=floor(L[i]);
1353 int last=floor(H[i]);
1354 float sumR,sumM;
1355
1356 if(first==last){
1357 float del=H[i]-L[i];
1358 sumR=r[first]*del;
1359 sumM=m[first]*del;
1360 }else{
1361 float del=1.-(L[i]-first);
1362 sumR=r[first]*del;
1363 sumM=m[first]*del;
1364
1365 for(j=first+1;j<last;j++){
1366 sumR+=r[j];
1367 sumM+=m[j];
1368 }
1369
1370 del=(H[i]-last);
1371 sumR+=r[last]*del;
1372 sumM+=m[last]*del;
1373 }
1374
1375 if(sumR>rn[i] && sumM>mn[i]){
1376 mag[i] = todB_a(&sumR)*.5;
1377 sumM /= sumR;
1378 om[i] = todB_a(&sumM)*.5;
1379 }else{
1380 om[i] = NAN;
1381 }
1382 }
1383 }
1384
1385 /* phase */
1386 if(active[ch+1]){
1387 for(i=0;i<width;i++){
1388 int first=floor(L[i]);
1389 int last=floor(H[i]);
1390 float sumR,sumI;
1391
1392 if(first==last){
1393 float del=H[i]-L[i];
1394 sumR=p[(first<<1)]*del;
1395 sumI=p[(first<<1)+1]*del;
1396 }else{
1397 float del=1.-(L[i]-first);
1398 sumR=p[(first<<1)]*del;
1399 sumI=p[(first<<1)+1]*del;
1400
1401 for(j=first+1;j<last;j++){
1402 sumR+=p[(j<<1)];
1403 sumI+=p[(j<<1)+1];
1404 }
1405
1406 del=(H[i]-last);
1407 sumR+=p[(last<<1)]*del;
1408 sumI+=p[(last<<1)+1]*del;
1409 }
1410
1411 if(!isnan(om[i])){
1412 op[i] = atan2(sumI,sumR)*57.29;
1413 }else{
1414 op[i]=NAN;
1415 }
1416 }
1417 }
1418
1419 /* scan the resulting buffers for marginal data that would
1420 produce spurious output. Specifically we look for sharp
1421 falloffs of > 40dB or an original test magnitude under
1422 -70dB. */
1423 if(active[ch] || active[ch+1]){
1424 for(i=0;i<width;i++){
1425 if(isnan(om[i])){
1426 int j=i-binspan;
1427 if(j<0)j=0;
1428 for(;j<i;j++){
1429 om[j]=NAN;
1430 op[j]=NAN;
1431 }
1432 for(;j<width;j++){
1433 if(!isnan(om[j]))break;
1434 om[j]=NAN;
1435 op[j]=NAN;
1436 }
1437 i=j+3;
1438 for(;j<i && j<width;j++){
1439 om[j]=NAN;
1440 op[j]=NAN;
1441 }
1442 }
1443 if(om[i]>*ymax)*ymax = om[i];
1444 if(op[i]>*pmax)*pmax = op[i];
1445 if(op[i]<*pmin)*pmin = op[i];
1446
1447 }
1448 }
1449 }
1450 }
1451 break;
1452
1453 case LINK_THD: /* THD */
1454 case LINK_THD2: /* THD-2 */
1455 case LINK_THDN: /* THD+N */
1456 case LINK_THDN2: /* THD+N-2 */
1457
1458
1459 break;
1460
1461 }
1462 ch+=channels[fi];
1463 }
1464
1465 return plot_data;
1466 }
1467
Mirror of Xiph's Subversion repository