xref: /dports/audio/snd/snd-21.2/snd-fft.c

#include "snd.h"
#include "clm2xen.h"
#include "sndlib2xen.h"


/* -------------------------------- WAVELET TRANSFORM -------------------------------- */
/*
 * taken from wavelets.cl in clm which is taken from
 * M. J. Shensa Naval Ocean Systems Center,
 * Wickerhauser "Adapted Wavelet Analysis",
 * and the UBC Imager Wavelet Package by Bob Lewis
 * later stuff from Numerical Recipes
 */

static void wavelet_transform(mus_float_t *data, mus_long_t num, mus_float_t *cc, mus_long_t cc_size)
{
  mus_float_t *data1 = NULL;
  mus_float_t sig = -1.0;
  mus_float_t *cr;
  mus_long_t i, j, n, ii, ni, k, jf;

  cr = (mus_float_t *)malloc(cc_size * sizeof(mus_float_t));

  for (i = 0, j = cc_size - 1; i < cc_size; i++, j--)
    {
      cr[j] = sig * cc[i];
      sig = -sig;
    }

  for (n = num; n >= 4; n /= 2)
    {
      mus_long_t n1, nh, nmod, joff;
      if (data1) free(data1);
      data1 = (mus_float_t *)calloc(n, sizeof(mus_float_t));
      n1 = n - 1;
      nmod = cc_size * n;
      nh = n >> 1;
      joff = -(cc_size >> 1);

      for (ii = 0, i = 1; i <= n; i += 2, ii++)
	{
	  ni = i + nmod + joff;
	  for (k = 0; k < cc_size; k++)
	    {
	      jf = n1 & (ni + k + 1);
	      data1[ii] += cc[k] * data[jf];
	      data1[ii + nh] += cr[k] * data[jf];
	    }
	}
      mus_copy_floats(data, data1, n);
    }

  if (data1) free(data1);
  if (cr) free(cr);
}


