Fix for autoreconf failures.
[flac.git] / src / flac / analyze.c
1 /* flac - Command-line FLAC encoder/decoder
2  * Copyright (C) 2000,2001,2002,2003,2004,2005,2006,2007,2008,2009  Josh Coalson
3  *
4  * This program is free software; you can redistribute it and/or
5  * modify it under the terms of the GNU General Public License
6  * as published by the Free Software Foundation; either version 2
7  * of the License, or (at your option) any later version.
8  *
9  * This program is distributed in the hope that it will be useful,
10  * but WITHOUT ANY WARRANTY; without even the implied warranty of
11  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
12  * GNU General Public License for more details.
13  *
14  * You should have received a copy of the GNU General Public License
15  * along with this program; if not, write to the Free Software
16  * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.
17  */
18
19 #if HAVE_CONFIG_H
20 #  include <config.h>
21 #endif
22
23 #include <errno.h>
24 #include <math.h>
25 #include <stdio.h>
26 #include <stdlib.h>
27 #include <string.h>
28 #include <inttypes.h>
29 #include "FLAC/all.h"
30 #include "analyze.h"
31
32 typedef struct {
33         FLAC__int32 residual;
34         unsigned count;
35 } pair_t;
36
37 typedef struct {
38         pair_t buckets[FLAC__MAX_BLOCK_SIZE];
39         int peak_index;
40         unsigned nbuckets;
41         unsigned nsamples;
42         double sum, sos;
43         double variance;
44         double mean;
45         double stddev;
46 } subframe_stats_t;
47
48 static subframe_stats_t all_;
49
50 static void init_stats(subframe_stats_t *stats);
51 static void update_stats(subframe_stats_t *stats, FLAC__int32 residual, unsigned incr);
52 static void compute_stats(subframe_stats_t *stats);
53 static FLAC__bool dump_stats(const subframe_stats_t *stats, const char *filename);
54
55 void flac__analyze_init(analysis_options aopts)
56 {
57         if(aopts.do_residual_gnuplot) {
58                 init_stats(&all_);
59         }
60 }
61
62 void flac__analyze_frame(const FLAC__Frame *frame, unsigned frame_number, FLAC__uint64 frame_offset, unsigned frame_bytes, analysis_options aopts, FILE *fout)
63 {
64         const unsigned channels = frame->header.channels;
65         char outfilename[1024];
66         subframe_stats_t stats;
67         unsigned i, channel, partitions;
68
69         /* do the human-readable part first */
70         fprintf(fout, "frame=%u\toffset=%" PRIu64 "\tbits=%u\tblocksize=%u\tsample_rate=%u\tchannels=%u\tchannel_assignment=%s\n", frame_number, frame_offset, frame_bytes*8, frame->header.blocksize, frame->header.sample_rate, channels, FLAC__ChannelAssignmentString[frame->header.channel_assignment]);
71         for(channel = 0; channel < channels; channel++) {
72                 const FLAC__Subframe *subframe = frame->subframes+channel;
73                 const FLAC__bool is_rice2 = subframe->data.fixed.entropy_coding_method.type == FLAC__ENTROPY_CODING_METHOD_PARTITIONED_RICE2;
74                 const unsigned pesc = is_rice2? FLAC__ENTROPY_CODING_METHOD_PARTITIONED_RICE2_ESCAPE_PARAMETER : FLAC__ENTROPY_CODING_METHOD_PARTITIONED_RICE_ESCAPE_PARAMETER;
75                 fprintf(fout, "\tsubframe=%u\twasted_bits=%u\ttype=%s", channel, subframe->wasted_bits, FLAC__SubframeTypeString[subframe->type]);
76                 switch(subframe->type) {
77                         case FLAC__SUBFRAME_TYPE_CONSTANT:
78                                 fprintf(fout, "\tvalue=%d\n", subframe->data.constant.value);
79                                 break;
80                         case FLAC__SUBFRAME_TYPE_FIXED:
81                                 FLAC__ASSERT(subframe->data.fixed.entropy_coding_method.type <= FLAC__ENTROPY_CODING_METHOD_PARTITIONED_RICE2);
82                                 fprintf(fout, "\torder=%u\tresidual_type=%s\tpartition_order=%u\n", subframe->data.fixed.order, is_rice2? "RICE2":"RICE", subframe->data.fixed.entropy_coding_method.data.partitioned_rice.order);
83                                 for(i = 0; i < subframe->data.fixed.order; i++)
84                                         fprintf(fout, "\t\twarmup[%u]=%d\n", i, subframe->data.fixed.warmup[i]);
85                                 partitions = (1u << subframe->data.fixed.entropy_coding_method.data.partitioned_rice.order);
86                                 for(i = 0; i < partitions; i++) {
87                                         unsigned parameter = subframe->data.fixed.entropy_coding_method.data.partitioned_rice.contents->parameters[i];
88                                         if(parameter == pesc)
89                                                 fprintf(fout, "\t\tparameter[%u]=ESCAPE, raw_bits=%u\n", i, subframe->data.fixed.entropy_coding_method.data.partitioned_rice.contents->raw_bits[i]);
90                                         else
91                                                 fprintf(fout, "\t\tparameter[%u]=%u\n", i, parameter);
92                                 }
93                                 if(aopts.do_residual_text) {
94                                         for(i = 0; i < frame->header.blocksize-subframe->data.fixed.order; i++)
95                                                 fprintf(fout, "\t\tresidual[%u]=%d\n", i, subframe->data.fixed.residual[i]);
96                                 }
97                                 break;
98                         case FLAC__SUBFRAME_TYPE_LPC:
99                                 FLAC__ASSERT(subframe->data.lpc.entropy_coding_method.type <= FLAC__ENTROPY_CODING_METHOD_PARTITIONED_RICE2);
100                                 fprintf(fout, "\torder=%u\tqlp_coeff_precision=%u\tquantization_level=%d\tresidual_type=%s\tpartition_order=%u\n", subframe->data.lpc.order, subframe->data.lpc.qlp_coeff_precision, subframe->data.lpc.quantization_level, is_rice2? "RICE2":"RICE", subframe->data.lpc.entropy_coding_method.data.partitioned_rice.order);
101                                 for(i = 0; i < subframe->data.lpc.order; i++)
102                                         fprintf(fout, "\t\tqlp_coeff[%u]=%d\n", i, subframe->data.lpc.qlp_coeff[i]);
103                                 for(i = 0; i < subframe->data.lpc.order; i++)
104                                         fprintf(fout, "\t\twarmup[%u]=%d\n", i, subframe->data.lpc.warmup[i]);
105                                 partitions = (1u << subframe->data.lpc.entropy_coding_method.data.partitioned_rice.order);
106                                 for(i = 0; i < partitions; i++) {
107                                         unsigned parameter = subframe->data.lpc.entropy_coding_method.data.partitioned_rice.contents->parameters[i];
108                                         if(parameter == pesc)
109                                                 fprintf(fout, "\t\tparameter[%u]=ESCAPE, raw_bits=%u\n", i, subframe->data.lpc.entropy_coding_method.data.partitioned_rice.contents->raw_bits[i]);
110                                         else
111                                                 fprintf(fout, "\t\tparameter[%u]=%u\n", i, parameter);
112                                 }
113                                 if(aopts.do_residual_text) {
114                                         for(i = 0; i < frame->header.blocksize-subframe->data.lpc.order; i++)
115                                                 fprintf(fout, "\t\tresidual[%u]=%d\n", i, subframe->data.lpc.residual[i]);
116                                 }
117                                 break;
118                         case FLAC__SUBFRAME_TYPE_VERBATIM:
119                                 fprintf(fout, "\n");
120                                 break;
121                 }
122         }
123
124         /* now do the residual distributions if requested */
125         if(aopts.do_residual_gnuplot) {
126                 for(channel = 0; channel < channels; channel++) {
127                         const FLAC__Subframe *subframe = frame->subframes+channel;
128                         unsigned residual_samples;
129
130                         init_stats(&stats);
131
132                         switch(subframe->type) {
133                                 case FLAC__SUBFRAME_TYPE_FIXED:
134                                         residual_samples = frame->header.blocksize - subframe->data.fixed.order;
135                                         for(i = 0; i < residual_samples; i++)
136                                                 update_stats(&stats, subframe->data.fixed.residual[i], 1);
137                                         break;
138                                 case FLAC__SUBFRAME_TYPE_LPC:
139                                         residual_samples = frame->header.blocksize - subframe->data.lpc.order;
140                                         for(i = 0; i < residual_samples; i++)
141                                                 update_stats(&stats, subframe->data.lpc.residual[i], 1);
142                                         break;
143                                 default:
144                                         break;
145                         }
146
147                         /* update all_ */
148                         for(i = 0; i < stats.nbuckets; i++) {
149                                 update_stats(&all_, stats.buckets[i].residual, stats.buckets[i].count);
150                         }
151
152                         /* write the subframe */
153                         sprintf(outfilename, "f%06u.s%u.gp", frame_number, channel);
154                         compute_stats(&stats);
155
156                         (void)dump_stats(&stats, outfilename);
157                 }
158         }
159 }
160
161 void flac__analyze_finish(analysis_options aopts)
162 {
163         if(aopts.do_residual_gnuplot) {
164                 compute_stats(&all_);
165                 (void)dump_stats(&all_, "all");
166         }
167 }
168
169 void init_stats(subframe_stats_t *stats)
170 {
171         stats->peak_index = -1;
172         stats->nbuckets = 0;
173         stats->nsamples = 0;
174         stats->sum = 0.0;
175         stats->sos = 0.0;
176 }
177
178 void update_stats(subframe_stats_t *stats, FLAC__int32 residual, unsigned incr)
179 {
180         unsigned i;
181         const double r = (double)residual, a = r*incr;
182
183         stats->nsamples += incr;
184         stats->sum += a;
185         stats->sos += (a*r);
186
187         for(i = 0; i < stats->nbuckets; i++) {
188                 if(stats->buckets[i].residual == residual) {
189                         stats->buckets[i].count += incr;
190                         goto find_peak;
191                 }
192         }
193         /* not found, make a new bucket */
194         i = stats->nbuckets;
195         stats->buckets[i].residual = residual;
196         stats->buckets[i].count = incr;
197         stats->nbuckets++;
198 find_peak:
199         if(stats->peak_index < 0 || stats->buckets[i].count > stats->buckets[stats->peak_index].count)
200                 stats->peak_index = i;
201 }
202
203 void compute_stats(subframe_stats_t *stats)
204 {
205         stats->mean = stats->sum / (double)stats->nsamples;
206         stats->variance = (stats->sos - (stats->sum * stats->sum / stats->nsamples)) / stats->nsamples;
207         stats->stddev = sqrt(stats->variance);
208 }
209
210 FLAC__bool dump_stats(const subframe_stats_t *stats, const char *filename)
211 {
212         FILE *outfile;
213         unsigned i;
214         const double m = stats->mean;
215         const double s1 = stats->stddev, s2 = s1*2, s3 = s1*3, s4 = s1*4, s5 = s1*5, s6 = s1*6;
216         const double p = stats->buckets[stats->peak_index].count;
217
218         outfile = fopen(filename, "w");
219
220         if(0 == outfile) {
221                 fprintf(stderr, "ERROR opening %s: %s\n", filename, strerror(errno));
222                 return false;
223         }
224
225         fprintf(outfile, "plot '-' title 'PDF', '-' title 'mean' with impulses, '-' title '1-stddev' with histeps, '-' title '2-stddev' with histeps, '-' title '3-stddev' with histeps, '-' title '4-stddev' with histeps, '-' title '5-stddev' with histeps, '-' title '6-stddev' with histeps\n");
226
227         for(i = 0; i < stats->nbuckets; i++) {
228                 fprintf(outfile, "%d %u\n", stats->buckets[i].residual, stats->buckets[i].count);
229         }
230         fprintf(outfile, "e\n");
231
232         fprintf(outfile, "%f %f\ne\n", stats->mean, p);
233         fprintf(outfile, "%f %f\n%f %f\ne\n", m-s1, p*0.8, m+s1, p*0.8);
234         fprintf(outfile, "%f %f\n%f %f\ne\n", m-s2, p*0.7, m+s2, p*0.7);
235         fprintf(outfile, "%f %f\n%f %f\ne\n", m-s3, p*0.6, m+s3, p*0.6);
236         fprintf(outfile, "%f %f\n%f %f\ne\n", m-s4, p*0.5, m+s4, p*0.5);
237         fprintf(outfile, "%f %f\n%f %f\ne\n", m-s5, p*0.4, m+s5, p*0.4);
238         fprintf(outfile, "%f %f\n%f %f\ne\n", m-s6, p*0.3, m+s6, p*0.3);
239
240         fprintf(outfile, "pause -1 'waiting...'\n");
241
242         fclose(outfile);
243         return true;
244 }