make audio more dynamic with low pass filter on amplitude scale
[maximus:dr1.git] / stringthing / stringthing.c
1 /*
2     stringthing -- pseudo-generative spinning string drone physical model
3     Copyright (C) 2011-10-18 Claude Heiland-Allen <claudiusmaximus@goto10.org>
4
5     This program is free software: you can redistribute it and/or modify
6     it under the terms of the GNU Affero General Public License as
7     published by the Free Software Foundation, either version 3 of the
8     License, or (at your option) any later version.
9
10     This program is distributed in the hope that it will be useful,
11     but WITHOUT ANY WARRANTY; without even the implied warranty of
12     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
13     GNU Affero General Public License for more details.
14
15     You should have received a copy of the GNU Affero General Public License
16     along with this program.  If not, see <http://www.gnu.org/licenses/>.
17
18 */
19
20 #include <ctype.h>
21 #include <math.h>
22 #include <stdio.h>
23 #include <stdlib.h>
24 #include <string.h>
25 #include <jack/jack.h>
26 #include <GL/glew.h>
27 #include <GL/glut.h>
28
29 #define PI 3.141592653589793
30 #define DT 0.0625
31 #define FRICTION 0.9999
32 #define GRAVITY 2
33
34 static double SR = 48000;
35 static double HZ;
36 #define MINN 32
37 #define MAXN 2048
38 static int N = 1024;
39
40 struct elem {
41   double x, y, z, dx, dy, dz;
42 };
43 static struct elem string[2][MAXN];
44 static volatile int w = 0;
45 static jack_default_audio_sample_t buf[MAXN][2];
46
47 static jack_client_t *client;
48 static jack_port_t *port[2];
49
50 static int phase = -1;
51 static int t = 0;
52 static double f1 = 0;
53 static double f2 = 0;
54 static double p1 = 0;
55 static double p2 = 0;
56 static double lo = 0;
57 static double lh = 0;
58 static double lb = 0;
59 static double ro = 0;
60 static double rh = 0;
61 static double rb = 0;
62
63 static float vstring[MAXN][4];
64 static int ufudge = 0;
65 static int uangle = 0;
66
67 static int width = 512;
68 static int height = 288;
69 static int fullscreen = 0;
70
71 static void displaycb(void) {
72   int s = w;
73   float tot = 0;
74   for (int i = 0; i < N; ++i) {
75     vstring[i][0] = string[s][i].dx; tot += string[s][i].dx * string[s][i].dx;
76     vstring[i][1] = string[s][i].dy; tot += string[s][i].dy * string[s][i].dy;
77     vstring[i][2] = string[s][i].dz; tot += string[s][i].dz * string[s][i].dz;
78     vstring[i][3] = 0.0;
79   }
80   tot = sqrt(tot / (3 * N));
81   for (int i = 0; i < N; ++i) {
82     vstring[i][0] /= tot;
83     vstring[i][1] /= tot;
84     vstring[i][2] /= tot;
85     vstring[i][3] /= tot;
86   }
87   glTexSubImage1D(GL_TEXTURE_1D, 0, 0, N, GL_RGBA, GL_FLOAT, vstring);
88   glUniform1fARB(ufudge, (N - 1) / (float) MAXN);
89   glUniform1fARB(uangle, atan2(vstring[N-1][1], vstring[N-1][0]));
90   glBegin(GL_QUADS); {
91     float a = width * 1.0 / height;
92     float f = N * 1.25 /  MAXN;
93     float x = a * f;
94     float y =     f;
95     glTexCoord2f( x, -y); glVertex2f(1, 0);
96     glTexCoord2f( x,  y); glVertex2f(1, 1);
97     glTexCoord2f(-x,  y); glVertex2f(0, 1);
98     glTexCoord2f(-x, -y); glVertex2f(0, 0);
99   } glEnd();
100   glutSwapBuffers();
101   glutReportErrors();
102 }
103
104 static void reshapecb(int w, int h) {
105   width = w;
106   height = h;
107   glViewport(0, 0, width, height);
108   glLoadIdentity();
109   gluOrtho2D(0, 1, 1, 0);
110   glutPostRedisplay();
111 }
112
113 static void timercb(int v) {
114   glutPostRedisplay();
115   glutTimerFunc(20, timercb, v);
116 }
117
118 static double f1s[256];
119 static double f2s[256];
120 static void keyboard1cb(unsigned char key, int x, int y) {
121   if (isalpha(key)) {
122     if (isupper(key)) {
123       f1s[tolower(key)] = f1;
124       f2s[tolower(key)] = f2;
125     } else {
126       f1 = f1s[key];
127       f2 = f2s[key];
128     }
129   } else {
130     switch (key) {
131     case ' ': { double temp = f1; f1 = f2; f2 = -temp; break; }
132     }
133   }
134 }
135
136 static void keyboard2cb(int key, int x, int y) {
137   switch (key) {
138   case GLUT_KEY_F11:
139     fullscreen = !fullscreen;
140     if (fullscreen) {
141       glutFullScreen();
142       glutSetCursor(GLUT_CURSOR_NONE);
143     } else {
144       glutReshapeWindow(512, 288);
145       glutSetCursor(GLUT_CURSOR_INHERIT);
146     }
147     break;
148   }
149 }
150
151 static void mousecb(int x, int y) {
152   double f = (x - width  / 2) * 2.0 / width;
153   double d = (height / 2 - y) * 4.0 / height;
154   double m = (pow(2, fabs(f)) - 1) * (f < 0 ? -1 : 1);
155   double r = pow(2, d);
156   f1 = (m * r) / 8.0;
157   f2 = (m / r) / 8.0;
158 }
159
160 static double ls = 0;
161
162 static int processcb(jack_nframes_t nframes, void *arg) {
163   jack_default_audio_sample_t *out[2];
164   out[0] = (jack_default_audio_sample_t *)
165     jack_port_get_buffer(port[0], nframes);
166   out[1] = (jack_default_audio_sample_t *)
167     jack_port_get_buffer(port[1], nframes);
168   for (unsigned int k = 0; k < nframes; ++k, ++phase) {
169     if (phase >= 2 * N) {
170       phase = 0;
171     }
172     if (phase == 0) {
173       double s = 0;
174       for (int i = 0; i < N; ++i) {
175         s += string[w][i].dx * string[w][i].dx
176           +  string[w][i].dy * string[w][i].dy
177           +  string[w][i].dz * string[w][i].dz;
178       }
179       s = 4 * sqrt(s / N);
180       ls = 0.99 * ls + 0.01 * s;
181       if (!(ls > 0)) { ls = 1; }
182       double a0 = -atan2(string[w][N-1].dy, string[w][N-1].dx);
183       double s0 = sin(a0);
184       double c0 = cos(a0);
185       for (int j = 0; j < 2 * N; ++j) {
186         int i = j >= N ? 2 * N - 1 - j : j;
187         double x =  c0 * string[w][i].dx + s0 * string[w][i].dy;
188         double y = -s0 * string[w][i].dx + c0 * string[w][i].dy;
189         double lv = x * x + string[w][i].dz * string[w][i].dz;
190         double rv = y * y + string[w][i].dz * string[w][i].dz;
191         double wd = (1 - cos(j * PI / N)) / 2;
192         lv = sqrt(lv) / ls;
193         rv = sqrt(rv) / ls;
194         lo = 0.999 * lo + 0.001 * lv;
195         ro = 0.999 * ro + 0.001 * rv;
196         lh = lv - lo;
197         rh = rv - ro;
198         lb = tanh(0.5 * lb + 0.5 * wd * lh);
199         rb = tanh(0.5 * rb + 0.5 * wd * rh);
200         buf[j][0] = tanh(2 * lb);
201         buf[j][1] = tanh(2 * rb);
202       }
203       double a1 = 1 / (1 + f1 * f2);
204       double a2 = 1 / (1 + f1 * f2);
205       for (int tt = 0; tt < 256; ++tt) {
206         p1 += f1 * N / SR;
207         p2 += f2 * N / SR;
208         while (p1 > 1) { p1 -= 1; } while (p1 < 0) { p1 += 1; }
209         while (p2 > 1) { p2 -= 1; } while (p2 < 0) { p2 += 1; }
210         {
211           double fx = a1 * cos(2 * PI * p1) - a2 * sin(2 * PI * p2);
212           double fy = a1 * sin(2 * PI * p1) + a2 * cos(2 * PI * p2);
213           double fz = 0;
214           string[1-w][0].dx = FRICTION * string[w][0].dx + fx * DT;
215           string[1-w][0].dy = FRICTION * string[w][0].dy + fy * DT;
216           string[1-w][0].dz = FRICTION * string[w][0].dz + fz * DT;
217           string[1-w][0].x = string[w][0].x + string[1-w][0].dx * DT;
218           string[1-w][0].y = string[w][0].y + string[1-w][0].dy * DT;
219           string[1-w][0].z = string[w][0].z + string[1-w][0].dz * DT;
220         }
221         for (int i = 1; i < N - 1; ++i) {
222           double fx = string[w][i-1].x + string[w][i+1].x - 2 * string[w][i].x;
223           double fy = string[w][i-1].y + string[w][i+1].y - 2 * string[w][i].y;
224           double fz = string[w][i-1].z + string[w][i+1].z - 2 * string[w][i].z
225                     - GRAVITY;
226           string[1-w][i].dx = FRICTION * string[w][i].dx + fx * DT;
227           string[1-w][i].dy = FRICTION * string[w][i].dy + fy * DT;
228           string[1-w][i].dz = FRICTION * string[w][i].dz + fz * DT;
229           string[1-w][i].x = string[w][i].x + string[1-w][i].dx * DT;
230           string[1-w][i].y = string[w][i].y + string[1-w][i].dy * DT;
231           string[1-w][i].z = string[w][i].z + string[1-w][i].dz * DT;
232         }
233         for (int i = N - 1; i < N; ++i) {
234           double fx = string[w][i-1].x - string[w][i].x;
235           double fy = string[w][i-1].y - string[w][i].y;
236           double fz = string[w][i-1].z - string[w][i].z - GRAVITY;
237           string[1-w][i].dx = FRICTION * string[w][i].dx + fx * DT;
238           string[1-w][i].dy = FRICTION * string[w][i].dy + fy * DT;
239           string[1-w][i].dz = FRICTION * string[w][i].dz + fz * DT;
240           string[1-w][i].x = string[w][i].x + string[1-w][i].dx * DT;
241           string[1-w][i].y = string[w][i].y + string[1-w][i].dy * DT;
242           string[1-w][i].z = string[w][i].z + string[1-w][i].dz * DT;
243         }
244         w = 1 - w;
245       }
246       ++t;
247     }
248     out[0][k] = buf[phase][0];
249     out[1][k] = buf[phase][1];
250   }
251   return 0;
252 }
253
254 static int sratecb(jack_nframes_t nframes, void *arg) {
255   SR = nframes;
256   N = fmin(fmax(SR / HZ, MINN), MAXN);
257   return 0;
258 }
259
260 static void errorcb(const char *desc) {
261   fprintf(stderr, "JACK error: %s\n", desc);
262 }
263
264 static void shutdowncb(void *arg) {
265   exit(1);
266 }
267
268 static void exitcb(void) {
269   jack_client_close(client);
270 }
271
272 static const char *src =
273 "uniform sampler1D tex;\n"
274 "uniform sampler2D palette;\n"
275 "uniform float fudge;\n"
276 "uniform float angle;\n"
277 "void main() {\n"
278 "  vec2 q = gl_TexCoord[0].xy;\n"
279 "  float l = length(q);\n"
280 "  float a = atan(q.y, q.x) + angle;\n"
281 "  float s = sin(a);\n"
282 "  float c = cos(a);\n"
283 "  mat2 m = mat2(c, s, -s, c);\n"
284 "  float r = clamp(l, 0.0, fudge);\n"
285 "  vec2 p = texture1D(tex, r).rg;\n"
286 "  vec4 rgba = texture2D(palette, p * m * 0.125 + vec2(0.5));\n"
287 "  rgba.rgb *= pow(clamp(r / l, 0.0, 1.0), 2.0);\n"
288 "  gl_FragColor = rgba;\n"
289 "}\n";
290
291 int main(int argc, char **argv) {
292   HZ = 0;
293   if (argc > 1) {
294     HZ = atof(argv[1]);
295   }
296   if (! HZ) { HZ = 96000.0 / 1024.0; }
297   memset(string, 0, 2 * MAXN * sizeof(struct elem));
298   memset(buf, 0, 2 * MAXN * sizeof(jack_default_audio_sample_t));
299   for (int i = 0; i < N; ++i) {
300     string[w][i].z = -i;
301   }
302   memset(f1s, 0, 256 * sizeof(double));
303   memset(f2s, 0, 256 * sizeof(double));
304   /* initialize OpenGL */
305   glutInitWindowSize(width, height);
306   glutInit(&argc, argv);
307   glutInitDisplayMode(GLUT_RGBA | GLUT_DOUBLE);
308   glutCreateWindow("stringthing");
309   glewInit();
310   glShadeModel(GL_FLAT);
311   GLuint texp;
312   glActiveTexture(GL_TEXTURE1);
313   glGenTextures(1, &texp);
314   glEnable(GL_TEXTURE_2D);
315   glBindTexture(GL_TEXTURE_2D, texp);
316   unsigned char *bufp = malloc(512 * 512 * 4);
317   double gamma = sqrt(2);
318   for (int j = 0; j < 512; ++j) {
319     for (int i = 0; i < 512; ++i) {
320       double x = (i - 256) / 64.0;
321       double y = (j - 256) / 64.0;
322       double h = atan2(y, x) * 0.15915494309189535;
323       double s = fmin(0.125 * sqrt(x * x + y * y), 1.0);
324       double v = 0.0625 * (x * x + y * y);
325       h -= floor(h);
326       h *= 6;
327       int k = floor(h);
328       double f = h - k;
329       double p = v * (1 - s);
330       double q = v * (1 - s * f);
331       double t = v * (1 - s * (1 - f));
332       double r, g, b;
333       switch (k) {
334       case  0: r = v; g = t; b = p; break;
335       case  1: r = q; g = v; b = p; break;
336       case  2: r = p; g = v; b = t; break;
337       case  3: r = p; g = q; b = v; break;
338       case  4: r = t; g = p; b = v; break;
339       default: r = v; g = p; b = q; break;
340       }
341       double z = s / (0.001 + 0.299 * r + 0.587 * g + 0.114 * b);
342       double w = pow(256, 1.0/gamma) * 256 * z / (1 + pow(fmax(x * x + y * y, 16), 2));
343       bufp[(512 * j + i)*4 + 0] = fmin(fmax(pow(w * r, gamma), 0), 255);
344       bufp[(512 * j + i)*4 + 1] = fmin(fmax(pow(w * g, gamma), 0), 255);
345       bufp[(512 * j + i)*4 + 2] = fmin(fmax(pow(w * b, gamma), 0), 255);
346       bufp[(512 * j + i)*4 + 3] = 1;
347     }
348   }
349   glTexImage2D(GL_TEXTURE_2D, 0, GL_RGBA, 512, 512, 0, GL_RGBA, GL_UNSIGNED_BYTE, bufp);
350   free(bufp);
351   glGenerateMipmap(GL_TEXTURE_2D);
352   glTexParameteri(GL_TEXTURE_2D, GL_TEXTURE_WRAP_S, GL_CLAMP_TO_EDGE);
353   glTexParameteri(GL_TEXTURE_2D, GL_TEXTURE_WRAP_T, GL_CLAMP_TO_EDGE);
354   glActiveTexture(GL_TEXTURE0);
355   GLuint tex;
356   glGenTextures(1, &tex);
357   glEnable(GL_TEXTURE_1D);
358   glBindTexture(GL_TEXTURE_1D, tex);
359   glTexImage1D(GL_TEXTURE_1D, 0, GL_RGBA32F_ARB, MAXN, 0, GL_RGBA, GL_FLOAT, 0);
360   glTexParameteri(GL_TEXTURE_1D, GL_TEXTURE_MIN_FILTER, GL_NEAREST);
361   glTexParameteri(GL_TEXTURE_1D, GL_TEXTURE_MAG_FILTER, GL_NEAREST);
362
363   /* compile shaders */
364   GLint success;
365   int prog = glCreateProgramObjectARB();
366   int frag = glCreateShaderObjectARB(GL_FRAGMENT_SHADER_ARB);
367   glShaderSourceARB(frag, 1, (const GLcharARB **) &src, 0);
368   glCompileShaderARB(frag);
369   glAttachObjectARB(prog, frag);
370   glLinkProgramARB(prog);
371   glGetObjectParameterivARB(prog, GL_OBJECT_LINK_STATUS_ARB, &success);
372   if (!success) exit(1);
373   glUseProgramObjectARB(prog);
374   glUniform1iARB(glGetUniformLocationARB(prog, "tex"), 0);
375   glUniform1iARB(glGetUniformLocationARB(prog, "palette"), 1);
376   ufudge = glGetUniformLocationARB(prog, "fudge");
377   uangle = glGetUniformLocationARB(prog, "angle");
378
379   /* set up callbacks */
380   glutDisplayFunc(displaycb);
381   glutReshapeFunc(reshapecb);
382   glutMotionFunc(mousecb);
383   glutPassiveMotionFunc(mousecb);
384   glutKeyboardFunc(keyboard1cb);
385   glutSpecialFunc(keyboard2cb);
386   glutTimerFunc(20, timercb, 0);
387
388   /* initialize JACK */
389   jack_set_error_function(errorcb);
390   if (!(client = jack_client_open("stringthing", 0, 0))) {
391     return 1;
392   }
393   atexit(exitcb);
394   jack_set_process_callback(client, processcb, 0);
395   jack_set_sample_rate_callback(client, sratecb, 0);
396   jack_on_shutdown(client, shutdowncb, 0);
397   port[0] = jack_port_register(
398     client, "output_1", JACK_DEFAULT_AUDIO_TYPE, JackPortIsOutput, 0
399   );
400   port[1] = jack_port_register(
401     client, "output_2", JACK_DEFAULT_AUDIO_TYPE, JackPortIsOutput, 0
402   );
403   if (jack_activate(client)) {
404     fprintf (stderr, "cannot activate JACK client");
405     return 1;
406   }
407   const char **ports;
408   if ((ports = jack_get_ports(
409     client, NULL, NULL, JackPortIsPhysical | JackPortIsInput
410   ))) {
411     int i = 0;
412     while (ports[i] && i < 2) {
413       if (jack_connect(
414         client, jack_port_name(port[i]), ports[i]
415       )) {
416         fprintf(stderr, "cannot connect output port\n");
417       }
418       i++;
419     }
420     free(ports);
421   }
422
423   glutMainLoop();
424   return 0;
425 }