/* these are from daub.h musicdsp.org, Computed by Kazuo Hatano, Aichi Institute of Technology. */
/* static mus_float_t Daub1[2] = { 7.07106781186e-01, 7.07106781186e-01}; */
static mus_float_t Daub2[4] = { 4.82962913144e-01, 8.36516303737e-01, 2.24143868042e-01,-1.29409522551e-01};
static mus_float_t Daub3[6] = { 3.32670552950e-01, 8.06891509311e-01, 4.59877502118e-01,-1.35011020010e-01,-8.54412738820e-02, 3.52262918857e-02};
static mus_float_t Daub4[8] = { 2.30377813308e-01, 7.14846570552e-01, 6.30880767929e-01,-2.79837694168e-02,-1.87034811719e-01, 3.08413818355e-02, 3.28830116668e-02,-1.05974017850e-02};
static mus_float_t Daub5[10] = { 1.60102397974e-01, 6.03829269797e-01, 7.24308528437e-01, 1.38428145901e-01,-2.42294887066e-01,-3.22448695846e-02, 7.75714938400e-02,-6.24149021279e-03,-1.25807519990e-02, 3.33572528547e-03};
static mus_float_t Daub6[12] = { 1.11540743350e-01, 4.94623890398e-01, 7.51133908021e-01, 3.15250351709e-01, -2.26264693965e-01,-1.29766867567e-01, 9.75016055873e-02, 2.75228655303e-02,-3.15820393174e-02, 5.53842201161e-04, 4.77725751094e-03,-1.07730108530e-03};
static mus_float_t Daub7[14] = { 7.78520540850e-02, 3.96539319481e-01, 7.29132090846e-01, 4.69782287405e-01,-1.43906003928e-01,-2.24036184993e-01, 7.13092192668e-02, 8.06126091510e-02,-3.80299369350e-02,-1.65745416306e-02, 1.25509985560e-02, 4.29577972921e-04,-1.80164070404e-03, 3.53713799974e-04};
static mus_float_t Daub8[16] = { 5.44158422431e-02, 3.12871590914e-01, 6.75630736297e-01, 5.85354683654e-01,-1.58291052563e-02,-2.84015542961e-01, 4.72484573913e-04, 1.28747426620e-01,-1.73693010018e-02,-4.40882539307e-02, 1.39810279173e-02, 8.74609404740e-03,-4.87035299345e-03,-3.91740373376e-04, 6.75449406450e-04,-1.17476784124e-04};
static mus_float_t Daub9[18] = { 3.80779473638e-02, 2.43834674612e-01, 6.04823123690e-01, 6.57288078051e-01, 1.33197385825e-01,-2.93273783279e-01,-9.68407832229e-02, 1.48540749338e-01, 3.07256814793e-02,-6.76328290613e-02, 2.50947114831e-04, 2.23616621236e-02,-4.72320475775e-03,-4.28150368246e-03, 1.84764688305e-03, 2.30385763523e-04,-2.51963188942e-04, 3.93473203162e-05};
static mus_float_t Daub10[20] = { 2.66700579005e-02, 1.88176800077e-01, 5.27201188931e-01, 6.88459039453e-01, 2.81172343660e-01,-2.49846424327e-01,-1.95946274377e-01, 1.27369340335e-01, 9.30573646035e-02,-7.13941471663e-02,-2.94575368218e-02, 3.32126740593e-02, 3.60655356695e-03,-1.07331754833e-02, 1.39535174705e-03, 1.99240529518e-03,-6.85856694959e-04,-1.16466855129e-04, 9.35886703200e-05,-1.32642028945e-05};
static mus_float_t Daub11[22] = { 1.86942977614e-02, 1.44067021150e-01, 4.49899764356e-01, 6.85686774916e-01, 4.11964368947e-01,-1.62275245027e-01,-2.74230846817e-01, 6.60435881966e-02, 1.49812012466e-01,-4.64799551166e-02,-6.64387856950e-02, 3.13350902190e-02, 2.08409043601e-02,-1.53648209062e-02,-3.34085887301e-03, 4.92841765605e-03,-3.08592858815e-04,-8.93023250666e-04, 2.49152523552e-04, 5.44390746993e-05,-3.46349841869e-05, 4.49427427723e-06};
static mus_float_t Daub12[24] = { 1.31122579572e-02, 1.09566272821e-01, 3.77355135214e-01, 6.57198722579e-01, 5.15886478427e-01,-4.47638856537e-02,-3.16178453752e-01,-2.37792572560e-02, 1.82478605927e-01, 5.35956967435e-03,-9.64321200965e-02, 1.08491302558e-02, 4.15462774950e-02,-1.22186490697e-02,-1.28408251983e-02, 6.71149900879e-03, 2.24860724099e-03,-2.17950361862e-03, 6.54512821250e-06, 3.88653062820e-04,-8.85041092082e-05,-2.42415457570e-05, 1.27769522193e-05,-1.52907175806e-06};
static mus_float_t Daub13[26] = { 9.20213353896e-03, 8.28612438729e-02, 3.11996322160e-01, 6.11055851158e-01, 5.88889570431e-01, 8.69857261796e-02,-3.14972907711e-01,-1.24576730750e-01, 1.79476079429e-01, 7.29489336567e-02,-1.05807618187e-01,-2.64884064753e-02, 5.61394771002e-02, 2.37997225405e-03,-2.38314207103e-02, 3.92394144879e-03, 7.25558940161e-03,-2.76191123465e-03,-1.31567391189e-03, 9.32326130867e-04, 4.92515251262e-05,-1.65128988556e-04, 3.06785375793e-05, 1.04419305714e-05,-4.70041647936e-06, 5.22003509845e-07};
static mus_float_t Daub14[28] = { 6.46115346008e-03, 6.23647588493e-02, 2.54850267792e-01, 5.54305617940e-01, 6.31187849104e-01, 2.18670687758e-01,-2.71688552278e-01,-2.18033529993e-01, 1.38395213864e-01, 1.39989016584e-01,-8.67484115681e-02,-7.15489555040e-02, 5.52371262592e-02, 2.69814083079e-02,-3.01853515403e-02,-5.61504953035e-03, 1.27894932663e-02,-7.46218989268e-04,-3.84963886802e-03, 1.06169108560e-03, 7.08021154235e-04,-3.86831947312e-04,-4.17772457703e-05, 6.87550425269e-05,-1.03372091845e-05,-4.38970490178e-06, 1.72499467536e-06,-1.78713996831e-07};
static mus_float_t Daub15[30] = { 4.53853736157e-03, 4.67433948927e-02, 2.06023863986e-01, 4.92631771708e-01, 6.45813140357e-01, 3.39002535454e-01,-1.93204139609e-01,-2.88882596566e-01, 6.52829528487e-02, 1.90146714007e-01,-3.96661765557e-02,-1.11120936037e-01, 3.38771439235e-02, 5.47805505845e-02,-2.57670073284e-02,-2.08100501696e-02, 1.50839180278e-02, 5.10100036040e-03,-6.48773456031e-03,-2.41756490761e-04, 1.94332398038e-03,-3.73482354137e-04,-3.59565244362e-04, 1.55896489920e-04, 2.57926991553e-05,-2.81332962660e-05, 3.36298718173e-06, 1.81127040794e-06,-6.31688232588e-07, 6.13335991330e-08};
static mus_float_t Daub16[32] = { 3.18922092534e-03, 3.49077143236e-02, 1.65064283488e-01, 4.30312722846e-01, 6.37356332083e-01, 4.40290256886e-01,-8.97510894024e-02,-3.27063310527e-01,-2.79182081330e-02, 2.11190693947e-01, 2.73402637527e-02,-1.32388305563e-01,-6.23972275247e-03, 7.59242360442e-02,-7.58897436885e-03,-3.68883976917e-02, 1.02976596409e-02, 1.39937688598e-02,-6.99001456341e-03,-3.64427962149e-03, 3.12802338120e-03, 4.07896980849e-04,-9.41021749359e-04, 1.14241520038e-04, 1.74787245225e-04,-6.10359662141e-05,-1.39456689882e-05, 1.13366086612e-05,-1.04357134231e-06,-7.36365678545e-07, 2.30878408685e-07,-2.10933963010e-08};
static mus_float_t Daub17[34] = { 2.24180700103e-03, 2.59853937036e-02, 1.31214903307e-01, 3.70350724152e-01, 6.10996615684e-01, 5.18315764056e-01, 2.73149704032e-02,-3.28320748363e-01,-1.26599752215e-01, 1.97310589565e-01, 1.01135489177e-01,-1.26815691778e-01,-5.70914196316e-02, 8.11059866541e-02, 2.23123361781e-02,-4.69224383892e-02,-3.27095553581e-03, 2.27336765839e-02,-3.04298998135e-03,-8.60292152032e-03, 2.96799669152e-03, 2.30120524215e-03,-1.43684530480e-03,-3.28132519409e-04, 4.39465427768e-04,-2.56101095665e-05,-8.20480320245e-05, 2.31868137987e-05, 6.99060098507e-06,-4.50594247722e-06, 3.01654960999e-07, 2.95770093331e-07,-8.42394844600e-08, 7.26749296856e-09};
static mus_float_t Daub18[36] = { 1.57631021844e-03, 1.92885317241e-02, 1.03588465822e-01, 3.14678941337e-01, 5.71826807766e-01, 5.71801654888e-01, 1.47223111969e-01,-2.93654040736e-01,-2.16480934005e-01, 1.49533975565e-01, 1.67081312763e-01,-9.23318841508e-02,-1.06752246659e-01, 6.48872162119e-02, 5.70512477385e-02,-4.45261419029e-02,-2.37332103958e-02, 2.66707059264e-02, 6.26216795430e-03,-1.30514809466e-02, 1.18630033858e-04, 4.94334360546e-03,-1.11873266699e-03,-1.34059629833e-03, 6.28465682965e-04, 2.13581561910e-04,-1.98648552311e-04,-1.53591712353e-07, 3.74123788074e-05,-8.52060253744e-06,-3.33263447888e-06, 1.76871298362e-06,-7.69163268988e-08,-1.17609876702e-07, 3.06883586304e-08,-2.50793445494e-09};
static mus_float_t Daub19[38] = { 1.10866976318e-03, 1.42810984507e-02, 8.12781132654e-02, 2.64388431740e-01, 5.24436377464e-01, 6.01704549127e-01, 2.60894952651e-01,-2.28091394215e-01,-2.85838631755e-01, 7.46522697081e-02, 2.12349743306e-01,-3.35185419023e-02,-1.42785695038e-01, 2.75843506256e-02, 8.69067555558e-02,-2.65012362501e-02,-4.56742262772e-02, 2.16237674095e-02, 1.93755498891e-02,-1.39883886785e-02,-5.86692228101e-03, 7.04074736710e-03, 7.68954359257e-04,-2.68755180070e-03, 3.41808653458e-04, 7.35802520505e-04,-2.60676135678e-04,-1.24600791734e-04, 8.71127046721e-05, 5.10595048707e-06,-1.66401762971e-05, 3.01096431629e-06, 1.53193147669e-06,-6.86275565776e-07, 1.44708829879e-08, 4.63693777578e-08,-1.11640206703e-08, 8.66684883899e-10};
static mus_float_t Daub20[40] = { 7.79953613666e-04, 1.05493946249e-02, 6.34237804590e-02, 2.19942113551e-01, 4.72696185310e-01, 6.10493238938e-01, 3.61502298739e-01,-1.39212088011e-01,-3.26786800434e-01,-1.67270883090e-02, 2.28291050819e-01, 3.98502464577e-02,-1.55458750707e-01,-2.47168273386e-02, 1.02291719174e-01, 5.63224685730e-03,-6.17228996246e-02, 5.87468181181e-03, 3.22942995307e-02,-8.78932492390e-03,-1.38105261371e-02, 6.72162730225e-03, 4.42054238704e-03,-3.58149425960e-03,-8.31562172822e-04, 1.39255961932e-03,-5.34975984399e-05,-3.85104748699e-04, 1.01532889736e-04, 6.77428082837e-05,-3.71058618339e-05,-4.37614386218e-06, 7.24124828767e-06,-1.01199401001e-06,-6.84707959700e-07, 2.63392422627e-07, 2.01432202355e-10,-1.81484324829e-08, 4.05612705555e-09,-2.99883648961e-10};
static mus_float_t Daub21[42] = { 5.48822509852e-04, 7.77663905235e-03, 4.92477715381e-02, 1.81359625440e-01, 4.19687944939e-01, 6.01506094935e-01, 4.44590451927e-01,-3.57229196172e-02,-3.35664089530e-01,-1.12397071568e-01, 2.11564527680e-01, 1.15233298439e-01,-1.39940424932e-01,-8.17759429808e-02, 9.66003903237e-02, 4.57234057492e-02,-6.49775048937e-02,-1.86538592021e-02, 3.97268354278e-02, 3.35775639033e-03,-2.08920536779e-02, 2.40347092080e-03, 8.98882438197e-03,-2.89133434858e-03,-2.95837403893e-03, 1.71660704063e-03, 6.39418500512e-04,-6.90671117082e-04,-3.19640627768e-05, 1.93664650416e-04,-3.63552025008e-05,-3.49966598498e-05, 1.53548250927e-05, 2.79033053981e-06,-3.09001716454e-06, 3.16609544236e-07, 2.99213663046e-07,-1.00040087903e-07,-2.25401497467e-09, 7.05803354123e-09,-1.47195419765e-09, 1.03880557102e-10};
static mus_float_t Daub22[44] = { 3.86263231491e-04, 5.72185463133e-03, 3.80699372364e-02, 1.48367540890e-01, 3.67728683446e-01, 5.78432731009e-01, 5.07901090622e-01, 7.37245011836e-02,-3.12726580428e-01,-2.00568406104e-01, 1.64093188106e-01, 1.79973187992e-01,-9.71107984091e-02,-1.31768137686e-01, 6.80763143927e-02, 8.45573763668e-02,-5.13642542974e-02,-4.65308118275e-02, 3.69708466206e-02, 2.05867076275e-02,-2.34800013444e-02,-6.21378284936e-03, 1.25647252183e-02, 3.00137398507e-04,-5.45569198615e-03, 1.04426073918e-03, 1.82701049565e-03,-7.70690988123e-04,-4.23787399839e-04, 3.28609414213e-04, 4.34589990453e-05,-9.40522363481e-05, 1.13743496621e-05, 1.73737569575e-05,-6.16672931646e-06,-1.56517913199e-06, 1.29518205731e-06,-8.77987987336e-08,-1.28333622875e-07, 3.76122874933e-08, 1.68017140492e-09,-2.72962314663e-09, 5.33593882166e-10,-3.60211348433e-11};
static mus_float_t Daub23[46] = { 2.71904194128e-04, 4.20274889318e-03, 2.93100036578e-02, 1.20515531783e-01, 3.18450813852e-01, 5.44931147873e-01, 5.51018517241e-01, 1.81392625363e-01,-2.61392148030e-01,-2.71402098607e-01, 9.21254070824e-02, 2.23573658242e-01,-3.30374470942e-02,-1.64011321531e-01, 2.02830745756e-02, 1.12297043618e-01,-2.11262123562e-02,-7.02073915749e-02, 2.17658568344e-02, 3.84953325225e-02,-1.85235136501e-02,-1.75371010030e-02, 1.27519439315e-02, 6.03184065002e-03,-7.07531927370e-03,-1.13486547335e-03, 3.12287644981e-03,-2.46501400516e-04,-1.06123122888e-03, 3.19420492709e-04, 2.56762452007e-04,-1.50021850349e-04,-3.37889483412e-05, 4.42607120310e-05,-2.63520788924e-06,-8.34787556785e-06, 2.39756954684e-06, 8.14757483477e-07,-5.33900540520e-07, 1.85309178563e-08, 5.41754917953e-08,-1.39993549543e-08,-9.47288590181e-10, 1.05044645369e-09,-1.93240511131e-10, 1.25020330235e-11};
static mus_float_t Daub24[48] = { 1.91435800947e-04, 3.08208171490e-03, 2.24823399497e-02, 9.72622358336e-02, 2.72908916067e-01, 5.04371040839e-01, 5.74939221095e-01, 2.80985553233e-01,-1.87271406885e-01,-3.17943078999e-01, 4.77661368434e-03, 2.39237388780e-01, 4.25287296414e-02,-1.71175351370e-01,-3.87771735779e-02, 1.21016303469e-01, 2.09801137091e-02,-8.21616542080e-02,-4.57843624181e-03, 5.13016200399e-02,-4.94470942812e-03,-2.82131070949e-02, 7.66172188164e-03, 1.30499708710e-02,-6.29143537001e-03,-4.74656878632e-03, 3.73604617828e-03, 1.15376493683e-03,-1.69645681897e-03,-4.41618485614e-05, 5.86127059318e-04,-1.18123323796e-04,-1.46007981776e-04, 6.55938863930e-05, 2.18324146046e-05,-2.02288829261e-05, 1.34115775080e-08, 3.90110033859e-06,-8.98025314393e-07,-4.03250775687e-07, 2.16633965327e-07,-5.05764541979e-10,-2.25574038817e-08, 5.15777678967e-09, 4.74837582425e-10,-4.02465864458e-10, 6.99180115763e-11,-4.34278250380e-12};
static mus_float_t Daub25[50] = { 1.34802979347e-04, 2.25695959185e-03, 1.71867412540e-02, 7.80358628721e-02, 2.31693507886e-01, 4.59683415146e-01, 5.81636896746e-01, 3.67885074802e-01,-9.71746409646e-02,-3.36473079641e-01,-8.75876145876e-02, 2.24537819745e-01, 1.18155286719e-01,-1.50560213750e-01,-9.85086152899e-02, 1.06633805018e-01, 6.67521644940e-02,-7.70841110565e-02,-3.71739628611e-02, 5.36179093987e-02, 1.55426059291e-02,-3.40423204606e-02,-3.07983679484e-03, 1.89228044766e-02,-1.98942578220e-03,-8.86070261804e-03, 2.72693625873e-03, 3.32270777397e-03,-1.84248429020e-03,-8.99977423746e-04, 8.77258193674e-04, 1.15321244046e-04,-3.09880099098e-04, 3.54371452327e-05, 7.90464000396e-05,-2.73304811996e-05,-1.27719529319e-05, 8.99066139306e-06, 5.23282770815e-07,-1.77920133265e-06, 3.21203751886e-07, 1.92280679014e-07,-8.65694173227e-08,-2.61159855611e-09, 9.27922448008e-09,-1.88041575506e-09,-2.22847491022e-10, 1.53590157016e-10,-2.52762516346e-11, 1.50969208282e-12};
static mus_float_t Daub26[52] = { 9.49379575071e-05, 1.65052023353e-03, 1.30975542925e-02, 6.22747440251e-02, 1.95039438716e-01, 4.13292962278e-01, 5.73669043034e-01, 4.39158311789e-01, 1.77407678098e-03,-3.26384593691e-01,-1.74839961289e-01, 1.81291832311e-01, 1.82755409589e-01,-1.04323900285e-01,-1.47977193275e-01, 6.98231861132e-02, 1.06482405249e-01,-5.34485616814e-02,-6.86547596040e-02, 4.22321857963e-02, 3.85357159711e-02,-3.13781103630e-02,-1.77609035683e-02, 2.07349201799e-02, 5.82958055531e-03,-1.17854979061e-02,-5.28738399262e-04, 5.60194723942e-03,-9.39058250473e-04,-2.14553028156e-03, 8.38348805654e-04, 6.16138220457e-04,-4.31955707426e-04,-1.06057474828e-04, 1.57479523860e-04,-5.27779549303e-06,-4.10967399639e-05, 1.07422154087e-05, 7.00007868296e-06,-3.88740016185e-06,-4.65046322064e-07, 7.93921063370e-07,-1.07900423757e-07,-8.90446637016e-08, 3.40779562129e-08, 2.16932825985e-09,-3.77601047853e-09, 6.78004724582e-10, 1.00230319104e-10,-5.84040818534e-11, 9.13051001637e-12,-5.25187122424e-13};
static mus_float_t Daub27[54] = { 6.68713138543e-05, 1.20553123167e-03, 9.95258878087e-03, 4.94525999829e-02, 1.62922027502e-01, 3.67110214125e-01, 5.53849860990e-01, 4.93406122677e-01, 1.02840855061e-01,-2.89716803314e-01,-2.48264581903e-01, 1.14823019517e-01, 2.27273288414e-01,-3.87864186318e-02,-1.78031740959e-01, 1.57993974602e-02, 1.31197971717e-01,-1.40627515558e-02,-9.10229065295e-02, 1.73110182654e-02, 5.79694057347e-02,-1.85124935619e-02,-3.27390666310e-02, 1.61469669223e-02, 1.56655956489e-02,-1.15771864589e-02,-5.86209634546e-03, 6.85663560968e-03, 1.34262687730e-03,-3.33285446952e-03, 1.45752962593e-04, 1.30117745024e-03,-3.41835122691e-04,-3.87901857410e-04, 2.01971987969e-04, 7.66005838706e-05,-7.71114551779e-05,-3.51748361490e-06, 2.06344264773e-05,-3.90116407063e-06,-3.65750090818e-06, 1.63436962472e-06, 3.05088068625e-07,-3.47246814739e-07, 3.28655896805e-08, 4.02625505286e-08,-1.32133227399e-08,-1.30946560685e-09, 1.52161498477e-09,-2.41552692801e-10,-4.37498622429e-11, 2.21366208806e-11,-3.29579012247e-12, 1.82818835288e-13};
static mus_float_t Daub28[56] = { 4.71080777501e-05, 8.79498515984e-04, 7.54265037764e-03, 3.90926081154e-02, 1.35137914253e-01, 3.22563361285e-01, 5.24998231630e-01, 5.30516293441e-01, 2.00176144045e-01,-2.30498954047e-01,-3.01327809532e-01, 3.28578791633e-02, 2.45808151373e-01, 3.69068853157e-02,-1.82877330732e-01,-4.68382337445e-02, 1.34627567910e-01, 3.44786312750e-02,-9.76853558056e-02,-1.73419228313e-02, 6.77478955019e-02, 3.44801895554e-03,-4.33333686160e-02, 4.43173291006e-03, 2.46880600101e-02,-6.81554976455e-03,-1.20635919682e-02, 5.83881662774e-03, 4.78486311245e-03,-3.72546124707e-03,-1.36037384563e-03, 1.87599866820e-03, 1.41567239314e-04,-7.48674955911e-04, 1.15465606365e-04, 2.29579098223e-04,-8.90390149004e-05,-4.90771341619e-05, 3.64140121105e-05, 4.63866498139e-06,-1.00432604133e-05, 1.24790031757e-06, 1.84036373451e-06,-6.67021547995e-07,-1.75746117320e-07, 1.49066001353e-07,-8.26238731562e-09,-1.78413869087e-08, 5.04404705638e-09, 6.94454032894e-10,-6.07704124722e-10, 8.49222001105e-11, 1.86736726378e-11,-8.36549047125e-12, 1.18885053340e-12,-6.36777235471e-14};
static mus_float_t Daub29[58] = { 3.31896627984e-05, 6.40951680304e-04, 5.70212651777e-03, 3.07735802214e-02, 1.11370116951e-01, 2.80653455970e-01, 4.89758804762e-01, 5.51374432758e-01, 2.89105238335e-01,-1.54028734459e-01,-3.30040948917e-01,-5.57068000729e-02, 2.36105236153e-01, 1.12419174873e-01,-1.60877988594e-01,-1.07845949938e-01, 1.14472295893e-01, 8.32207471624e-02,-8.51254926156e-02,-5.50274895253e-02, 6.34791645842e-02, 3.05315432727e-02,-4.51879812777e-02,-1.29171425542e-02, 2.94704318717e-02, 2.64832730767e-03,-1.70412245736e-02, 1.73788033272e-03, 8.46972549356e-03,-2.55080712778e-03,-3.47379898968e-03, 1.87712092572e-03, 1.08705394222e-03,-1.00077832708e-03,-2.00071136307e-04, 4.11128345474e-04,-2.29201804121e-05,-1.29304484008e-04, 3.64502606856e-05, 2.91334475016e-05,-1.65732839530e-05,-3.59364480402e-06, 4.75060924645e-06,-3.02905459205e-07,-8.97570175063e-07, 2.63389838699e-07, 9.38719741109e-08,-6.28615692201e-08, 1.07659190661e-09, 7.76897885477e-09,-1.89399538617e-09,-3.42680086326e-10, 2.40709945350e-10,-2.94058925076e-11,-7.83250973362e-12, 3.15276241337e-12,-4.28565487006e-13, 2.21919131158e-14};
static mus_float_t Daub30[60] = { 2.33861617273e-05, 4.66637950428e-04, 4.30079716504e-03, 2.41308326715e-02, 9.12383040670e-02, 2.42020670940e-01, 4.50487821853e-01, 5.57572232912e-01, 3.66242683371e-01,-6.61836707759e-02,-3.32966975020e-01,-1.41968513330e-01, 1.99462121580e-01, 1.77829873244e-01,-1.14558219432e-01,-1.57236817959e-01, 7.27786589703e-02, 1.22747746045e-01,-5.38064654582e-02,-8.76586900363e-02, 4.38016646714e-02, 5.67123657447e-02,-3.56733974967e-02,-3.22637589193e-02, 2.70786195952e-02, 1.52879607698e-02,-1.83997438681e-02,-5.29685966613e-03, 1.09156316583e-02, 6.19671756497e-04,-5.53073014819e-03, 8.43384586662e-04, 2.32452009406e-03,-8.60927696811e-04,-7.67878250438e-04, 5.05094823903e-04, 1.72482584235e-04,-2.16171830116e-04,-8.54830546758e-06, 6.98200837080e-05,-1.33971686329e-05,-1.63615247872e-05, 7.25214553589e-06, 2.32754909849e-06,-2.18726767699e-06, 1.09947433852e-08, 4.26166232601e-07,-1.00041468235e-07,-4.76437996513e-08, 2.60544275497e-08, 5.55339786139e-10,-3.33110568046e-09, 6.98486269183e-10, 1.61362297827e-10,-9.46138799727e-11, 1.00010513139e-11, 3.23942863853e-12,-1.18523759210e-12, 1.54399757084e-13,-7.73794263095e-15};
static mus_float_t Daub31[62] = { 1.64801338645e-05, 3.39412203776e-04, 3.23688406862e-03, 1.88536916129e-02, 7.43360930116e-02, 2.07012874485e-01, 4.09192200037e-01, 5.51139840914e-01, 4.29468808206e-01, 2.71692124973e-02,-3.10955118319e-01,-2.17978485523e-01, 1.40178288765e-01, 2.24966711473e-01,-4.99263491604e-02,-1.86962360895e-01, 1.54369884294e-02, 1.45089500931e-01,-8.13983227346e-03,-1.07612773323e-01, 1.09412974523e-02, 7.53536117432e-02,-1.48800266181e-02,-4.86190754648e-02, 1.61541715659e-02, 2.80476193667e-02,-1.42762752777e-02,-1.39005529392e-02, 1.05176394873e-02, 5.51616357331e-03,-6.52085237587e-03,-1.42826422321e-03, 3.39306677671e-03,-6.39790110601e-05,-1.45904174198e-03, 3.43139829690e-04, 4.99881617563e-04,-2.39658346940e-04,-1.24341161725e-04, 1.08958435041e-04, 1.50133572744e-05,-3.63125515786e-05, 4.03452023518e-06, 8.79530134269e-06,-3.03514236589e-06,-1.36906023094e-06, 9.81001542204e-07, 5.32725065697e-08,-1.97592512917e-07, 3.61682651733e-08, 2.32830971382e-08,-1.06152960215e-08,-6.47431168795e-10, 1.40856815102e-09,-2.52404395415e-10,-7.34893003248e-11, 3.69210880887e-11,-3.32700896712e-12,-1.32433491724e-12, 4.44546709629e-13,-5.55944205057e-14, 2.69938287976e-15};
static mus_float_t Daub32[64] = { 1.16146330213e-05, 2.46656690638e-04, 2.43126191957e-03, 1.46810463814e-02, 6.02574991203e-02, 1.75750783639e-01, 3.67509628597e-01, 5.34317919340e-01, 4.77809163733e-01, 1.20630538265e-01,-2.66698181476e-01,-2.77421581558e-01, 6.47133548055e-02, 2.48310642356e-01, 2.46624448396e-02,-1.92102344708e-01,-4.89951171846e-02, 1.45232079475e-01, 4.44049081999e-02,-1.09456113116e-01,-2.96278725084e-02, 8.08741406384e-02, 1.41061515161e-02,-5.69263140624e-02,-2.38026446493e-03, 3.70514579235e-02,-4.14590766082e-03,-2.16628228363e-02, 6.16752731068e-03, 1.10174007154e-02,-5.41156825727e-03,-4.64921675118e-03, 3.62722464068e-03, 1.46895510046e-03,-1.96474055582e-03,-2.21167872957e-04, 8.67305851845e-04,-1.02453731060e-04,-3.05965442382e-04, 1.05391546173e-04, 8.10367832913e-05,-5.25980928268e-05,-1.29404577940e-05, 1.82426840198e-05,-6.36178153226e-07,-4.55830957626e-06, 1.20288903632e-06, 7.56004762559e-07,-4.28597069315e-07,-5.00336186874e-08, 8.96596631195e-08,-1.21992435948e-08,-1.10438302172e-08, 4.25042231198e-09, 4.38438779994e-10,-5.88109146263e-10, 8.90472379622e-11, 3.26327074133e-11,-1.43091876516e-11, 1.07561065350e-12, 5.36148222961e-13,-1.66380048943e-13, 2.00071530381e-14,-9.42101913953e-16};
static mus_float_t Daub33[66] = { 8.18635831417e-06, 1.79101615370e-04, 1.82270943516e-03, 1.13959433745e-02, 4.86146665317e-02, 1.48186313180e-01, 3.26718130117e-01, 5.09376172514e-01, 5.11254770583e-01, 2.09582350713e-01,-2.04202622398e-01,-3.15997410766e-01,-1.92783394369e-02, 2.45420612119e-01, 9.98515586803e-02,-1.71428099051e-01,-1.10844133116e-01, 1.21967856403e-01, 9.47880880506e-02,-9.11469683513e-02,-7.03024850540e-02, 7.01911439409e-02, 4.57345618938e-02,-5.34712513358e-02,-2.52485829774e-02, 3.86870607602e-02, 1.07032658200e-02,-2.57287617547e-02,-2.16775861735e-03, 1.53169541158e-02,-1.59428878241e-03,-7.95354038705e-03, 2.38906240816e-03, 3.48080095340e-03,-1.86071821445e-03,-1.20430925760e-03, 1.07438069635e-03, 2.72730584733e-04,-4.90832900759e-04, 4.39316625176e-06, 1.78043189825e-04,-4.16043851627e-05,-4.92956442341e-05, 2.42333539881e-05, 9.07080575782e-06,-8.86612136675e-06,-3.60751610287e-07, 2.28837127614e-06,-4.42692340795e-07,-3.98579129198e-07, 1.82244333257e-07, 3.37797270373e-08,-3.98783819851e-08, 3.67286357683e-09, 5.11121185734e-09,-1.67139267725e-09,-2.49640210524e-10, 2.42683310230e-10,-3.04957445394e-11,-1.42023685988e-11, 5.50941472076e-12,-3.34348121895e-13,-2.15248838683e-13, 6.21474024717e-14,-7.19651054536e-15, 3.28937367841e-16};
static mus_float_t Daub34[68] = { 5.77051063273e-06, 1.29947620067e-04, 1.36406139005e-03, 8.81988940388e-03, 3.90488413517e-02, 1.24152482111e-01, 2.87765059233e-01, 4.78478746279e-01, 5.30555099656e-01, 2.90366329507e-01,-1.28246842174e-01,-3.31525301508e-01,-1.03891915515e-01, 2.16907220187e-01, 1.66601750412e-01,-1.27337358223e-01,-1.60924927177e-01, 7.79918469379e-02, 1.34125960271e-01,-5.44829680641e-02,-1.02947596992e-01, 4.35760946496e-02, 7.31852354367e-02,-3.70128384178e-02,-4.74385596452e-02, 3.07397465739e-02, 2.72283507563e-02,-2.36717379228e-02,-1.31439800166e-02, 1.64093741998e-02, 4.71364926099e-03,-1.00455067083e-02,-6.19474884515e-04, 5.33495076875e-03,-7.69212797506e-04,-2.39945394353e-03, 8.58995987436e-04, 8.75199906407e-04,-5.52735576214e-04,-2.32673214023e-04, 2.65077239755e-04, 2.66005001845e-05,-9.91469777078e-05, 1.35311722724e-05, 2.84495141969e-05,-1.05765749425e-05,-5.71082651099e-06, 4.16987175854e-06, 4.97971810142e-07,-1.11630653481e-06, 1.44819570833e-07, 2.02599066666e-07,-7.52670174041e-08,-1.99034650153e-08, 1.74042333293e-08,-8.66574426136e-10,-2.31650194699e-09, 6.44637821032e-10, 1.30041031860e-10,-9.90477453763e-11, 1.00420873546e-11, 6.08012535400e-12,-2.10787910891e-12, 9.79945115821e-14, 8.57919405179e-14,-2.31708370390e-14, 2.58733838193e-15,-1.14894475448e-16};
static mus_float_t Daub35[70] = { 4.06793406114e-06, 9.42146947557e-05, 1.01912268037e-03, 6.80729288431e-03, 3.12362885114e-02, 1.03404455861e-01, 2.51307378994e-01, 4.43592739224e-01, 5.37008427509e-01, 3.60345640518e-01,-4.38838818739e-02,-3.23822864912e-01,-1.81786976766e-01, 1.66041357490e-01, 2.17299289321e-01,-6.52628713106e-02,-1.91919589298e-01, 1.93095446660e-02, 1.55292480396e-01,-4.75268083411e-03,-1.20585522643e-01, 4.73422917264e-03, 8.99135475707e-02,-9.31855894990e-03,-6.33560374404e-02, 1.32285495850e-02, 4.12546930647e-02,-1.43668397842e-02,-2.41694978016e-02, 1.27664567156e-02, 1.22894360081e-02,-9.57779789923e-03,-5.08599164923e-03, 6.13775458674e-03, 1.42808879407e-03,-3.35764438092e-03, 7.61596943517e-06, 1.54963746970e-03,-3.34669216425e-04,-5.86481031899e-04, 2.64832881996e-04, 1.70001228366e-04,-1.36588307226e-04,-2.97699596284e-05, 5.30414312291e-05,-2.43700152682e-06,-1.57244207727e-05, 4.30804786171e-06, 3.35334586287e-06,-1.89592961769e-06,-3.90393173328e-07, 5.30236861690e-07,-3.70030837820e-08,-9.99039694453e-08, 3.00818865071e-08, 1.08490273378e-08,-7.45811655289e-09, 5.89795131038e-11, 1.03082334548e-09,-2.43354557375e-10,-6.40793825650e-11, 4.00053662725e-11,-3.12563935710e-12,-2.56706547615e-12, 8.01508853368e-13,-2.59795432889e-14,-3.39772085679e-14, 8.62403743472e-15,-9.29801252932e-16, 4.01462871233e-17};
static mus_float_t Daub36[72] = { 2.86792518275e-06, 6.82602867854e-05, 7.60215109966e-04, 5.24029737740e-03, 2.48905656448e-02, 8.56520925952e-02, 2.17756953097e-01, 4.06433697708e-01, 5.32266895260e-01, 4.17875335600e-01, 4.39751975293e-02,-2.94421039589e-01,-2.46807036978e-01, 9.81142041631e-02, 2.46537277608e-01, 7.27851509579e-03,-1.99337205608e-01,-4.58614007463e-02, 1.54106236627e-01, 5.02761800735e-02,-1.18803754310e-01,-3.98808535755e-02, 9.11567822580e-02, 2.50387214495e-02,-6.82090166368e-02,-1.13191003168e-02, 4.85130835478e-02, 1.42497266176e-03,-3.19807206776e-02, 3.98404019871e-03, 1.90635947806e-02,-5.65781324505e-03,-9.99026347328e-03, 5.02298910666e-03, 4.41348483535e-03,-3.48454144540e-03,-1.50307406629e-03, 1.99079377185e-03, 2.77681279571e-04,-9.46340382326e-04, 8.61456575899e-05, 3.69350728496e-04,-1.15511889584e-04,-1.13189946808e-04, 6.69474119693e-05, 2.37510668366e-05,-2.73139082465e-05,-1.18347105998e-06, 8.37221819816e-06,-1.58614578243e-06,-1.87081160285e-06, 8.31142127970e-07, 2.54842352255e-07,-2.45537765843e-07, 2.75324907333e-09, 4.79904346545e-08,-1.15609368881e-08,-5.61278434332e-09, 3.13884169578e-09, 1.09081555371e-10,-4.51254577856e-10, 8.96241820385e-11, 3.03742909811e-11,-1.59971668926e-11, 8.87684628721e-13, 1.07096935711e-12,-3.02928502697e-13, 5.54226318263e-15, 1.33807138629e-14,-3.20462854340e-15, 3.33997198481e-16,-1.40327417537e-17};
static mus_float_t Daub37[74] = { 2.02206086249e-06, 4.94234375062e-05, 5.66241837706e-04, 4.02414036825e-03, 1.97622861538e-02, 7.05848259771e-02, 1.87326331862e-01, 3.68440972400e-01, 5.18167040855e-01, 4.62207553661e-01, 1.30878963233e-01,-2.46180429761e-01,-2.94375915262e-01, 1.96715004523e-02, 2.51523254360e-01, 8.18060283872e-02,-1.81962291778e-01,-1.08451713823e-01, 1.29929646959e-01, 1.01780296838e-01,-9.66075406166e-02,-8.23302119065e-02, 7.50476199483e-02, 5.95674108715e-02,-5.92568156326e-02,-3.82538294793e-02, 4.58079441512e-02, 2.09728005925e-02,-3.35235840641e-02,-8.83349389041e-03, 2.26186515445e-02, 1.69047238348e-03,-1.37639819628e-02, 1.51930577883e-03, 7.38775745285e-03,-2.24805318700e-03,-3.39452327640e-03, 1.81687134380e-03, 1.26393425811e-03,-1.11148486531e-03,-3.28078847088e-04, 5.49053277337e-04, 1.53443902319e-05,-2.20894403245e-04, 4.33672612594e-05, 7.05513878206e-05,-3.09866292761e-05,-1.63916249616e-05, 1.35432771841e-05, 1.84994500311e-06,-4.30994155659e-06, 4.85473139699e-07, 1.00212139929e-06,-3.49494860344e-07,-1.50988538867e-07, 1.10903123221e-07, 5.35065751546e-09,-2.25219383672e-08, 4.22448570636e-09, 2.79397446595e-09,-1.29720500146e-09,-1.03141112909e-10, 1.94616489408e-10,-3.20339824412e-11,-1.39841571553e-11, 6.33495544097e-12,-2.09636319423e-13,-4.42161240987e-13, 1.13805283092e-13,-4.51888960746e-16,-5.24302569188e-15, 1.18901238750e-15,-1.19928033585e-16, 4.90661506493e-18};
static mus_float_t Daub38[76] = { 1.42577664167e-06, 3.57625199426e-05, 4.21170266472e-04, 3.08308811925e-03, 1.56372493475e-02, 5.78899436128e-02, 1.60071993564e-01, 3.30775781411e-01, 4.96591175311e-01, 4.93356078517e-01, 2.13050571355e-01,-1.82867667708e-01,-3.21675637808e-01,-6.22665060478e-02, 2.32125963835e-01, 1.49985119618e-01,-1.41795685973e-01,-1.59912565158e-01, 8.56381215561e-02, 1.41414734073e-01,-5.65864586307e-02,-1.14731170710e-01, 4.30958954330e-02, 8.72043982620e-02,-3.66051034028e-02,-6.17662087084e-02, 3.19898775315e-02, 4.00549811051e-02,-2.68914938808e-02,-2.31141340205e-02, 2.09046452556e-02, 1.12904972786e-02,-1.47018820653e-02,-4.13130665603e-03, 9.21478503219e-03, 5.62571574840e-04,-5.07131450921e-03, 7.16982182106e-04, 2.40069778189e-03,-8.44862666553e-04,-9.42461407722e-04, 5.81075975053e-04, 2.81763925038e-04,-3.03102046072e-04,-4.55568269666e-05, 1.26204335016e-04,-1.15540910383e-05,-4.17514164854e-05, 1.33417614992e-05, 1.03735918404e-05,-6.45673042846e-06,-1.55084435011e-06, 2.14996026993e-06,-8.48708758607e-08,-5.18773373887e-07, 1.39637754550e-07, 8.40035104689e-08,-4.88475793745e-08,-5.42427480028e-09, 1.03470453927e-08,-1.43632948779e-09,-1.34919775398e-09, 5.26113255735e-10, 6.73233649018e-11,-8.27825652253e-11, 1.10169293459e-11, 6.29153731703e-12,-2.48478923756e-12, 2.62649650406e-14, 1.80866123627e-13,-4.24981781957e-14,-4.56339716212e-16, 2.04509967678e-15,-4.40530704248e-16, 4.30459683955e-17,-1.71615245108e-18};


static mus_float_t Battle_Lemarie[24] = {-0.0028284274, -0.004242641, 0.008485282, 0.008485282, -0.018384777, -0.016970564, 0.042426407,
				   0.032526914, -0.11030867, -0.049497478, 0.4341636, 0.7665038, 0.4341636, -0.049497478, -0.11030867,
				   0.032526914, 0.042426407, -0.016970564, -0.018384777, 0.008485282, 0.008485282, -0.004242641, -0.0028284274, 0.0};
static mus_float_t Burt_Adelson[6] = {-0.07071068, 0.3535534, 0.8485282, 0.3535534, -0.07071068, 0.0};
static mus_float_t Beylkin[18] = {0.099305765374353, 0.424215360812961, 0.699825214056600, 0.449718251149468,
			    -.110927598348234, -.264497231446384, 0.026900308803690, 0.155538731877093,
			    -.017520746266529, -.088543630622924, 0.019679866044322, 0.042916387274192,
			    -.017460408696028, -.014365807968852, 0.010040411844631, .0014842347824723,
			    -.002736031626258, .0006404853285212};
static mus_float_t coif2[6] = {0.038580775, -0.12696913, -0.07716155, 0.6074917, 0.74568766, 0.2265843};
static mus_float_t coif4[12] = {0.0011945726958388, -0.01284557955324, 0.024804330519353, 0.050023519962135, -0.15535722285996,
			  -0.071638282295294, 0.57046500145033, 0.75033630585287, 0.28061165190244, -0.0074103835186718,
			  -0.014611552521451, -0.0013587990591632};
static mus_float_t coif6[18] = {-0.0016918510194918, -0.00348787621998426, 0.019191160680044, 0.021671094636352, -0.098507213321468,
			  -0.056997424478478, 0.45678712217269, 0.78931940900416, 0.38055713085151, -0.070438748794943,
			  -0.056514193868065, 0.036409962612716, 0.0087601307091635, -0.011194759273835, -0.0019213354141368,
			  0.0020413809772660, 0.00044583039753204, -0.00021625727664696};
static mus_float_t sym2[5] = {-0.1767767, 0.3535534, 1.0606602, 0.3535534, -0.1767767};
static mus_float_t sym3[4] = {0.1767767, 0.5303301, 0.5303301, 0.1767767};
static mus_float_t sym4[10] = {0.033145633, -0.066291265, -0.1767767, 0.4198447, 0.994369, 0.4198447, -0.1767767, -0.066291265, 0.033145633, 0.0};
static mus_float_t sym5[8] = {0.066291265, -0.19887379, -0.15467963, 0.994369, 0.994369, -0.15467963, -0.19887379, 0.066291265};
static mus_float_t sym6[16] = {-0.0030210863, -0.009063259, -0.016831767, 0.07466399, 0.03133298, -0.30115914, -0.026499243,
			 0.9516422, 0.9516422, -0.026499243, -0.30115914, 0.03133298, 0.07466399, -0.016831767, -0.009063259, -0.0030210863};

static const char *wavelet_names_1[NUM_WAVELETS] =
  {"daub4", "daub6", "daub8", "daub10", "daub12", "daub14", "daub16", "daub18", "daub20",
   "battle-lemarie", "burt-adelson", "beylkin", "coif2", "coif4", "coif6",
   "sym2", "sym3", "sym4", "sym5", "sym6",
   "daub22", "daub24", "daub26", "daub28", "daub30", "daub32", "daub34", "daub36", "daub38", "daub40",
   "daub42", "daub44", "daub46", "daub48", "daub50", "daub52", "daub54", "daub56", "daub58", "daub60",
   "daub62", "daub64", "daub66", "daub68", "daub70", "daub72", "daub74", "daub76"};

static int wavelet_sizes[NUM_WAVELETS] =
  {4, 6, 8, 10, 12, 14, 16, 18, 20,
   24, 6, 18, 6, 12, 18,
   5, 4, 10, 8, 16,
   22, 24, 26, 28, 30, 32, 34, 36, 38, 40, 42, 44, 46, 48, 50, 52, 54, 56, 58, 60, 62, 64, 66, 68, 70, 72, 74, 76};

static mus_float_t *wavelet_data[NUM_WAVELETS] =
  {Daub2, Daub3, Daub4, Daub5, Daub6, Daub7, Daub8, Daub9, Daub10,
   Battle_Lemarie, Burt_Adelson, Beylkin, coif2, coif4, coif6,
   sym2, sym3, sym4, sym5, sym6,
   Daub11, Daub12, Daub13, Daub14, Daub15, Daub16, Daub17, Daub18, Daub19, Daub20,
   Daub21, Daub22, Daub23, Daub24, Daub25, Daub26, Daub27, Daub28, Daub29, Daub30,
   Daub31, Daub32, Daub33, Daub34, Daub35, Daub36, Daub37, Daub38};

const char *wavelet_name(int i) {return(wavelet_names_1[i]);}

const char **wavelet_names(void) {return(wavelet_names_1);}



/* -------------------------------- HAAR TRANSFORM -------------------------------- */
/*
 * from fxt/haar/haar.cc
 */

static void haar_transform(mus_float_t *f, mus_long_t n)
{
  mus_long_t m, i, j, k;
  mus_float_t s2;
  mus_float_t v = 1.0;
  mus_float_t x, y;
  mus_float_t *g;

  s2 = sqrt(0.5);
  g = (mus_float_t *)calloc(n, sizeof(mus_float_t));

  for (m = n; m > 1; m >>= 1)
    {
      mus_long_t mh;
      v *= s2;
      mh = (m >> 1);
      for (j = 0, k = 0; j < m; j += 2, k++)
        {
	  x = f[j];
	  y = f[j + 1];
	  g[k] = x + y;
	  g[mh + k] = (x - y) * v;
        }
      for (i = m - 1; i >= 0; i--) f[i] = g[i];
    }

  f[0] *= v;
  free(g);
}



/* -------------------------------- WALSH TRANSFORM -------------------------------- */

/* borrowed from walsh/walshdit2.cc in the fxt package fxt970929.tgz written by (and copyright) Joerg Arndt
 * arndt@spektracom.de, arndt@jjj.de, http://www.spektracom.de/~arndt/joerg.html, http://www.jjj.de/fxt/
 * fxt.doc appears to say I can use it here (Snd is freeware and I've modified the original to some extent).
 */

static void walsh_transform(mus_float_t *data, mus_long_t n)
{
  int i, ipow;
  ipow = (int)((log(n) / log(2)) + .0001); /* added fudge factor 21-Sep-01 -- (int)3.0000 = 2 on PC */

  for (i = ipow; i >= 1; --i)
    {
      mus_long_t r, m, mh, t1, t2;
      m = (1 << i);
      mh = (m >> 1);
      for (r = 0; r < n; r += m)
        {
	  int j;
	  for (j = 0, t1 = r, t2 = r + mh; j < mh; ++j, ++t1, ++t2)
            {
	      mus_float_t u, v;
	      u = data[t1];
	      v = data[t2];
	      data[t1] = u + v;
	      data[t2] = u - v;
            }
        }
    }
}



static int compare_peaks(const void *pk1, const void *pk2)
{
  if (((fft_peak *)pk1)->freq > ((fft_peak *)pk2)->freq) return(1);
  if (((fft_peak *)pk1)->freq == ((fft_peak *)pk2)->freq) return(0);
  return(-1);
}

int find_and_sort_peaks(mus_float_t *buf, fft_peak *found, int num_peaks, mus_long_t losamp, mus_long_t hisamp)
{
  /* in the fft peak finder below we assume data between 0 and 1 */
  /* this procedure is for the list graph -- see below for fft */

  mus_long_t i, j;
  int pks, minpk;
  mus_float_t minval, ra, ca;
  mus_float_t *peaks;
  mus_long_t *inds;

  if (num_peaks <= 0) return(0);
  peaks = (mus_float_t *)calloc(num_peaks, sizeof(mus_float_t));
  inds = (mus_long_t *)calloc(num_peaks, sizeof(mus_long_t));

  pks = 0;
  ca = 0.0;
  ra = 0.0;
  minval = 0.00001;

  for (i = losamp; i < hisamp; i++)
    {
      mus_float_t la;
      la = ca;
      ca = ra;
      ra = buf[i];
      if ((ca > minval) && (ca > ra) && (ca > la))
	{
	  if (pks < num_peaks)
	    {
	      inds[pks] = i - 1;
	      peaks[pks++] = ca;
	    }
	  else
	    {
	      minval = peaks[0];
	      minpk = 0;
	      for (j = 1; j < num_peaks; j++)
		if (peaks[j] < minval)
		  {
		    minval = peaks[j];
		    minpk = j;
		  }
	      if (ca > minval)
		{
		  inds[minpk] = i - 1;
		  peaks[minpk] = ca;
		}
	    }
	}
    }

  for (i = 0; i < pks; i++)
    {
      j = inds[i];
      found[i].amp = buf[j];
      found[i].freq = (mus_float_t)j;
    }

  if (pks > 0) qsort((void *)found, pks, sizeof(fft_peak), compare_peaks);
  free(peaks);
  free(inds);
  return(pks);
}


#define MIN_CHECK 0.000001

int find_and_sort_transform_peaks(mus_float_t *buf, fft_peak *found, int num_peaks,
				  mus_long_t losamp, mus_long_t hisamp, mus_float_t samps_per_pixel, mus_float_t fft_scale)
{
  /* we want to reflect the graph as displayed, so each "bin" is samps_per_pixel wide */
  mus_long_t i, j, k, hop, pkj;
  int pks, minpk;
  mus_float_t minval, la, ra, ca, logca, logra, logla, offset, fscl, ascl, bscl, fftsize2;
  mus_float_t *peaks;
  mus_long_t *inds;

  fftsize2 = (mus_float_t)hisamp;
  peaks = (mus_float_t *)calloc(num_peaks, sizeof(mus_float_t));
  inds = (mus_long_t *)calloc(num_peaks, sizeof(mus_long_t));

  fscl = 1.0 / fftsize2;
  hop = (mus_long_t)(samps_per_pixel + 0.5);
  if (hop < 1) hop = 1;

  pks = 0;
  /* la = 0.0; */
  ca = 0.0;
  ra = 0.0;
  minval = 0.0;
  ascl = 0.0;
  pkj = 0;

  for (i = losamp; i < hisamp - hop; i += hop)
    {
      mus_long_t oldpkj;
      la = ca;
      ca = ra;
      oldpkj = pkj;
      ra = 0.0;

      for (k = 0; k < hop; k++)
	if (buf[i + k] > ra)
	  {
	    pkj = i + k;
	    ra = buf[pkj]; /* reflect user's view of the graph */
	  }

      if ((ca > minval) && (ca > ra) && (ca > la))
	{
          if (ascl < ca) ascl = ca;
	  if (pks < num_peaks)
	    {
	      inds[pks] = oldpkj;
	      peaks[pks] = ca;
	      pks++;
	    }
	  else
	    {
	      minval = peaks[0];
	      minpk = 0;
	      for (j = 1; j < num_peaks; j++)
		{
		  if (peaks[j] < minval)
		    {
		      minval = peaks[j];
		      minpk = j;
		    }
		}
	      if (ca > minval)
		{
		  inds[minpk] = oldpkj;
		  peaks[minpk] = ca;
		}
	    }
	}
    }

  /* now we have the peaks; turn these into interpolated peaks/amps, and sort */
  if (ascl > 0.0)
    ascl = 1.0 / ascl;
  else ascl = 1.0;

  if (fft_scale > 0.0)
    bscl = fft_scale / ascl;
  else bscl = 1.0;

  for (i = 0, k = 0; i < pks; i++)
    {
      j = inds[i];
      ca = buf[j] * ascl;
      if (j > 0)
	la = buf[j - 1] * ascl;
      else la = ca;
      ra = buf[j + 1] * ascl;

      if ((la < MIN_CHECK) || (ra < MIN_CHECK))
	{
	  found[k].amp = bscl * ca;
	  found[k].freq = fscl * j;
	}
      else
	{
	  logla = log10(la);
	  logca = log10(ca);
	  logra = log10(ra);
	  offset = (0.5 * (logla - logra)) / (logla + logra - 2 * logca); /* this assumes amps < 1.0 (from XJS sms code) */
	  found[k].amp = bscl * pow(10.0, logca - 0.25 * offset * (logla - logra));

	  if ((found[k].amp > 1.0) &&
	      (fft_scale > 0.0))
	    found[k].amp = 1.0;
	  found[k].freq = fscl * (j + offset);
	}
      if (found[k].freq < 0.0) found[k].freq = 0.0;
      if (found[k].amp > 0.0) k++;
    }

  for (i = k; i < num_peaks; i++) found[i].freq = 1.0; /* move blank case to end of sorted list */
  qsort((void *)found, pks, sizeof(fft_peak), compare_peaks);

  free(peaks);
  free(inds);
  return(k);
}


static mus_float_t beta_maxes[MUS_NUM_FFT_WINDOWS] = {1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0,
						1.0, 15.0 /* kaiser */, 10.0, 10.0, 10.0, 1.0, 18.0 /* dolph */, 10.0, 1.0, 18.0 /* samaraki */,
						18.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0,
						1.0, 1.0, 1.0, 1.0, 1.0, 0.24};

mus_float_t fft_beta_max(mus_fft_window_t win) {return(beta_maxes[(int)win]);}



static const char *transform_type_names[NUM_BUILTIN_TRANSFORM_TYPES] = {"Fourier", "Wavelet", "Walsh", "Autocorrelate", "Cepstrum", "Haar"};

static const char *transform_type_program_names[NUM_BUILTIN_TRANSFORM_TYPES] = {
  S_fourier_transform, S_wavelet_transform, S_walsh_transform, S_autocorrelation, S_cepstrum, S_haar_transform};


typedef struct {
  char *name, *xlabel;
  mus_float_t lo, hi;
  Xen proc;
  int type, row, gc_loc;
} added_transform;

static added_transform **added_transforms = NULL;
static int added_transforms_size = 0;


static added_transform *new_transform(void)
{
  int loc = -1;
  if (!added_transforms)
    {
      loc = 0;
      added_transforms_size = 4;
      added_transforms = (added_transform **)calloc(added_transforms_size, sizeof(added_transform *));
    }
  else
    {
      int i;
      for (i = 0; i < added_transforms_size; i++)
	if (!added_transforms[i])
	  {
	    loc = i;
	    break;
	  }
      if (loc == -1)
	{
	  int i;
	  loc = added_transforms_size;
	  added_transforms_size += 4;
	  added_transforms = (added_transform **)realloc(added_transforms, added_transforms_size * sizeof(added_transform *));
	  for (i = loc; i < added_transforms_size; i++) added_transforms[i] = NULL;
	}
    }
  added_transforms[loc] = (added_transform *)calloc(1, sizeof(added_transform));
  added_transforms[loc]->type = loc + NUM_BUILTIN_TRANSFORM_TYPES;
  return(added_transforms[loc]);
}


static int add_transform(const char *name, const char *xlabel, mus_float_t lo, mus_float_t hi, Xen proc)
{
  added_transform *af;
  af = new_transform();
  af->name = mus_strdup(name);
  af->xlabel = mus_strdup(xlabel);
  af->lo = lo;
  af->hi = hi;
  af->proc = proc;
  af->gc_loc = snd_protect(proc);
  make_transform_type_list();
  return(af->type);
}


static added_transform *type_to_transform(int type)
{
  int loc;
  loc = type - NUM_BUILTIN_TRANSFORM_TYPES;
  if ((loc >= 0) && (loc < added_transforms_size))
    return(added_transforms[loc]);
  return(NULL);
}


bool is_transform(int type)
{
  if (type < 0) return(false);
  if (type < NUM_BUILTIN_TRANSFORM_TYPES) return(true);
  return(type_to_transform(type));
}


/* delete-transform would also need to remove its name from the various UI lists */

const char *transform_name(int type)
{
  if (is_transform(type))
    {
      added_transform *af;
      if (type < NUM_BUILTIN_TRANSFORM_TYPES)
	return(transform_type_names[type]);
      af = type_to_transform(type);
      return(af->name);
    }
  return("unknown");
}


const char *transform_program_name(int type)
{
  if (is_transform(type))
    {
      added_transform *af;
      if (type < NUM_BUILTIN_TRANSFORM_TYPES)
	return(transform_type_program_names[type]);
      af = type_to_transform(type);
      return(af->name);
    }
  return("unknown");
}


static const char *added_transform_xlabel(int type)
{
  added_transform *af;
  af = type_to_transform(type);
  if (af) return(af->xlabel);
  return("unknown");
}


static mus_float_t added_transform_lo(int type)
{
  added_transform *af;
  af = type_to_transform(type);
  if (af) return(af->lo);
  return(0.0);
}


static mus_float_t added_transform_hi(int type)
{
  added_transform *af;
  af = type_to_transform(type);
  if (af) return(af->hi);
  return(1.0);
}


static Xen added_transform_proc(int type)
{
  added_transform *af;
  af = type_to_transform(type);
  if (af) return(af->proc);
  return(Xen_false);
}


void set_transform_position(int i, int j)
{
  if (i >= NUM_BUILTIN_TRANSFORM_TYPES)
    {
      added_transform *af;
      af = type_to_transform(i);
      af->row = j;
    }
}


int max_transform_type(void)
{
  return(NUM_BUILTIN_TRANSFORM_TYPES + added_transforms_size);
}


int transform_type_to_position(int type)
{
  added_transform *af;
  if (type < NUM_BUILTIN_TRANSFORM_TYPES)
    return(type);
  af = type_to_transform(type);
  if (af) return(af->row);
  return(-1);
}


int transform_position_to_type(int pos)
{
  int i;
  if (pos < NUM_BUILTIN_TRANSFORM_TYPES)
    return(pos);
  for (i = 0; i < added_transforms_size; i++)
    {
      added_transform *af;
      af = added_transforms[i];
      if ((af) && (af->row = pos))
	return(af->type);
    }
  return(-1);
}


static const char *spectro_xlabel(chan_info *cp)
{
  switch (cp->transform_type)
    {
    case FOURIER:
      if (cp->fft_log_frequency)
	return("log freq");
      else return("frequency");

    case WAVELET:         return(wavelet_names_1[cp->wavelet_type]);
    case HAAR:            return("Haar spectrum");
    case CEPSTRUM:        return("cepstrum");
    case WALSH:           return("Sequency");
    case AUTOCORRELATION: return("Lag time");
    default:              return(added_transform_xlabel(cp->transform_type));
    }
  return(NULL);
}


static void make_sonogram_axes(chan_info *cp)
{
  fft_info *fp;
  fp = cp->fft;
  if (fp)
    {
      axis_info *ap;
      mus_float_t max_freq, min_freq, yang;
      ap = cp->axis;
      if (cp->transform_type == FOURIER)
	{
	  max_freq = cp->spectrum_end * (mus_float_t)snd_srate(cp->sound) * 0.5;
	  if ((cp->fft_log_frequency) && ((snd_srate(cp->sound) * 0.5 * cp->spectrum_start) < log_freq_start(ss)))
	    min_freq = log_freq_start(ss);
	  else min_freq = cp->spectrum_start * (mus_float_t)snd_srate(cp->sound) * 0.5;
	}
      else
	{
          if ((cp->transform_type == AUTOCORRELATION) ||
              (cp->transform_type == CEPSTRUM))
	    {
	      max_freq = fp->current_size * cp->spectrum_end / 2;
	      min_freq = fp->current_size * cp->spectrum_start / 2;
	    }
	  else
	    {
	      max_freq = fp->current_size * cp->spectrum_end;
	      min_freq = fp->current_size * cp->spectrum_start;
	    }
	}
      yang = fmod(cp->spectro_y_angle, 360.0);
      if (yang < 0.0) yang += 360.0;
      if (cp->transform_graph_type == GRAPH_AS_SPECTROGRAM)
	{
	  const char *xlabel;
	  if (cp->transform_type == FOURIER)
	    {
	      if (yang < 45.0) xlabel = "frequency";
	      else if (yang < 135.0) xlabel = "time";
	      else if (yang < 225.0) xlabel = "wavelength?";
	      else if (yang < 315.0) xlabel = "reversed time?";
	      else xlabel = "frequency";
	    }
	  else xlabel = spectro_xlabel(cp);
	  fp->axis = make_axis_info(cp,
				    min_freq, max_freq,
				    ap->x0, ap->x1,
				    xlabel,
				    min_freq, max_freq,
				    ap->x0, ap->x1,
				    fp->axis);
	}
      else
	{
	  if (!fp->xlabel)
	    fp->xlabel = mus_strdup("time");
	  fp->axis = make_axis_info(cp,
				    ap->x0, ap->x1,
				    min_freq, max_freq,
				    fp->xlabel,
				    ap->x0, ap->x1,
				    min_freq, max_freq,
				    fp->axis);
	}
    }
}


typedef struct fft_state {
  mus_long_t size;
  mus_fft_window_t wintype;
  graph_type_t graph_type;
  bool done;
  chan_info *cp;
  mus_float_t *window;
  mus_float_t *data;
  mus_float_t alpha, beta;
  int wavelet_choice, transform_type;
  mus_long_t beg, databeg, datalen;
  mus_long_t losamp;
  int edit_ctr;
  bool dBing, lfreq, with_phases;
  int pad_zero;
  mus_float_t cutoff;
} fft_state;


static mus_float_t *fs_idata = NULL;
static mus_long_t fs_idata_size = 0;

void fourier_spectrum(snd_fd *sf, mus_float_t *fft_data, mus_long_t fft_size, mus_long_t data_len, mus_float_t *window, chan_info *cp)
{
  mus_long_t i, j, lim;

  if (window)
    {
      for (i = 0; i < data_len; i++)
	fft_data[i] = window[i] * read_sample(sf);
    }
  else
    {
      for (i = 0; i < data_len; i++)
	fft_data[i] = read_sample(sf);
    }

  if (data_len < fft_size)
    mus_clear_floats(fft_data + data_len, fft_size - data_len);
  if (fft_size <= fs_idata_size)
    mus_clear_floats(fs_idata, fft_size);
  else
    {
      if (!fs_idata)
	fs_idata = (mus_float_t *)malloc(fft_size * sizeof(mus_float_t));
      else fs_idata = (mus_float_t *)realloc(fs_idata, fft_size * sizeof(mus_float_t));
      mus_clear_floats(fs_idata, fft_size);
      fs_idata_size = fft_size;
    }

  mus_fft(fft_data, fs_idata, fft_size, 1);

  lim = fft_size / 2;
  if ((cp) && (cp->fft_with_phases))
    {
      fft_data[0] = hypot(fft_data[0], fs_idata[0]);
      cp->fft->phases[0] = -atan2(fs_idata[0], fft_data[0]);
      for (i = 1, j = fft_size - 1; i < lim; i++, j--)
	{
	  fft_data[i] = hypot(fft_data[i], fs_idata[i]);
	  fft_data[j] = fft_data[i];
	  cp->fft->phases[i] = -atan2(fs_idata[i], fft_data[i]);
	  cp->fft->phases[j] = -cp->fft->phases[i];
	}
    }
  else
    {
      fft_data[0] = hypot(fft_data[0], fs_idata[0]);
      for (i = 1, j = fft_size - 1; i < lim; i++, j--)
	{
	  fft_data[i] = hypot(fft_data[i], fs_idata[i]);
	  fft_data[j] = fft_data[i];
	}
    }
  if (fs_idata_size >= 4194304)
    {
      free(fs_idata);
      fs_idata = NULL;
      fs_idata_size = 0;
    }
}


static Xen before_transform_hook;

static void apply_fft(fft_state *fs)
{
  mus_long_t i, ind0, data_len;
  mus_float_t *fft_data;
  snd_fd *sf;
  chan_info *cp;

  cp = fs->cp;
  fft_data = fs->data;
  data_len = cp->transform_size;

  if ((show_selection_transform(ss)) &&
      (selection_is_active_in_channel(cp)) &&
      (fs->datalen > 0))
    {
      ind0 = fs->databeg;
      if (cp->transform_graph_type == GRAPH_ONCE)
	data_len = fs->datalen;
    }
  else
    {
      if (Xen_hook_has_list(before_transform_hook))
	{
	  Xen res;
	  res = run_progn_hook(before_transform_hook,
			       Xen_list_2(C_int_to_Xen_sound(cp->sound->index),
					  C_int_to_Xen_integer(cp->chan)),
			       S_before_transform_hook);
	  if (Xen_is_llong(res))
	    ind0 = Xen_llong_to_C_llong(res) + fs->beg;
	  else ind0 = cp->axis->losamp + fs->beg;
	}
      else
	ind0 = cp->axis->losamp + fs->beg;
    }

  if (data_len > fs->size) data_len = fs->size;
  sf = init_sample_read(ind0, cp, READ_FORWARD);
  if (!sf) return;

  switch (cp->transform_type)
    {
    case FOURIER:
      fourier_spectrum(sf, fft_data, fs->size, data_len, fs->window, cp);
      break;

    case WAVELET:
      for (i = 0; i < data_len; i++) fft_data[i] = read_sample(sf);
      if (data_len < fs->size)
	mus_clear_floats(fft_data + data_len, fs->size - data_len);
      wavelet_transform(fft_data, fs->size, wavelet_data[cp->wavelet_type], wavelet_sizes[cp->wavelet_type]);
      break;

    case HAAR:
      for (i = 0; i < data_len; i++) fft_data[i] = read_sample(sf);
      if (data_len < fs->size)
	mus_clear_floats(fft_data + data_len, fs->size - data_len);
      haar_transform(fft_data, fs->size);
      break;

    case CEPSTRUM:
      for (i = 0; i < data_len; i++) fft_data[i] = read_sample(sf);
      if (data_len < fs->size)
	mus_clear_floats(fft_data + data_len, fs->size - data_len);
      mus_cepstrum(fft_data, fs->size);
      break;

    case WALSH:
      for (i = 0; i < data_len; i++) fft_data[i] = read_sample(sf);
      if (data_len < fs->size)
	mus_clear_floats(fft_data + data_len, fs->size - data_len);
      walsh_transform(fft_data, fs->size);
      break;

    case AUTOCORRELATION:
      for (i = 0; i < data_len; i++) fft_data[i] = read_sample(sf);
      if (data_len < fs->size)
	mus_clear_floats(fft_data + data_len, fs->size - data_len);
      mus_autocorrelate(fft_data, fs->size);
      break;

    default:
      {
	Xen res, sfd;
	int gc_loc, sf_loc;
	sfd = g_c_make_sampler(sf);
	sf_loc = snd_protect(sfd);
	res = Xen_call_with_2_args(added_transform_proc(cp->transform_type),
			 C_llong_to_Xen_llong(data_len),
			 sfd,
			 "added transform func");
	gc_loc = snd_protect(res);
	if (mus_is_vct(res))
	  {
	    vct *v;
	    mus_long_t len;
	    v = Xen_to_vct(res);
	    len = mus_vct_length(v);
	    mus_copy_floats(fft_data, mus_vct_data(v), len);
	  }
	snd_unprotect_at(gc_loc);
	snd_unprotect_at(sf_loc);
	free_snd_fd_almost(sf);
	return;
      }
    }

  free_snd_fd(sf);
}


static void display_fft(fft_state *fs)
{
  fft_info *fp;
  chan_info *cp;
  mus_float_t max_freq = 0.0, min_freq = 0.0, max_val, min_val, data_max = 0.0, scale = 1.0;
  const char *xlabel;
  fft_info *nfp;
  mus_float_t *data, *tdata;
  chan_info *ncp;
  snd_info *sp;
  mus_long_t i, lo, hi;

  cp = fs->cp;
  if ((!cp) || (cp->active < CHANNEL_HAS_AXES)) return;
  if (cp->transform_graph_type != GRAPH_ONCE) return;

  fp = cp->fft;
  if (!fp) return; /* can happen if selection transform set, but no selection */

  data = fp->data;
  if (!data) return;

  sp = cp->sound;
  if (!fp->xlabel)
    xlabel = spectro_xlabel(cp);
  else xlabel = fp->xlabel;
  /* this only works until the fft data is remade in some way, then the label reverts to "frequency" */

  switch (cp->transform_type)
    {
    case FOURIER:
      max_freq = ((mus_float_t)(snd_srate(sp)) * 0.5 * cp->spectrum_end);
      if ((cp->fft_log_frequency) && ((snd_srate(sp) * 0.5 * cp->spectrum_start) < log_freq_start(ss)))
	min_freq = log_freq_start(ss);
      else min_freq = ((mus_float_t)(snd_srate(sp)) * 0.5 * cp->spectrum_start);
      break;

    case WAVELET: case WALSH: case HAAR:
      max_freq = fs->size * cp->spectrum_end;
      min_freq = fs->size * cp->spectrum_start;
      break;

    case AUTOCORRELATION: case CEPSTRUM:
      max_freq = fs->size * cp->spectrum_end / 2;
      min_freq = fs->size * cp->spectrum_start / 2;
      break;

    default:
      min_freq = added_transform_lo(cp->transform_type) * fs->size * cp->spectrum_end;
      max_freq = added_transform_hi(cp->transform_type) * fs->size * cp->spectrum_end;
      break;
    }

  if (cp->transform_normalization == DONT_NORMALIZE)
    {
      lo = 0;
      hi = (int)(fp->current_size / 2);
    }
  else
    {
      if (cp->transform_type == FOURIER)
	{
	  hi = (int)(fs->size * cp->spectrum_end / 2);
	  lo = (int)(fs->size * cp->spectrum_start / 2);
	}
      else
	{
	  hi = (int)(fs->size * cp->spectrum_end);
	  lo = (int)(fs->size * cp->spectrum_start);
	}
    }

  data_max = 0.0;
  if ((cp->transform_normalization == NORMALIZE_BY_SOUND) ||
      ((cp->transform_normalization == DONT_NORMALIZE) &&
       (sp->nchans > 1) &&
       (sp->channel_style == CHANNELS_SUPERIMPOSED)))
    {
      uint32_t j;
      for (j = 0; j < sp->nchans; j++)
	{
	  ncp = sp->chans[j];
	  if ((ncp->graph_transform_on) && (ncp->fft)) /* normalize-by-sound but not ffting all chans? */
	    {
	      nfp = ncp->fft;
	      tdata = nfp->data;
	      for (i = lo; i < hi; i++)
		if (tdata[i] > data_max)
		  data_max = tdata[i];
	    }
	}
    }
  else
    {
      if (cp->transform_type == FOURIER)
	{
	  for (i = lo; i < hi; i++)
	    if (data[i] > data_max)
	      data_max = data[i];
	}
      else
	{
	  for (i = lo; i < hi; i++)
	    {
	      if (data[i] > data_max)
		data_max = data[i];
	      else
		if (data[i] < -data_max)
		  data_max = -data[i];
	    }
	}
    }

  if (data_max == 0.0) data_max = 1.0;
  if (cp->transform_normalization != DONT_NORMALIZE)
    scale = 1.0 / data_max;
  else
    {
      if (cp->transform_type == FOURIER)
	{
	  int di;
	  scale = 2.0 / (mus_float_t)(fs->size);
	  di = (int)(10 * data_max * scale + 1);
	  if (di == 1)
	    {
	      di = (int)(100 * data_max * scale + 1);
	      if (di == 1)
		{
		  di = (int)(1000 * data_max * scale + 1);
		  data_max = (mus_float_t)di / 1000.0;
		}
	      else data_max = (mus_float_t)di / 100.0;
	    }
	  else data_max = (mus_float_t)di / 10.0;
	}
      else
	{
	  scale = 1.0;
	}
    }

  if (cp->fft_log_magnitude)
    {
      /* the DONT_NORMALIZE option is ignored because it really looks dumb in this context */
      max_val = 0.0;
      min_val = cp->min_dB;
    }
  else
    {
      if (cp->transform_normalization == DONT_NORMALIZE)
	{
	  if (cp->transform_type == FOURIER)
	    min_val = 0.0;
	  else min_val = -data_max;
	  max_val = data_max;
	}
      else
	{
	  if (cp->transform_type == FOURIER)
	    min_val = 0.0;
	  else min_val = -1.0;
	  max_val = 1.0;
	}
    }

  fp->scale = scale;
  fp->axis = make_axis_info(cp,
			    min_freq, max_freq,
			    min_val, max_val,
			    xlabel,
			    min_freq, max_freq,
			    min_val, max_val,
			    fp->axis);
}


static fft_state *free_fft_state(fft_state *fs)
{
  if (fs)
    {
      if (fs->window) {free(fs->window); fs->window = NULL;}
      free(fs);
    }
  return(NULL);
}


void cp_free_fft_state(chan_info *cp)
{
  if (cp->fft_data)
    cp->fft_data = (fft_state *)free_fft_state(cp->fft_data);
}


bool fft_window_beta_in_use(mus_fft_window_t win)
{
  return((win == MUS_KAISER_WINDOW) ||
	 (win == MUS_CAUCHY_WINDOW) ||
	 (win == MUS_POISSON_WINDOW) ||
	 (win == MUS_GAUSSIAN_WINDOW) ||
	 (win == MUS_TUKEY_WINDOW) ||
	 (win == MUS_DOLPH_CHEBYSHEV_WINDOW) ||
	 (win == MUS_HANN_POISSON_WINDOW) ||
	 (win == MUS_SAMARAKI_WINDOW) ||
	 (win == MUS_ULTRASPHERICAL_WINDOW) ||
	 (win == MUS_DPSS_WINDOW));
}


bool fft_window_alpha_in_use(mus_fft_window_t win)
{
  return(win == MUS_ULTRASPHERICAL_WINDOW);
}


static fft_state *make_fft_state(chan_info *cp, bool force_recalc)
{
  fft_state *fs = NULL;
  axis_info *ap;
  bool reuse_old = false;
  mus_long_t fftsize;
  mus_long_t dbeg = 0, dlen = 0;

  ap = cp->axis;

  if ((show_selection_transform(ss)) &&
      (cp->transform_graph_type == GRAPH_ONCE) &&
      (selection_is_active_in_channel(cp)))
    {
      /* override transform_size(ss) in this case (sonograms cover selection but use preset size) */
      dbeg = selection_beg(cp);
      dlen = selection_len();
      /* these need to be handled at the same time, and not re-examined until the next call */
      /* if we're sweeping the mouse defining the selection, by the time we get to apply_fft, selection_len() can change */
      fftsize = snd_mus_long_t_pow2((int)ceil(log(dlen * (1 + cp->zero_pad) + 1) / log(2.0)));
      if (fftsize < 2) fftsize = 2;
      cp->selection_transform_size = fftsize;
    }
  else
    {
      if ((cp->zero_pad == 0) && (is_power_of_2(cp->transform_size)))
	fftsize = cp->transform_size;
      else fftsize = snd_mus_long_t_pow2((int)(ceil(log((mus_float_t)(cp->transform_size * (1 + cp->zero_pad))) / log(2.0))));
      if (fftsize < 2) fftsize = 2;
      cp->selection_transform_size = 0;
    }

  if ((!force_recalc) && (cp->fft_data) && (cp->selection_transform_size == 0))
    {
      fs = cp->fft_data;
      if ((fs->losamp == ap->losamp) &&
	  (!(Xen_hook_has_list(before_transform_hook))) &&
	  (fs->size == fftsize) &&
	  (fs->transform_type == cp->transform_type) &&
	  (fs->wintype == cp->fft_window) &&
	  ((!(fft_window_alpha_in_use(fs->wintype))) || (fs->alpha == cp->fft_window_alpha)) &&
	  ((!(fft_window_beta_in_use(fs->wintype))) || (fs->beta == cp->fft_window_beta)) &&
	  (fs->dBing == cp->fft_log_magnitude) &&
	  (fs->lfreq == cp->fft_log_frequency) &&
	  (fs->with_phases == cp->fft_with_phases) &&
	  (fs->pad_zero == cp->zero_pad) &&
	  (fs->cutoff == cp->spectrum_end) &&
	  (fs->graph_type == cp->transform_graph_type) &&
	  (fs->wavelet_choice == cp->wavelet_type) &&
	  (fs->edit_ctr == cp->edit_ctr))
	reuse_old = true;
    }

  if (!reuse_old)
    {
      cp_free_fft_state(cp);
      fs = (fft_state *)calloc(1, sizeof(fft_state));
      fs->cp = cp;
      fs->cutoff = cp->spectrum_end;
      fs->size = fftsize;
      fs->pad_zero = cp->zero_pad;
      fs->wintype = cp->fft_window;
      fs->dBing = cp->fft_log_magnitude;
      fs->lfreq = cp->fft_log_frequency;
      fs->with_phases = cp->fft_with_phases;
      fs->window = NULL;
      fs->losamp = ap->losamp;
      fs->edit_ctr = cp->edit_ctr;
      fs->wavelet_choice = cp->wavelet_type;
      fs->transform_type = cp->transform_type;
      fs->graph_type = cp->transform_graph_type;
      fs->alpha = cp->fft_window_alpha;
      fs->beta = cp->fft_window_beta;
    }

  fs->done = reuse_old;
  fs->beg = 0; /* offset from losamp */
  fs->databeg = dbeg;
  fs->datalen = dlen;
  return(fs);
}


static fft_info *make_fft_info(mus_long_t size, mus_fft_window_t window, mus_float_t alpha, mus_float_t beta)
{
  fft_info *fp;
  fp = (fft_info *)calloc(1, sizeof(fft_info));
  fp->size = size;
  fp->window = window;
  fp->alpha = alpha;
  fp->beta = beta;
  fp->xlabel = NULL;
  fp->data = (mus_float_t *)calloc(size + 1, sizeof(mus_float_t)); /*  + 1 for complex storage or starts at 1 or something */
  fp->phases = NULL;
  return(fp);
}


void set_fft_info_xlabel(chan_info *cp, const char *new_label)
{
  if ((cp) && (cp->fft))
    {
      if (cp->fft->xlabel) free(cp->fft->xlabel);
      if (new_label)
	cp->fft->xlabel = mus_strdup(new_label);
      else cp->fft->xlabel = NULL;
    }
}


fft_info *free_fft_info(fft_info *fp)
{
  fp->chan = NULL;
  if (fp->data) free(fp->data);
  if (fp->phases) free(fp->phases);
  if (fp->axis) free_axis_info(fp->axis);
  if (fp->xlabel) free(fp->xlabel);
  free(fp);
  return(NULL);
}


static void one_fft(fft_state *fs)
{
  if (!fs->done)
    {
      /* allocate arrays if needed */
      fft_info *fp;
      chan_info *cp;
      cp = fs->cp;
      fp = cp->fft;
      if (!fp)                              /* associated channel hasn't done any ffts yet, so there's no struct */
	{
	  cp->fft = make_fft_info(fs->size, fs->wintype, fs->alpha, fs->beta);
	  fp = cp->fft;
	}
      else
	{
	  if ((!fp->data) ||
	      (fs->size > fp->size))
	    {
	      fp->size = fs->size;
	      if (fp->data) free(fp->data);
	      fp->data = (mus_float_t *)calloc(fp->size + 1, sizeof(mus_float_t));
	      if (fp->phases) free(fp->phases);
	      fp->phases = NULL;
	    }
	}

      if ((cp->fft_with_phases) &&
	  (!(fp->phases)))
	fp->phases = (mus_float_t *)calloc(fp->size + 1, sizeof(mus_float_t));

      fp->current_size = fs->size; /* protect against parallel size change via fft size menu */
      fs->data = fp->data;
      if (!fs->window)
	{
	  static mus_long_t last_size = 0;
	  static int last_zero = 0;
	  static mus_fft_window_t last_wintype = MUS_RECTANGULAR_WINDOW;
	  static mus_float_t last_beta = 0.0;
	  static mus_float_t last_alpha = 0.0;
	  static mus_float_t *last_window = NULL;

	  fs->window = (mus_float_t *)calloc(fs->size, sizeof(mus_float_t));

	  if ((fs->wintype != last_wintype) ||
	      (fs->size != last_size) ||
	      (fs->beta != last_beta) ||
	      (fs->alpha != last_alpha) ||
	      (fs->pad_zero != last_zero))
	    {
	      if (last_window) free(last_window);
	      last_window = (mus_float_t *)calloc(fs->size, sizeof(mus_float_t));
	      if (cp->selection_transform_size > 0)
		mus_make_fft_window_with_window(fs->wintype, cp->selection_transform_size,
						fs->beta * beta_maxes[fs->wintype],
						fs->alpha, last_window);
	      else mus_make_fft_window_with_window(fs->wintype, cp->transform_size,
						   fs->beta * beta_maxes[fs->wintype],
						   fs->alpha, last_window);
	      last_size = fs->size;
	      last_alpha = fs->alpha;
	      last_beta = fs->beta;
	      last_wintype = fs->wintype;
	      last_zero = fs->pad_zero;
	    }
	  mus_copy_floats(fs->window, last_window, fs->size);
	}
      apply_fft(fs);
    }
  display_fft(fs);
}


void single_fft(chan_info *cp, bool update_display, bool force_recalc)
{
  if (cp->transform_size < 2) return;
  cp->fft_data = make_fft_state(cp, force_recalc);
  one_fft(cp->fft_data);
  if (update_display) display_channel_fft_data(cp);
}



/* -------------------------------- GRAPH_AS_SONOGRAM -------------------------------- */
/*
 * calls calculate_fft for each slice, each bin being a gray-scaled rectangle in the display
 */

/* as we run the ffts, we need to save the fft data for subsequent redisplay/printing etc */
/* many of these can be running in parallel so the pointers can't be global */

/* this work proc calls a loop by pixels (hop depends on pixel/samps decision)
   each pixel(group) sets up the fft_state pointer with beg reflecting hop
   then loops, each time called, calling fft_in_slices until it returns true.
   then grab that data, update the channel display, look to see if we're
   behind the times, if so cleanup and exit, else jump back to outer loop.
   */

typedef enum {SONO_INIT, SONO_RUN, SONO_QUIT, SONO_DONE} sono_slice_t;

typedef struct sonogram_state {
  sono_slice_t slice;
  int outlim, outer;
  fft_state *fs;
  chan_info *cp;
  int spectrum_size;
  sono_info *scp;
  mus_long_t beg, losamp, hisamp;
  bool done;
  double acc_hop, hop;
  mus_fft_window_t window;
  int msg_ctr;
  int edit_ctr;
  mus_float_t old_scale, alpha, beta;
  graph_type_t graph_type;
  int transform_type, w_choice;
  bool old_logxing;
  bool status_needs_to_be_cleared;
  bool force_recalc;
} sonogram_state;


void clear_transform_edit_ctrs(chan_info *cp)
{
  if (cp->fft_data)
    {
      fft_state *fs;
      fs = cp->fft_data;
      fs->edit_ctr = -1;
    }
  if (cp->last_sonogram)
    {
      sonogram_state *lsg;
      lsg = (sonogram_state *)(cp->last_sonogram);
      lsg->edit_ctr = -1;
    }
}


void *make_sonogram_state(chan_info *cp, bool force_recalc)
{
  sonogram_state *sg;
  fft_state *fs;
  sg = (sonogram_state *)calloc(1, sizeof(sonogram_state));
  sg->cp = cp;
  sg->done = false;
  sg->force_recalc = force_recalc;
  fs = make_fft_state(cp, true);
  cp->fft_data = NULL;
  sg->fs = fs;
  sg->msg_ctr = 8;
  sg->transform_type = cp->transform_type;
  sg->w_choice = cp->wavelet_type;
  sg->status_needs_to_be_cleared = false;
  if (cp->temp_sonogram)
    {
      sonogram_state *temp_sg;
      /* we must have restarted fft process without letting the previous run at all */
      temp_sg = (sonogram_state *)(cp->temp_sonogram);
      if (temp_sg->fs) temp_sg->fs = free_fft_state(temp_sg->fs);
      free(temp_sg);
      /* cp->last_sonogram = NULL; */
    }
  cp->temp_sonogram = sg; /* background process may never run, so we need a way to find this pointer at cleanup time */
  return((void *)sg);
}


void free_sonogram_fft_state(void *ptr)
{
  sonogram_state *sg = (sonogram_state *)ptr;
  if (sg->fs) free_fft_state(sg->fs);
  sg->fs = NULL;
}


void free_sono_info(chan_info *cp)
{
  sono_info *si;
  si = cp->sonogram_data;
  if (si)
    {
      if (si->begs) free(si->begs);
      if (si->data)
	{
	  int i;
	  for (i = 0; i < si->total_slices; i++)
	    if (si->data[i])
	      free(si->data[i]);
	  free(si->data);
	}
      free(si);
      cp->sonogram_data = NULL;
    }
}


static bool memory_is_available(mus_long_t slices, mus_long_t bins)
{
  /* as far as I can tell, the approved way to make sure we can allocate enough memory is to allocate it... */
  mus_long_t bytes_needed; /* "mus_long_t" throughout is vital here */

  bytes_needed = (mus_long_t)(slices * (sizeof(mus_long_t) + sizeof(mus_float_t *))) +  /* begs + pointer to each slice data */
                 (mus_long_t)(bins * slices * sizeof(mus_float_t)) +               /* data for all the slices */
                 (mus_long_t)(bins * 2 * sizeof(double));                    /* FFT space */

  /* apparently Linux always returns a pointer, even when it will fail later, so we have to
   *   touch the memory to find out if it's there; but this works only if we then try another alloc!
   *   this is ridiculous...
   */

  if (bytes_needed > 100000000) /* assume 100 MBytes is ok */
    {
      char *check_alloc[10];
      int i;
      mus_long_t bytes_needed_10;

      bytes_needed_10 = bytes_needed / 10;
      for (i = 0; i < 10; i++)
	{
	  check_alloc[i] = (char *)malloc(bytes_needed_10);
	  if (!check_alloc[i])
	    {
	      int j;
	      snd_warning("can't allocate enough memory to run this set of FFTS: %" print_mus_long " bytes needed", bytes_needed);
	      for (j = 0; j < i; j++)
		free(check_alloc[j]);
	      return(false);
	    }
	  check_alloc[i][0] = ' '; /* touch it! */
	}

      for (i = 0; i < 10; i++)
	free(check_alloc[i]);
    }
  return(true);
}


static sono_slice_t set_up_sonogram(sonogram_state *sg)
{
  /* return SONO_RUN to go on, SONO_QUIT to quit early */
  sono_info *si;
  axis_info *ap;
  chan_info *cp;
  sonogram_state *lsg = NULL;
  int i, dpys = 1;

  cp = sg->cp;
  if (cp->fft_changed != FFT_CHANGE_LOCKED)
    cp->fft_changed = FFT_UNCHANGED;
  else cp->fft_changed = FFT_CHANGED;

  if ((!(cp->graph_transform_on)) || (cp->transform_size <= 2)) return(SONO_QUIT);
  ap = cp->axis;
  sg->slice = SONO_INIT;
  sg->outer = 0;
  sg->beg = ap->losamp;
  sg->losamp = ap->losamp;
  sg->hisamp = ap->hisamp;
  sg->window = cp->fft_window;
  sg->alpha = cp->fft_window_alpha;
  sg->beta = cp->fft_window_beta;

  if (cp->graph_time_on) dpys++;
  if (cp->graph_lisp_on) dpys++;

  if (cp->transform_graph_type == GRAPH_AS_SPECTROGRAM)
    sg->outlim = ap->height / cp->spectro_hop;
  else sg->outlim = ap->window_width / dpys;
  if (sg->outlim <= 1) return(SONO_QUIT);

  /* sg->hop = (int)(ceil((double)(ap->hisamp - ap->losamp + 1) / (double)(sg->outlim))); */
  /*   this ^ old form causes a noticeable step in the zoom sequence when the hop value changes */
  sg->hop = (double)(ap->hisamp - ap->losamp + 1) / (double)(sg->outlim);
  sg->acc_hop = 0.0;

  /* if fewer samps than pixels, draw rectangles */
  if ((cp->transform_type == FOURIER) ||
      (cp->transform_type == AUTOCORRELATION) ||
      (cp->transform_type == CEPSTRUM))
    sg->spectrum_size = sg->fs->size / 2;
  else sg->spectrum_size = sg->fs->size;
   if (sg->spectrum_size <= 0) return(SONO_QUIT);

  sg->edit_ctr = cp->edit_ctr;
  si = cp->sonogram_data;

  if (!si)
    {
      si = (sono_info *)calloc(1, sizeof(sono_info));
      si->total_bins = sg->spectrum_size;
      si->total_slices = snd_to_int_pow2(sg->outlim);

      if (!memory_is_available((mus_long_t)(si->total_slices), (mus_long_t)(si->total_bins)))
	{
	  free(si);
	  return(SONO_QUIT);
	}

      cp->sonogram_data = si;
      si->begs = (mus_long_t *)calloc(si->total_slices, sizeof(mus_long_t));
      si->data = (mus_float_t **)calloc(si->total_slices, sizeof(mus_float_t *));
      for (i = 0; i < si->total_slices; i++) si->data[i] = (mus_float_t *)calloc(si->total_bins, sizeof(mus_float_t));
    }
  else
    if ((si->total_slices < sg->outlim) ||
	(si->total_bins < sg->spectrum_size))
      {
	int tempsize;

	tempsize = snd_to_int_pow2(sg->outlim);
	if (!memory_is_available((mus_long_t)tempsize, (mus_long_t)(sg->spectrum_size)))
	  return(SONO_QUIT);

	for (i = 0; i < si->total_slices; i++)
	  if (si->data[i])
	    {
	      free(si->data[i]);
	      si->data[i] = NULL;
	    }

	if (si->total_slices < tempsize)
	  {
	    si->total_slices = tempsize;
	    free(si->data);
	    si->begs = (mus_long_t *)realloc(si->begs, si->total_slices * sizeof(mus_long_t));
	    si->data = (mus_float_t **)calloc(si->total_slices, sizeof(mus_float_t *));
	  }

	if (si->total_bins < sg->spectrum_size) si->total_bins = sg->spectrum_size;
	for (i = 0; i < si->total_slices; i++) si->data[i] = (mus_float_t *)calloc(si->total_bins, sizeof(mus_float_t));
      }
  sg->scp = si;
  si->target_bins = sg->spectrum_size;
  si->active_slices = 0;
  si->target_slices = sg->outlim;
  si->scale = 0.0;
  if ((!(sg->force_recalc)) && (cp->last_sonogram))      /* there was a previous run */
    {
      lsg = (sonogram_state *)(cp->last_sonogram);
      if ((lsg->done) &&                                 /* it completed all ffts */
	  (lsg->outlim == sg->outlim) &&                 /* the number of ffts is the same */
	  (lsg->spectrum_size == sg->spectrum_size) &&   /* ditto fft sizes [fs->size] */
	  (lsg->losamp == sg->losamp) &&                 /* begins are same */
	  (lsg->hisamp == sg->hisamp) &&                 /* ends are same */
	  (lsg->window == sg->window) &&                 /* data windows are same [fs->wintype, cp->fft_window] */
	  (lsg->transform_type == sg->transform_type) && /* transform types are the same */
	  (lsg->w_choice == sg->w_choice) &&             /* wavelets are the same [cp->wavelet_type] */
	  (lsg->edit_ctr == sg->edit_ctr) &&             /* underlying data is the same */
	  ((!(fft_window_alpha_in_use(sg->window))) || (lsg->alpha == sg->alpha)) &&
	  ((!(fft_window_beta_in_use(sg->window))) || (lsg->beta == sg->beta)))
	{
	  sg->outer = sg->outlim;                        /* fake up the run */
	  si->active_slices = si->target_slices;
	  sg->old_scale = lsg->old_scale;
	  si->scale = sg->old_scale;
	  if ((lsg->graph_type != cp->transform_graph_type) ||
	      (lsg->old_logxing != cp->fft_log_frequency))
	    make_sonogram_axes(cp);                      /* may need to fixup frequency axis labels */
	  sg->graph_type = cp->transform_graph_type;
	  sg->old_logxing = cp->fft_log_frequency;
	  return(SONO_QUIT);                             /* so skip the ffts! */
	}
    }

  cp->fft_changed = FFT_CHANGED;
  start_progress_report(cp);
  sg->status_needs_to_be_cleared = true;
  return(SONO_RUN);
}


static sono_slice_t run_all_ffts(sonogram_state *sg)
{
  fft_state *fs;
  sono_info *si;
  chan_info *cp;
  axis_info *ap;
  /* check for losamp/hisamp change? */

  one_fft((fft_state *)(sg->fs));

  fs = sg->fs;
  cp = sg->cp;
  si = cp->sonogram_data;
  if (si->active_slices < si->total_slices)
    si->begs[si->active_slices] = sg->beg + fs->beg;
  sg->msg_ctr--;

  if (sg->msg_ctr == 0)
    {
      progress_report(cp, ((mus_float_t)(si->active_slices) / (mus_float_t)(si->target_slices)));
      sg->status_needs_to_be_cleared = true;
      sg->msg_ctr = 8;
      if ((!(cp->graph_transform_on)) ||
	  (cp->active < CHANNEL_HAS_AXES))
	return(SONO_QUIT);
    }

  if (si->active_slices < si->total_slices)
    {
      int i;
      mus_float_t val;
      if (cp->transform_type == FOURIER)
	{
	  for (i = 0; i < sg->spectrum_size; i++)
	    {
	      val = fs->data[i];
	      if (val > si->scale) si->scale = val;
	      si->data[si->active_slices][i] = val;
	    }
	}
      else
	{
	  for (i = 0; i < sg->spectrum_size; i++)
	    {
	      val = fs->data[i];
	      if (val < 0.0) val = -val;  /* kinda dubious but I can't think of a good alternative */
	      if (val > si->scale) si->scale = val;
	      si->data[si->active_slices][i] = val;
	    }
	}
      si->active_slices++;
    }
  sg->outer++;
  if ((sg->outer == sg->outlim) || (!(cp->graph_transform_on)) || (cp->transform_graph_type == GRAPH_ONCE)) return(SONO_QUIT);
  sg->acc_hop += sg->hop;
  fs->beg = (mus_long_t)(sg->acc_hop);
  ap = cp->axis;
  if ((sg->losamp != ap->losamp) || (sg->hisamp != ap->hisamp))
    {
      fs->beg = 0;
      return(SONO_INIT);
    }
  return(SONO_RUN);
}


static void finish_sonogram(sonogram_state *sg)
{
  if (sg)
    {
      chan_info *cp;
      cp = sg->cp;
      if ((cp->active < CHANNEL_HAS_AXES) ||
	  (!(cp->graph_transform_on)))
	{
	  if (sg->fs)
	    sg->fs = free_fft_state(sg->fs);
	}
      else
	{
	  if ((sg->scp) && (sg->outlim > 1))
	    make_sonogram_axes(cp);
	  if (sg->fs)
	    sg->fs = free_fft_state(sg->fs);
	  cp->fft_data = NULL;
	  set_chan_fft_in_progress(cp, 0); /* i.e. clear it */
	  if ((sg->scp) && (sg->outlim > 1))
	    {
	      display_channel_fft_data(cp);
	      if (sg->outer == sg->outlim) sg->done = true;
	      sg->old_scale = (sg->scp)->scale;
	    }
	  else sg->done = true;
	  if ((cp->last_sonogram) && (cp->last_sonogram != sg))
	    free(cp->last_sonogram);
	  cp->last_sonogram = sg;
	  if (sg->status_needs_to_be_cleared)
	    {
	      finish_progress_report(cp);
	      sg->status_needs_to_be_cleared = false;
	    }
	}
    }
}


idle_func_t sonogram_in_slices(void *sono)
{
  sonogram_state *sg = (sonogram_state *)sono;
  chan_info *cp;
  cp = sg->cp;
  cp->temp_sonogram = NULL;
  if ((cp->active < CHANNEL_HAS_AXES) ||
      (!(cp->graph_transform_on)))
    {
      if ((sg) && (sg->fs)) sg->fs = free_fft_state(sg->fs);
      return(BACKGROUND_QUIT);
    }
  switch (sg->slice)
    {
    case SONO_INIT:
      sg->slice = set_up_sonogram(sg);
      break;

    case SONO_RUN:
      sg->slice = run_all_ffts(sg);
      break;

    case SONO_QUIT:
    default:
      finish_sonogram(sg);
      return(BACKGROUND_QUIT);
    }
  return(BACKGROUND_CONTINUE);
}


void sono_update(chan_info *cp)
{
  if (cp->transform_graph_type != GRAPH_ONCE) make_sonogram_axes(cp);
  update_graph(cp);
}


static void spectral_multiply(mus_float_t *rl1, mus_float_t *rl2, mus_long_t n)
{
  mus_long_t j, n2;
  mus_float_t invn;

  n2 = (int)(n * 0.5);
  invn = 0.25 / n;
  rl1[0] = ((rl1[0] * rl2[0]) / n);
  rl2[0] = 0.0;

  for (j = 1; j <= n2; j++)
    {
      mus_long_t nn2;
      mus_float_t rem, rep, aim, aip;
      nn2 = n - j;
      rep = (rl1[j] + rl1[nn2]);
      rem = (rl1[j] - rl1[nn2]);
      aip = (rl2[j] + rl2[nn2]);
      aim = (rl2[j] - rl2[nn2]);
      rl1[j] = invn * (rep * aip + aim * rem);
      rl1[nn2] = rl1[j];
      rl2[j] = invn * (aim * aip - rep * rem);
      rl2[nn2] = -rl2[j];
    }
}


void c_convolve(const char *fname, mus_float_t amp, int filec, mus_long_t filehdr, int filterc, mus_long_t filterhdr, mus_long_t filtersize,
		mus_long_t fftsize, int filter_chans, int filter_chan, mus_long_t data_size, snd_info *gsp)
{
  int err;

  /* need file to hold convolution output */
  err = mus_write_header(fname, MUS_NEXT, DEFAULT_OUTPUT_SRATE, 1, data_size * mus_bytes_per_sample(MUS_OUT_SAMPLE_TYPE), MUS_OUT_SAMPLE_TYPE, "c_convolve temp");
  if (err != MUS_NO_ERROR)
    snd_error("can't open convolution temp file %s: %s", fname, snd_io_strerror());
  else
    {
      mus_float_t *rl0, *rl1 = NULL;
      mus_float_t **pbuffer = NULL, **fbuffer = NULL;
      mus_long_t oloc;
      oloc = mus_header_data_location();

      /* get to start point in the two sound files and allocate space */
      lseek(filec, filehdr, SEEK_SET);
      lseek(filterc, filterhdr, SEEK_SET);

      rl0 = (mus_float_t *)calloc(fftsize, sizeof(mus_float_t));
      if (rl0) rl1 = (mus_float_t *)calloc(fftsize, sizeof(mus_float_t));
      if (rl1) pbuffer = (mus_float_t **)calloc(1, sizeof(mus_float_t *));
      if (pbuffer) pbuffer[0] = (mus_float_t *)calloc(data_size, sizeof(mus_float_t));
      fbuffer = (mus_float_t **)calloc(filter_chans, sizeof(mus_float_t *));
      if (fbuffer) fbuffer[filter_chan] = (mus_float_t *)calloc(filtersize, sizeof(mus_float_t));

      if ((!rl0) || (!rl1) ||
	  (!pbuffer) || (!pbuffer[0]) ||
	  (!fbuffer) || (!fbuffer[filter_chan]))
	{
	  snd_error("not enough memory for convolve of %s (filter size: %" print_mus_long ", fft size: %" print_mus_long ")",
		    fname, filtersize, fftsize);
	}
      else
	{
	  mus_float_t *pbuf;
	  mus_long_t i;
	  int tempfile;
	  chan_info *gcp;
	  gcp = gsp->chans[0];
	  start_progress_report(gcp);

	  tempfile = snd_reopen_write(fname);
	  snd_file_open_descriptors(tempfile, fname, MUS_OUT_SAMPLE_TYPE, oloc, 1, MUS_NEXT);
	  lseek(tempfile, oloc, SEEK_SET);

	  /* read in the "impulse response" */
	  pbuf = pbuffer[0];
	  mus_file_read_any(filterc, 0, filter_chans, filtersize, fbuffer, fbuffer);
	  for (i = 0; i < filtersize; i++)
	    rl1[i] = fbuffer[filter_chan][i];
	  progress_report(gcp, .1);

	  /* get the convolution data */
	  mus_file_read_any(filec, 0, 1, data_size, pbuffer, pbuffer);
	  for (i = 0; i < data_size; i++) rl0[i] = pbuf[i];

	  progress_report(gcp, .3);
	  mus_fft(rl0, rl1, fftsize, 1);
	  progress_report(gcp, .5);
	  spectral_multiply(rl0, rl1, fftsize);
	  progress_report(gcp, .6);
	  mus_fft(rl0, rl1, fftsize, -1);
	  progress_report(gcp, .8);

	  if (amp != 0.0)
	    {
	      mus_float_t scl;
	      /* normalize the results */
	      scl = 0.0;
	      for (i = 0; i < data_size; i++)
		{
		  mus_float_t val;
		  val = fabs(rl0[i]);
		  if (val > scl) scl = val;
		}
	      if (scl != 0.0) scl = amp / scl;
	      for (i = 0; i < data_size; i++)
		pbuf[i] = (scl * rl0[i]);
	    }
	  else
	    {
	      /* amp == 0.0 means un-normalized output */
	      mus_copy_floats(pbuf, rl0, data_size);
	    }
	  progress_report(gcp, .9);

	  /* and save as temp file */
	  mus_file_write(tempfile, 0, data_size - 1, 1, &(pbuf));
	  finish_progress_report(gcp);

	  if (mus_file_close(tempfile) != 0)
	    snd_error("convolve: can't close temp file %s!", fname);
	}
      if (rl0) free(rl0);
      if (rl1) free(rl1);
      if (pbuffer)
	{
	  if (pbuffer[0]) free(pbuffer[0]);
	  free(pbuffer);
	}
      if (fbuffer)
	{
	  if (fbuffer[filter_chan]) free(fbuffer[filter_chan]);
	  free(fbuffer);
	}
    }
}


/* -------------------------------------------------------------------------------- */

static void update_log_freq_fft_graph(chan_info *cp)
{
  if ((cp->active < CHANNEL_HAS_AXES) ||
      (!cp->sounds) ||
      (!cp->sounds[cp->sound_ctr]) ||
      (!(cp->graph_transform_on)) ||
      (!(cp->fft_log_frequency)) ||
      (chan_fft_in_progress(cp)))
    return;
  calculate_fft(cp);
}


void set_log_freq_start(mus_float_t base)
{
  in_set_log_freq_start(base);
  for_each_chan(update_log_freq_fft_graph);
}


static Xen g_log_freq_start(void) {return(C_double_to_Xen_real(log_freq_start(ss)));}

static Xen g_set_log_freq_start(Xen val)
{
  mus_float_t base;
  #define H_log_freq_start "(" S_log_freq_start "): log freq base (default: 25.0)"

  Xen_check_type(Xen_is_number(val), val, 1, S_set S_log_freq_start, "a number");
  base = Xen_real_to_C_double(val);
  if ((base < 0.0) ||
      (base > 100000.0))
    Xen_out_of_range_error(S_log_freq_start, 1, val, "base should be between 0 and srate/2");

  set_log_freq_start(base);
  reflect_log_freq_start_in_transform_dialog();

  return(C_double_to_Xen_real(log_freq_start(ss)));
}


static Xen g_show_selection_transform(void) {return(C_bool_to_Xen_boolean(show_selection_transform(ss)));}

static Xen g_set_show_selection_transform(Xen val)
{
  #define H_show_selection_transform "(" S_show_selection_transform "): " PROC_TRUE " if transform display reflects selection, not time-domain window"
  Xen_check_type(Xen_is_boolean(val), val, 1, S_set S_show_selection_transform, "a boolean");
  set_show_selection_transform(Xen_boolean_to_C_bool(val));
  return(C_bool_to_Xen_boolean(show_selection_transform(ss)));
}


static Xen g_transform_framples(Xen snd, Xen chn)
{
  #define H_transform_framples "(" S_transform_framples " :optional snd chn): \
return a description of transform graph data in snd's channel chn, based on " S_transform_graph_type ".\
If there is no transform graph, return 0; if " S_graph_once ", return " S_transform_size ",\
and otherwise return a list (spectrum-cutoff time-slices fft-bins)"

  chan_info *cp;
  sono_info *si;

  Snd_assert_channel(S_transform_framples, snd, chn, 1);
  cp = get_cp(snd, chn, S_transform_framples);
  if (!cp) return(Xen_false);
  if (!(cp->graph_transform_on))
    return(Xen_integer_zero);

  if (cp->transform_graph_type == GRAPH_ONCE)
    return(C_int_to_Xen_integer(cp->transform_size));
  si = cp->sonogram_data;

  if (si) return(Xen_list_3(C_double_to_Xen_real(cp->spectrum_end),
			    C_llong_to_Xen_llong(si->active_slices),
			    C_llong_to_Xen_llong(si->target_bins)));
  return(Xen_integer_zero);
}


static Xen g_transform_sample(Xen bin, Xen slice, Xen snd, Xen chn_n)
{
  #define H_transform_sample "(" S_transform_sample " :optional (bin 0) (slice 0) snd chn): \
return the current transform sample at bin and slice in snd channel chn (assuming sonogram or spectrogram)"

  chan_info *cp;

  Xen_check_type(Xen_is_llong_or_unbound(bin), bin, 1, S_transform_sample, "an integer");
  Xen_check_type(Xen_is_llong_or_unbound(slice), slice, 2, S_transform_sample, "an integer");

  Snd_assert_channel(S_transform_sample, snd, chn_n, 3);
  cp = get_cp(snd, chn_n, S_transform_sample);
  if (!cp) return(Xen_false);

  if (cp->graph_transform_on)
    {
      fft_info *fp;
      fp = cp->fft;
      if (fp)
	{
	  mus_long_t fbin = 0;
	  if (Xen_is_llong(bin)) fbin = Xen_llong_to_C_llong(bin);

	  if (fbin < fp->current_size)
	    {
	      if (cp->transform_graph_type == GRAPH_ONCE)
		return(C_double_to_Xen_real(fp->data[fbin]));
	      else
		{
		  mus_long_t fslice = 0;
		  sono_info *si;

		  if (Xen_is_llong(slice)) fslice = Xen_llong_to_C_llong(slice);
		  si = cp->sonogram_data;

		  if ((si) &&
		      (fbin < si->target_bins) &&
		      (fslice < si->active_slices))
		    return(C_double_to_Xen_real(si->data[fslice][fbin]));
		  else Xen_error(Xen_make_error_type("no-such-sample"),
				 Xen_list_8(C_string_to_Xen_string(S_transform_sample ": no such sample, bin: ~A, max bin: ~A, slice: ~A, max slice: ~A, sound index: ~A (~A), chan: ~A"),
					    bin,
					    C_int_to_Xen_integer((si) ? si->target_bins : 0),
					    slice,
					    C_int_to_Xen_integer((si) ? si->active_slices : 0),
					    snd,
					    C_string_to_Xen_string(cp->sound->short_filename),
					    chn_n));
		}
	    }
	  else Xen_error(Xen_make_error_type("no-such-sample"),
			 Xen_list_6(C_string_to_Xen_string(S_transform_sample ": no such sample, bin: ~A, max bin: ~A, sound index: ~A (~A), chan: ~A"),
				    bin,
				    C_int_to_Xen_integer(fp->current_size),
				    snd,
				    C_string_to_Xen_string(cp->sound->short_filename),
				    chn_n));
	}
    }
  return(Xen_false);
}


static Xen g_transform_to_vct(Xen snd, Xen chn_n, Xen v)
{
  #define H_transform_to_vct "(" S_transform_to_vct " :optional snd chn obj): \
return a " S_vct " (obj if it's passed), with the current transform data from snd's channel chn"

  chan_info *cp;
  vct *v1;
  v1 = xen_to_vct(v);

  Snd_assert_channel(S_transform_to_vct, snd, chn_n, 1);
  cp = get_cp(snd, chn_n, S_transform_to_vct);
  if (!cp) return(Xen_false);

  if ((cp->graph_transform_on) &&
      (cp->fft))
    {
      mus_long_t len;
      mus_float_t *fvals;

      if (cp->transform_graph_type == GRAPH_ONCE)
	{
	  fft_info *fp;
	  fp = cp->fft;
	  len = fp->current_size;
	  if (v1)
	    fvals = mus_vct_data(v1);
	  else fvals = (mus_float_t *)malloc(len * sizeof(mus_float_t));
	  mus_copy_floats(fvals, fp->data, len);
	  if (v1)
	    return(v);
	  else return(xen_make_vct(len, fvals));
	}
      else
	{
	  sono_info *si;
	  si = cp->sonogram_data;

	  if (si)
	    {
	      int i, j, k, bins, slices;
	      slices = si->active_slices;
	      bins = si->target_bins;
	      len = bins * slices;
	      if (v1)
		fvals = mus_vct_data(v1);
	      else fvals = (mus_float_t *)calloc(len, sizeof(mus_float_t));
	      for (i = 0, k = 0; i < slices; i++)
		for (j = 0; j < bins; j++, k++)
		  fvals[k] = si->data[i][j];
	      if (v1)
		return(v);
	      else return(xen_make_vct(len, fvals));
	    }
	}
    }
  return(Xen_false);
}



/* ---------------------------------------- transform objects ---------------------------------------- */

typedef struct {
  int n;
} xen_transform;


#define Xen_to_xen_transform(arg) ((xen_transform *)Xen_object_ref(arg))

int xen_transform_to_int(Xen n)
{
  xen_transform *col;
  col = Xen_to_xen_transform(n);
  return(col->n);
}


static Xen_object_type_t xen_transform_tag;

bool xen_is_transform(Xen obj)
{
  return(Xen_c_object_is_type(obj, xen_transform_tag));
}

#if (!HAVE_SCHEME)
static void xen_transform_free(xen_transform *v) {if (v) free(v);}

Xen_wrap_free(xen_transform, free_xen_transform, xen_transform_free)
#else
static s7_pointer s7_xen_transform_free(s7_scheme *sc, s7_pointer obj)
{
  xen_transform *v;
  v = (xen_transform *)s7_c_object_value(obj);
  if (v) free(v);
  return(NULL);
}
#endif

static char *xen_transform_to_string(xen_transform *v)
{
  #define TRANSFORM_PRINT_BUFFER_SIZE 64
  char *buf;
  if (!v) return(NULL);
  buf = (char *)calloc(TRANSFORM_PRINT_BUFFER_SIZE, sizeof(char));
  snprintf(buf, TRANSFORM_PRINT_BUFFER_SIZE, "#<transform %s>", transform_name(v->n));
  return(buf);
}


#if HAVE_FORTH || HAVE_RUBY
Xen_wrap_print(xen_transform, print_xen_transform, xen_transform_to_string)

static Xen g_xen_transform_to_string(Xen obj)
{
  char *vstr;
  Xen result;
  #define S_xen_transform_to_string "transform->string"
  Xen_check_type(xen_is_transform(obj), obj, 1, S_xen_transform_to_string, "a transform");
  vstr = xen_transform_to_string(Xen_to_xen_transform(obj));
  result = C_string_to_Xen_string(vstr);
  free(vstr);
  return(result);
}
#else
#if HAVE_SCHEME
static s7_pointer g_xen_transform_to_string(s7_scheme *sc, s7_pointer args)
{
  char *vstr;
  s7_pointer result;
  vstr = xen_transform_to_string(Xen_to_xen_transform(s7_car(args)));
  result = C_string_to_Xen_string(vstr);
  free(vstr);
  return(result);
}
#endif
#endif


#if (!HAVE_SCHEME)
static bool xen_transform_equalp(xen_transform *v1, xen_transform *v2)
{
  return((v1 == v2) ||
	 (v1->n == v2->n));
}

static Xen equalp_xen_transform(Xen obj1, Xen obj2)
{
  if ((!(xen_is_transform(obj1))) || (!(xen_is_transform(obj2)))) return(Xen_false);
  return(C_bool_to_Xen_boolean(xen_transform_equalp(Xen_to_xen_transform(obj1), Xen_to_xen_transform(obj2))));
}
#endif


static xen_transform *xen_transform_make(int n)
{
  xen_transform *new_v;
  new_v = (xen_transform *)malloc(sizeof(xen_transform));
  new_v->n = n;
  return(new_v);
}


Xen new_xen_transform(int n)
{
  xen_transform *mx;
  if (n < 0)
    return(Xen_false);

  mx = xen_transform_make(n);
  return(Xen_make_object(xen_transform_tag, mx, 0, free_xen_transform));
}

#define C_int_to_Xen_transform(Val) new_xen_transform(Val)


#if HAVE_SCHEME
static s7_pointer s7_xen_transform_is_equal(s7_scheme *sc, s7_pointer args)
{
  s7_pointer p1, p2;
  p1 = s7_car(args);
  p2 = s7_cadr(args);
  if (p1 == p2) return(s7_t(sc));
  if (s7_c_object_type(p2) == xen_transform_tag)
    return(s7_make_boolean(sc, (((xen_transform *)s7_c_object_value(p1))->n == ((xen_transform *)s7_c_object_value(p2))->n)));
  return(s7_f(sc));
}

static Xen s7_xen_transform_length(s7_scheme *sc, Xen args)
{
  return(C_int_to_Xen_integer(transform_size(ss)));
}
#endif


static void init_xen_transform(void)
{
#if HAVE_SCHEME
  xen_transform_tag = s7_make_c_type(s7, "<transform>");
  s7_c_type_set_gc_free(s7, xen_transform_tag, s7_xen_transform_free);
  s7_c_type_set_is_equal(s7, xen_transform_tag, s7_xen_transform_is_equal);
  s7_c_type_set_length(s7, xen_transform_tag, s7_xen_transform_length);
  s7_c_type_set_to_string(s7, xen_transform_tag, g_xen_transform_to_string);
#else
#if HAVE_RUBY
  xen_transform_tag = Xen_make_object_type("XenTransform", sizeof(xen_transform));
#else
  xen_transform_tag = Xen_make_object_type("Transform", sizeof(xen_transform));
#endif
#endif

#if HAVE_FORTH
  fth_set_object_inspect(xen_transform_tag,   print_xen_transform);
  fth_set_object_dump(xen_transform_tag,      g_xen_transform_to_string);
  fth_set_object_equal(xen_transform_tag,     equalp_xen_transform);
  fth_set_object_free(xen_transform_tag,      free_xen_transform);
#endif

#if HAVE_RUBY
  rb_define_method(xen_transform_tag, "to_s",     Xen_procedure_cast print_xen_transform, 0);
  rb_define_method(xen_transform_tag, "eql?",     Xen_procedure_cast equalp_xen_transform, 1);
  rb_define_method(xen_transform_tag, "==",       Xen_procedure_cast equalp_xen_transform, 1);
  rb_define_method(xen_transform_tag, "to_str",   Xen_procedure_cast g_xen_transform_to_string, 0);
#endif
}


static Xen g_integer_to_transform(Xen n)
{
  #define H_integer_to_transform "(" S_integer_to_transform " n) returns a transform object corresponding to the given integer"
  int type;
  Xen_check_type(Xen_is_integer(n), n, 1, S_integer_to_transform, "an integer");
  type = Xen_integer_to_C_int(n);
  if (is_transform(type))
    return(new_xen_transform(type));
  return(Xen_false);
}


static Xen g_transform_to_integer(Xen n)
{
  #define H_transform_to_integer "(" S_transform_to_integer " id) returns the integer corresponding to the given transform"
  Xen_check_type(xen_is_transform(n), n, 1, S_transform_to_integer, "a transform");
  return(C_int_to_Xen_integer(xen_transform_to_int(n)));
}


static Xen g_is_transform(Xen type)
{
  #define H_is_transform "(" S_is_transform " obj): " PROC_TRUE " if 'obj' is a transform object."
  return(C_bool_to_Xen_boolean(xen_is_transform(type) &&
			  is_transform(Xen_transform_to_C_int(type))));
}


static Xen g_snd_transform(Xen type, Xen data, Xen hint)
{
  #define H_snd_transform "(snd-transform type data choice) calls whatever FFT is being used by the \
display.  'type' is a transform object such as " S_fourier_transform "; 'data' is a " S_vct ". In the wavelet case, \
'choice' is the wavelet to use."

  int trf, hnt;
  mus_long_t i, j, n2, vlen;
  vct *v;
  mus_float_t *dat, *vdata;

  Xen_check_type(xen_is_transform(type), type, 1, "snd-transform", "a transform object");
  Xen_check_type(mus_is_vct(data), data, 2, "snd-transform", "a " S_vct);

  trf = Xen_transform_to_C_int(type);
  if ((trf < 0) || (trf > HAAR))
    Xen_out_of_range_error("snd-transform", 1, type, "invalid transform choice");

  v = Xen_to_vct(data);
  vlen = mus_vct_length(v);
  vdata = mus_vct_data(v);

  switch (trf)
    {
    case FOURIER:
      n2 = vlen / 2;
      dat = (mus_float_t *)calloc(vlen, sizeof(mus_float_t));
      mus_fft(vdata, dat, vlen, 1);
      vdata[0] *= vdata[0];
      vdata[n2] *= vdata[n2];
      for (i = 1, j = vlen - 1; i < n2; i++, j--)
	{
	  vdata[i] = vdata[i] * vdata[i] + dat[i] * dat[i];
	  vdata[j] = vdata[i];
	}
      free(dat);
      break;

    case WAVELET:
      hnt = Xen_integer_to_C_int(hint);
      if (hnt < NUM_WAVELETS)
	wavelet_transform(vdata, vlen, wavelet_data[hnt], wavelet_sizes[hnt]);
      break;

    case HAAR:
      haar_transform(vdata, vlen);
      break;

    case CEPSTRUM:
      mus_cepstrum(vdata, vlen);
      break;

    case WALSH:
      walsh_transform(vdata, vlen);
      break;

    case AUTOCORRELATION:
      mus_autocorrelate(vdata, vlen);
      break;
    }
  return(data);
}


static Xen g_add_transform(Xen name, Xen xlabel, Xen lo, Xen hi, Xen proc)
{
  #define H_add_transform "(" S_add_transform " name x-label low high func): add the transform func \
to the transform lists; func should be a function of two arguments, the length of the transform \
and a sampler to get the data, and should return a " S_vct " containing the transform results. \
name is the transform's name, x-label is its x-axis label, and the relevant returned data \
to be displayed goes from low to high (normally 0.0 to 1.0).  " S_add_transform " returns the new transform object."

  char *errmsg;
  errmsg = procedure_ok(proc, 2, S_add_transform, "transform", 5);
  if (errmsg)
    {
      Xen errstr;
      errstr = C_string_to_Xen_string(errmsg);
      free(errmsg);
      return(snd_bad_arity_error(S_add_transform, errstr, proc));
    }

#if HAVE_SCHEME
  if (mus_is_xen(proc)) /* happens a lot in snd-test.scm, so add a check */
    Xen_wrong_type_arg_error(S_add_transform, 5, proc, "a procedure");
#endif

  Xen_check_type(Xen_is_string(name), name, 1, S_add_transform, "a string");
  Xen_check_type(Xen_is_string(xlabel), xlabel, 2, S_add_transform, "a string");
  Xen_check_type(Xen_is_number(lo), lo, 3, S_add_transform, "a number");
  Xen_check_type(Xen_is_number(hi), hi, 4, S_add_transform, "a number");
  Xen_check_type(Xen_is_procedure(proc), proc, 5, S_add_transform, "a procedure");

  return(C_int_to_Xen_transform(add_transform(Xen_string_to_C_string(name),
					      Xen_string_to_C_string(xlabel),
					      Xen_real_to_C_double(lo),
					      Xen_real_to_C_double(hi),
					      proc)));
}


static int deleted_type = 0;
static void unset_deleted_transform_type(chan_info *cp)
{
  if (cp->transform_type == deleted_type)
    cp->transform_type = DEFAULT_TRANSFORM_TYPE;
}


static Xen g_delete_transform(Xen type)
{
  int typ;
  added_transform *af;
  #define H_delete_transform "(" S_delete_transform " obj) deletes the specified transform if it was created via " S_add_transform "."

  Xen_check_type(xen_is_transform(type), type, 1, S_delete_transform, "a transform");

  typ = Xen_transform_to_C_int(type);
  if ((typ < NUM_BUILTIN_TRANSFORM_TYPES) || (!is_transform(typ)))
    Xen_out_of_range_error(S_delete_transform, 1, type, "an integer (an active added transform)");

  af = type_to_transform(typ);
  if (af)
    {
      added_transforms[af->type - NUM_BUILTIN_TRANSFORM_TYPES] = NULL;
      free(af->xlabel);
      free(af->name);
      snd_unprotect_at(af->gc_loc);
      af->proc = Xen_false;
      free(af);
      /* now make sure nobody expects to use that transform */
      if (transform_type(ss) == typ) set_transform_type(DEFAULT_TRANSFORM_TYPE);
      deleted_type = typ;
      for_each_chan(unset_deleted_transform_type);
      return(Xen_true);
    }
  return(Xen_false);
}


Xen_wrap_2_optional_args(g_transform_framples_w, g_transform_framples)
Xen_wrap_4_optional_args(g_transform_sample_w, g_transform_sample)
Xen_wrap_3_optional_args(g_transform_to_vct_w, g_transform_to_vct)
Xen_wrap_5_args(g_add_transform_w, g_add_transform)
Xen_wrap_3_optional_args(g_snd_transform_w, g_snd_transform)
Xen_wrap_1_arg(g_is_transform_w, g_is_transform)
Xen_wrap_1_arg(g_delete_transform_w, g_delete_transform)
Xen_wrap_no_args(g_log_freq_start_w, g_log_freq_start)
Xen_wrap_1_arg(g_set_log_freq_start_w, g_set_log_freq_start)
Xen_wrap_no_args(g_show_selection_transform_w, g_show_selection_transform)
Xen_wrap_1_arg(g_set_show_selection_transform_w, g_set_show_selection_transform)
Xen_wrap_1_arg(g_integer_to_transform_w, g_integer_to_transform)
Xen_wrap_1_arg(g_transform_to_integer_w, g_transform_to_integer)

#if HAVE_SCHEME
static s7_pointer acc_log_freq_start(s7_scheme *sc, s7_pointer args) {return(g_set_log_freq_start(s7_cadr(args)));}
static s7_pointer acc_show_selection_transform(s7_scheme *sc, s7_pointer args) {return(g_set_show_selection_transform(s7_cadr(args)));}
#endif

#if (!HAVE_SCHEME)
static Xen transform_temp[6]; /* static for Ruby's sake */
#endif

void g_init_fft(void)
{
#if HAVE_SCHEME
  s7_pointer i, b, p, t, fv, r, s, tr, pr;
  i = s7_make_symbol(s7, "integer?");
  b = s7_make_symbol(s7, "boolean?");
  p = s7_make_symbol(s7, "pair?");
  fv = s7_make_symbol(s7, "float-vector?");
  r = s7_make_symbol(s7, "real?");
  s = s7_make_symbol(s7, "string?");
  tr = s7_make_symbol(s7, "transform?");
  pr = s7_make_symbol(s7, "procedure?");
  t = s7_t(s7);
#endif

#if HAVE_SCHEME
  #define H_before_transform_hook S_before_transform_hook " (snd chn): called just before a transform is calculated.  If it returns \
an integer, it is used as the starting point of the transform.  The following \
somewhat brute-force code shows a way to have the transform reflect the position \
of a moving mark:\n\
  (define transform-position #f)\n\
  (hook-push " S_before_transform_hook "\n\
    (lambda (snd chn) transform-position))\n\
  (hook-push " S_mark_drag_hook "\n\
    (lambda (id)\n\
      (set! transform-position (" S_mark_sample " id))\n\
      (" S_update_transform_graph ")))"
#endif

#if HAVE_RUBY
  #define H_before_transform_hook S_before_transform_hook " (snd chn): called just before a transform is calculated.  If it returns \
an integer, it is used as the starting point of the transform.  The following \
somewhat brute-force code shows a way to have the transform reflect the position \
of a moving mark:\n\
$transform_position = false\n\
$before_transform_hook.add_hook!(\"snd-fft\") |snd, chn|\n\
  $transform_position\n\
end\n\
$mark_drag_hook.add_hook!(\"snd-fft\") |id|\n\
  $transform_position = " S_mark_sample "(id)\n\
  " S_update_transform_graph "\n\
end"
#endif

#if HAVE_FORTH
  #define H_before_transform_hook S_before_transform_hook " (snd chn): called just before a transform is calculated.  If it returns \
an integer, it is used as the starting point of the transform.  The following \
somewhat brute-force code shows a way to have the transform reflect the position \
of a moving mark:\n\
#f value transform-position\n\
" S_before_transform_hook " lambda: <{ snd chn }> transform-position ; add-hook!\n\
" S_mark_drag_hook " lambda: <{ id }>\n\
  id " S_mark_sample " to transform-position\n\
  " S_update_transform_graph "\n\
; add-hook!"
#endif

  init_xen_transform();

  before_transform_hook = Xen_define_hook(S_before_transform_hook, "(make-hook 'snd 'chn)", 2, H_before_transform_hook);

#if HAVE_SCHEME
  s7_define_constant(s7, S_fourier_transform, C_int_to_Xen_transform(FOURIER));
  s7_define_constant(s7, S_wavelet_transform, C_int_to_Xen_transform(WAVELET));
  s7_define_constant(s7, S_haar_transform,    C_int_to_Xen_transform(HAAR));
  s7_define_constant(s7, S_cepstrum,          C_int_to_Xen_transform(CEPSTRUM));
  s7_define_constant(s7, S_walsh_transform,   C_int_to_Xen_transform(WALSH));
  s7_define_constant(s7, S_autocorrelation,   C_int_to_Xen_transform(AUTOCORRELATION));
  /* *transform-type* is #<transform Fourier> by default */
#else
  Xen_define_variable(S_fourier_transform, transform_temp[0], C_int_to_Xen_transform(FOURIER));
  Xen_define_variable(S_wavelet_transform, transform_temp[1], C_int_to_Xen_transform(WAVELET));
  Xen_define_variable(S_haar_transform,    transform_temp[2], C_int_to_Xen_transform(HAAR));
  Xen_define_variable(S_cepstrum,          transform_temp[3], C_int_to_Xen_transform(CEPSTRUM));
  Xen_define_variable(S_walsh_transform,   transform_temp[4], C_int_to_Xen_transform(WALSH));
  Xen_define_variable(S_autocorrelation,   transform_temp[5], C_int_to_Xen_transform(AUTOCORRELATION));
#endif

  #define H_dont_normalize "The value for " S_transform_normalization " that causes the transform to display raw data"
  #define H_normalize_by_channel "The value for " S_transform_normalization " that causes the transform to be normalized in each channel independently"
  #define H_normalize_by_sound "The value for " S_transform_normalization " that causes the transform to be normalized across a sound's channels"
  #define H_normalize_globally "The value for " S_transform_normalization " that causes the transform to be normalized across all sounds"

  Xen_define_constant(S_dont_normalize,        DONT_NORMALIZE,       H_dont_normalize);
  Xen_define_constant(S_normalize_by_channel,  NORMALIZE_BY_CHANNEL, H_normalize_by_channel);
  Xen_define_constant(S_normalize_by_sound,    NORMALIZE_BY_SOUND,   H_normalize_by_sound);
  Xen_define_constant(S_normalize_globally,    NORMALIZE_GLOBALLY,   H_normalize_globally);

  Xen_define_typed_procedure(S_transform_framples,   g_transform_framples_w, 0, 2, 0, H_transform_framples,
			     s7_make_signature(s7, 3, s7_make_signature(s7, 2, i, p), t, t));
  Xen_define_typed_procedure(S_transform_sample,     g_transform_sample_w, 0, 4, 0, H_transform_sample, s7_make_signature(s7, 5, r, i, i, t, t));
  Xen_define_typed_procedure(S_transform_to_vct,     g_transform_to_vct_w, 0, 3, 0, H_transform_to_vct, s7_make_signature(s7, 4, fv, t, t, fv));
  Xen_define_typed_procedure(S_add_transform,        g_add_transform_w,    5, 0, 0, H_add_transform,    s7_make_signature(s7, 6, tr, s, s, r, r, pr));
  Xen_define_typed_procedure(S_is_transform,         g_is_transform_w,     1, 0, 0, H_is_transform,     s7_make_signature(s7, 2, b, t));
  Xen_define_typed_procedure(S_delete_transform,     g_delete_transform_w, 1, 0, 0, H_delete_transform, s7_make_signature(s7, 2, b, tr));
  Xen_define_typed_procedure("snd-transform",        g_snd_transform_w,    2, 1, 0, H_snd_transform,    s7_make_signature(s7, 4, fv, tr, fv, i));

  Xen_define_typed_dilambda(S_log_freq_start, g_log_freq_start_w, H_log_freq_start,
			    S_set S_log_freq_start, g_set_log_freq_start_w, 0, 0, 1, 0,
			    s7_make_signature(s7, 1, r), s7_make_signature(s7, 2, r, r));

  Xen_define_typed_dilambda(S_show_selection_transform, g_show_selection_transform_w, H_show_selection_transform,
			    S_set S_show_selection_transform, g_set_show_selection_transform_w, 0, 0, 1, 0,
			    s7_make_signature(s7, 1, b), s7_make_signature(s7, 2, b, b));

  Xen_define_typed_procedure(S_integer_to_transform, g_integer_to_transform_w, 1, 0, 0, H_integer_to_transform, s7_make_signature(s7, 2, tr, i));
  Xen_define_typed_procedure(S_transform_to_integer, g_transform_to_integer_w, 1, 0, 0, H_transform_to_integer, s7_make_signature(s7, 2, i, tr));

#if HAVE_SCHEME
  s7_set_setter(s7, ss->log_freq_start_symbol, s7_make_function(s7, "[acc-" S_log_freq_start "]", acc_log_freq_start, 2, 0, false, "accessor"));
  s7_set_setter(s7, ss->show_selection_transform_symbol, s7_make_function(s7, "[acc-" S_show_selection_transform "]", acc_show_selection_transform, 2, 0, false, "accessor"));

  s7_set_documentation(s7, ss->log_freq_start_symbol, "*log-freq-start*: log freq base (25.0)");
  s7_set_documentation(s7, ss->show_selection_transform_symbol, "*show-selection-transform*: #t if transform display reflects selection, not time-domain window");
#endif
}

/* display by wavelength is not so useful in the context of sound because
 *  the frequencies span say 8-10 octaves.  Even if we place the start point
 *  at 20Hz, that still puts the (linear) middle at 40 Hz (17 meters->1 inch
 *  is basically the range); if we goof around with log freq axis, why bother
 *  with the inverse?
 */