From 1dcb134a763b4598c98c86d77e71e4f026bcc2e1 Mon Sep 17 00:00:00 2001 From: Scott Wilson Date: Mon, 10 Jun 2024 13:58:26 +0100 Subject: [PATCH] Update channel limits and VBAP 1.032 Adjustments ported from new PD code. LS amounts allocated dynamically. VBAP 1.0.3.2 Remove ls hardcoded limit --- source/VBAPUGens/VBAP.cpp | 213 +++++++++++++++++++++--------------- source/VBAPUGens/sc/vbap.sc | 154 ++++++++++++++------------ 2 files changed, 209 insertions(+), 158 deletions(-) diff --git a/source/VBAPUGens/VBAP.cpp b/source/VBAPUGens/VBAP.cpp index b4b14d0c8a..4059b7c09a 100644 --- a/source/VBAPUGens/VBAP.cpp +++ b/source/VBAPUGens/VBAP.cpp @@ -87,19 +87,19 @@ using nova::wrap_argument; #endif -#define RES_ID 9171 /* resource ID for assistance (we'll add that later) */ -#define MAX_LS_SETS 100 /* maximum number of loudspeaker sets (triplets or pairs) allowed */ -#define MAX_LS_AMOUNT 55 /* maximum amount of loudspeakers, can be increased */ static InterfaceTable *ft; +static float rad2ang = 360.0 / ( 2.0f * pi ); +static float atorad = (2.0f * pi) / 360.0f; + struct VBAP : Unit { float x_azi; /* panning direction azimuth */ float x_ele; /* panning direction elevation */ - float x_set_inv_matx[MAX_LS_SETS][9]; /* inverse matrice for each loudspeaker set */ - float x_set_matx[MAX_LS_SETS][9]; /* matrice for each loudspeaker set */ - int x_lsset[MAX_LS_SETS][3]; /* channel numbers of loudspeakers in each LS set */ + float* *x_set_inv_matx; /* inverse matrice for each loudspeaker set */ + float* *x_set_matx; /* matrice for each loudspeaker set */ + int* *x_lsset; /* channel numbers of loudspeakers in each LS set */ int x_lsset_available; /* have loudspeaker sets been defined with define_loudspeakers */ int x_lsset_amount; /* amount of loudspeaker sets */ int x_ls_amount; /* amount of loudspeakers */ @@ -108,7 +108,7 @@ struct VBAP : Unit float x_spread_base[3]; /* used to create uniform spreading */ float *final_gs; - float m_chanamp[MAX_LS_AMOUNT]; // for smoothing amp changes max channels 55 at the moment + float *m_chanamp; }; // for circular smoothing @@ -140,7 +140,6 @@ extern "C" static void angle_to_cart(float azi, float ele, float res[3]) /* converts angular coordinates to cartesian */ { - float atorad = (2 * 3.1415927 / 360) ; res[0] = cos((float) azi * atorad) * cos((float) ele * atorad); res[1] = sin((float) azi * atorad) * cos((float) ele * atorad); res[2] = sin((float) ele * atorad); @@ -149,10 +148,6 @@ static void angle_to_cart(float azi, float ele, float res[3]) static void cart_to_angle(float cvec[3], float avec[3]) /* converts cartesian coordinates to angular */ { - // float tmp, tmp2, tmp3, tmp4; /* warning: unused variable */ - float atorad = (2 * 3.1415927 / 360) ; - float pi = 3.1415927; - // float power; /* warning: unused variable */ float dist, atan_y_per_x, atan_x_pl_y_per_z; float azi, ele; @@ -187,14 +182,13 @@ static void new_spread_dir(VBAP *x, float spreaddir[3], float vscartdir[3], floa { float beta,gamma; float a,b; - float pi = 3.1415927; float power; gamma = acos(vscartdir[0] * spread_base[0] + vscartdir[1] * spread_base[1] + vscartdir[2] * spread_base[2])/pi*180; - if(fabs(gamma) < 1){ - angle_to_cart(x->x_azi+90, 0, spread_base); + if(fabs(gamma) < 1 || fabs(gamma) > 179 || fabs(gamma) < -179){ + angle_to_cart(x->x_azi+90.0, 0, spread_base); gamma = acos(vscartdir[0] * spread_base[0] + vscartdir[1] * spread_base[1] + vscartdir[2] * spread_base[2])/pi*180; @@ -217,7 +211,6 @@ static void new_spread_base(VBAP *x, float spreaddir[3], float vscartdir[3]) /* subroutine for spreading */ { float d; - float pi = 3.1415927; float power; d = cos(x->x_spread/180*pi); @@ -254,11 +247,9 @@ static void additive_vbap(float *final_gs, float cartdir[3], VBAP *x) float small_g; float big_sm_g, gtmp[3]; int winner_set; - // float new_cartdir[3]; /* warning: unused variable */ - // float new_angle_dir[3]; /* warning: unused variable */ int dim = x->x_dimension; int neg_g_am, best_neg_g_am; - float g[3]; + float g[3] = { 0, 0, 0 }; int ls[3] = { 0, 0, 0 }; big_sm_g = -100000.0; @@ -294,33 +285,27 @@ static void additive_vbap(float *final_gs, float cartdir[3], VBAP *x) g[2]=0.0; ls[2]=0; } - } + } } - + gains_modified=0; for(i=0;i */ - if(dim==3) - final_gs[ls[2]-1] += g[2]; - } + + // new from PD + if(gains_modified != 1){ + power=sqrt(g[0]*g[0] + g[1]*g[1] + g[2]*g[2]); + g[0] /= power; + g[1] /= power; + g[2] /= power; + + final_gs[ls[0]-1] += g[0]; + final_gs[ls[1]-1] += g[1]; + if (dim==3) + final_gs[ls[2]-1] += g[2]; + } } static void spread_it(VBAP *x, float *final_gs) @@ -413,23 +398,28 @@ static void vbap(float g[3], int ls[3], VBAP *x) float new_angle_dir[3]; int dim = x->x_dimension; int neg_g_am, best_neg_g_am; - - /* transfering the azimuth angle to a decent value */ - while(x->x_azi > 180) - x->x_azi -= 360; - while(x->x_azi < -179) - x->x_azi += 360; - - /* transferring the elevation to a decent value */ - if(dim == 3){ - while(x->x_ele > 180) - x->x_ele -= 360; - while(x->x_ele < -179) - x->x_ele += 360; - } else - x->x_ele = 0; - - + + // new from PD + // transfering the azimuth angle to a decent value + if(x->x_azi > 360.0 || x->x_azi < -360.0) + x->x_azi = fmod(x->x_azi, 360.0); + if(x->x_azi > 180.0) + x->x_azi -= 360.0; + if(x->x_azi < -179.0) + x->x_azi += 360.0; + + + // transferring the elevation to a decent value + if(dim == 3){ + if(x->x_ele > 360.0 || x->x_ele < -360.0) + x->x_ele = fmod(x->x_ele, 360.0); + if(x->x_ele > 180.0) + x->x_ele -= 360.0; + if(x->x_ele < -179.0) + x->x_ele += 360.0; + } else + x->x_ele = 0.0; + /* go through all defined loudspeaker sets and find the set which // has all positive values. If such is not found, set with largest // minimum value is chosen. If at least one of gain factors of one LS set is negative @@ -472,29 +462,30 @@ static void vbap(float g[3], int ls[3], VBAP *x) // calculate direction that corresponds to these new // gain values. This happens when the virtual source is outside of // all loudspeaker sets. */ - - if(dim==3){ - gains_modified=0; - for(i=0;ix_set_matx[winner_set][0] * g[0] - + x->x_set_matx[winner_set][1] * g[1] - + x->x_set_matx[winner_set][2] * g[2]; - new_cartdir[1] = x->x_set_matx[winner_set][3] * g[0] - + x->x_set_matx[winner_set][4] * g[1] - + x->x_set_matx[winner_set][5] * g[2]; - new_cartdir[2] = x->x_set_matx[winner_set][6] * g[0] - + x->x_set_matx[winner_set][7] * g[1] - + x->x_set_matx[winner_set][8] * g[2]; - cart_to_angle(new_cartdir,new_angle_dir); - x->x_azi = (float) (new_angle_dir[0] + 0.5); - x->x_ele = (float) (new_angle_dir[1] + 0.5); - } - } + + gains_modified=0; + for(i=0;ix_set_matx[winner_set][0] * g[0] + + x->x_set_matx[winner_set][1] * g[1] + + x->x_set_matx[winner_set][2] * g[2]; + new_cartdir[1] = x->x_set_matx[winner_set][3] * g[0] + + x->x_set_matx[winner_set][4] * g[1] + + x->x_set_matx[winner_set][5] * g[2]; + if (dim == 3) { + new_cartdir[2] = x->x_set_matx[winner_set][6] * g[0] + + x->x_set_matx[winner_set][7] * g[1] + + x->x_set_matx[winner_set][8] * g[2]; + } else new_cartdir[2] = 0; + cart_to_angle(new_cartdir,new_angle_dir); + x->x_azi = (new_angle_dir[0]); + x->x_ele = (new_angle_dir[1]); + } power=sqrt(g[0]*g[0] + g[1]*g[1] + g[2]*g[2]); g[0] /= power; @@ -605,10 +596,14 @@ static inline_functions void VBAP_next_simd(VBAP *unit, int inNumSamples) } #endif +// needs to check that numOutputs and x_ls_amount match!! static void VBAP_Ctor(VBAP* unit) { + //printf("VBAP-1.0.3.2\n"); int numOutputs = unit->mNumOutputs, counter = 0, datapointer=0, setpointer=0, i; - + + unit->m_chanamp = (float*)RTAlloc(unit->mWorld, numOutputs * sizeof(float)); + // initialise interpolation levels and outputs for (int i=0; im_chanamp[i] = 0; @@ -660,6 +655,15 @@ static void VBAP_Ctor(VBAP* unit) // return; } + unit->x_set_inv_matx = (float**)RTAlloc(unit->mWorld, counter * sizeof(float*)); + unit->x_set_matx = (float**)RTAlloc(unit->mWorld, counter * sizeof(float*)); + unit->x_lsset = (int**)RTAlloc(unit->mWorld, counter * sizeof(int*)); + + for(i=0; ix_set_inv_matx[i] = (float*)RTAlloc(unit->mWorld, 9 * sizeof(float)); + unit->x_set_matx[i] = (float*)RTAlloc(unit->mWorld, 9 * sizeof(float)); + unit->x_lsset[i] = (int*)RTAlloc(unit->mWorld, 3 * sizeof(int)); + } // probably sets should be created with rtalloc while(counter-- > 0){ for(i=0; i < unit->x_dimension; i++){ @@ -688,23 +692,52 @@ static void VBAP_Ctor(VBAP* unit) #endif SETCALC(VBAP_next); - if (unit->x_lsset_available == 1) { - unit->x_spread_base[0] = 0.0; - unit->x_spread_base[1] = 1.0; - unit->x_spread_base[2] = 0.0; - VBAP_next(unit, 1); // calculate initial gain factors && compute initial sample - } else { - ZOUT0(0) = 0; + ZOUT0(0) = ZIN0(0); + unit->x_azi = ZIN0(2); + unit->x_ele = ZIN0(3); + unit->x_spread_base[0] = 0.0; + unit->x_spread_base[1] = 1.0; + unit->x_spread_base[2] = 0.0; + unit->x_spread = ZIN0(4); + + // calculate initial gain factors + float g[3]; + int ls[3]; + float *final_gs = unit->final_gs; + + if(unit->x_lsset_available ==1){ + vbap(g,ls, unit); + for(i=0;ix_ls_amount;i++) + final_gs[i]=0.0; + for(i=0;ix_dimension;i++){ + final_gs[ls[i]-1]=g[i]; + } + if(unit->x_spread != 0){ + spread_it(unit,final_gs); + } + + } else { // if the ls data was bad, just set every gain to 0 and bail for(i=0;ix_ls_amount;i++) - unit->final_gs[i]=0.f; - } + final_gs[i]=0.f; + } } static void VBAP_Dtor(VBAP* unit) { + int counter = unit->x_lsset_amount; RTFree(unit->mWorld, unit->final_gs); + + for(int i=0; imWorld, unit->x_set_inv_matx[i]); + RTFree(unit->mWorld, unit->x_set_matx[i]); + RTFree(unit->mWorld, unit->x_lsset[i]); + } + + RTFree(unit->mWorld, unit->x_set_inv_matx); + RTFree(unit->mWorld, unit->x_set_matx); + RTFree(unit->mWorld, unit->x_lsset); } ////////////////////////////////////////////////////////////////////////////////////////////////// diff --git a/source/VBAPUGens/sc/vbap.sc b/source/VBAPUGens/sc/vbap.sc index 0f1e9a2fad..00438128dc 100644 --- a/source/VBAPUGens/sc/vbap.sc +++ b/source/VBAPUGens/sc/vbap.sc @@ -1,6 +1,6 @@ /* VBAP created by Ville Pukki -This version ported from ver 0.99 PD code by Scott Wilson +This version ported from ver 1.0.3.2 PD code by Scott Wilson Development funded in part by the AHRC http://www.ahrc.ac.uk Copyright @@ -27,7 +27,7 @@ use or for distribution: VBAPSpeakerArray { - classvar <>maxNumSpeakers = 55, minSideLength = 0.01; + classvar minSideLength = 0.01; var 0.00001, { ^(volper / lgth)}, { ^0.0 }); } //unq_cross_prod(t_ls v1,t_ls v2, t_ls *res) /* vector cross product */ unq_cross_prod { |v1, v2| - var length, result; - result = VBAPSpeaker.new; - result.x = (v1.y * v2.z ) - (v1.z * v2.y); - result.y = (v1.z * v2.x ) - (v1.x * v2.z); + var length, result; + result = VBAPSpeaker.new; + result.x = (v1.y * v2.z ) - (v1.z * v2.y); + result.y = (v1.z * v2.x ) - (v1.x * v2.z); result.z = (v1.x * v2.y ) - (v1.y * v2.x); length = this.vec_length(result); result.x = result.x / length; @@ -262,17 +267,17 @@ VBAPSpeakerArray { vec_angle{ |v1, v2| /* angle between two loudspeakers */ - var inner; - inner = ((v1.x*v2.x) + (v1.y*v2.y) + (v1.z*v2.z)) / - (this.vec_length(v1) * this.vec_length(v2)); + var inner; + inner = ((v1.x*v2.x) + (v1.y*v2.y) + (v1.z*v2.z)) / + (this.vec_length(v1) * this.vec_length(v2)); if(inner > 1.0, {inner = 1.0}); if (inner < -1.0, {inner = -1.0}); ^abs(acos(inner)); } any_ls_inside_triplet { |a, b, c| // speakers, numSpeakers - /* returns true if there is loudspeaker(s) inside given ls triplet */ - var invdet; + /* returns true if there is loudspeaker(s) inside given ls triplet */ + var invdet; var lp1, lp2, lp3; var invmx; var tmp; @@ -309,9 +314,9 @@ VBAPSpeakerArray { tmp = speakers[i].z * invmx[2 + (j*3)] + tmp; if(tmp < -0.001, {this_inside = false;}); }); - if(this_inside, {any_ls_inside = true}); - }); - }); + if(this_inside, {any_ls_inside = true}); + }); + }); ^any_ls_inside; } @@ -325,12 +330,12 @@ VBAPSpeakerArray { var result; if(sets.isNil, { - postln("define-loudspeakers: Not valid 3-D configuration"); - ^nil; + postln("define-loudspeakers: Not valid 3-D configuration"); + ^nil; }); - triplet_amount = sets.size; - //"triplet_amount: %\n".postf(triplet_amount); + triplet_amount = sets.size; + //"triplet_amount: %\n".postf(triplet_amount); list_length = triplet_amount * 21 + 2; result = FloatArray.newClear(list_length); @@ -339,7 +344,7 @@ VBAPSpeakerArray { pointer=2; sets.do({|set| - lp1 = speakers[set.chanOffsets[0]]; + lp1 = speakers[set.chanOffsets[0]]; lp2 = speakers[set.chanOffsets[1]]; lp3 = speakers[set.chanOffsets[2]]; @@ -389,22 +394,25 @@ VBAPSpeakerArray { choose_ls_tuplets { /* selects the loudspeaker pairs, calculates the inversion matrices and stores the data to a global array*/ - var atorad = (2 * 3.1415927 / 360); - var w1,w2; - var p1,p2; + var atorad = (2 * pi / 360); + //var w1,w2; + //var p1,p2; var sorted_lss; var exist; var amount=0; var inv_mat; - var ls_table; + var mat; + //var ls_table; var list_length; var result; var pointer; - exist = Array.newClear(maxNumSpeakers); - inv_mat = Array.fill(maxNumSpeakers, {Array.newClear(4)}); + exist = Array.newClear(numSpeakers); + inv_mat = Array.fill(numSpeakers, {Array.newClear(4)}); + // not sure + mat = Array.fill(numSpeakers, {Array.newClear(4)}); - for(0, maxNumSpeakers - 1, {|i| + for(0, numSpeakers - 1, {|i| exist[i]=0; }); @@ -415,25 +423,25 @@ VBAPSpeakerArray { for(0, numSpeakers -2, {|i| if((speakers[sorted_lss[i+1]].azi - speakers[sorted_lss[i]].azi) <= (180 - 10), { if(this.calc_2D_inv_tmatrix(speakers[sorted_lss[i]].azi, - speakers[sorted_lss[i+1]].azi, inv_mat[i]),{ + speakers[sorted_lss[i+1]].azi, inv_mat[i], mat[i]),{ exist[i]=1; amount = amount + 1; }); }); }); - if(((6.283 - speakers[sorted_lss[numSpeakers-1]].azi) + if(((360 - speakers[sorted_lss[numSpeakers-1]].azi) + speakers[sorted_lss[0]].azi) <= (180 - 10), { if(this.calc_2D_inv_tmatrix(speakers[sorted_lss[numSpeakers-1]].azi, speakers[sorted_lss[0]].azi, - inv_mat[numSpeakers-1]), { + inv_mat[numSpeakers-1], mat[numSpeakers-1]), { exist[numSpeakers-1]=1; amount = amount + 1; }); }); /* Output */ - list_length= amount * 6 + 2; + list_length= amount * 10 + 2; result = Array.newClear(list_length); result[0] = dim; @@ -450,6 +458,10 @@ VBAPSpeakerArray { result[pointer] = inv_mat[i][j]; pointer = pointer + 1; }); + for(0, 3, {|j| + result[pointer] = mat[i][j]; + pointer = pointer + 1; + }); }); }); if(exist[numSpeakers-1] == 1, { @@ -461,6 +473,10 @@ VBAPSpeakerArray { result[pointer] = inv_mat[numSpeakers-1][j]; pointer = pointer + 1; }); + for(0, 3, {|j| + result[pointer] = mat[numSpeakers-1][j]; + pointer = pointer + 1; + }); }); ^result; } @@ -468,14 +484,14 @@ VBAPSpeakerArray { sort_2D_lss { /* sort loudspeakers according to azimuth angle */ - var i,j,index; + var i,j; var tmp, tmp_azi; var rad2ang = 360.0 / ( 2 * pi ); - var x,y; + //var x,y; var sorted_lss; - sorted_lss = Array.newClear(maxNumSpeakers); + sorted_lss = Array.newClear(numSpeakers); /* Transforming angles between -180 and 180 */ for (0, numSpeakers - 1, {|i| @@ -488,7 +504,9 @@ VBAPSpeakerArray { speakers[i].azi = speakers[i].azi * tmp; }); for (0, numSpeakers - 1, {|i| + var index = 0; tmp = 2000; + for (0, numSpeakers - 1, {|j| if (speakers[j].azi <= tmp, { tmp = speakers[j].azi; @@ -506,18 +524,18 @@ VBAPSpeakerArray { ^sorted_lss; } - calc_2D_inv_tmatrix { |azi1, azi2, inv_mat| + calc_2D_inv_tmatrix { |azi1, azi2, inv_mat, mat| /* calculate inverse 2x2 matrix */ var x1,x2,x3,x4; - var y1,y2,y3,y4; + //var y1,y2,y3,y4; var det; - var rad2ang = 360.0 / ( 2 * 3.141592 ); + var rad2ang = 360.0 / ( 2 * pi ); - x1 = cos(azi1 / rad2ang); - x2 = sin(azi1 / rad2ang); - x3 = cos(azi2 / rad2ang); - x4 = sin(azi2 / rad2ang); + mat[0] = x1 = cos(azi1 / rad2ang); + mat[1] = x2 = sin(azi1 / rad2ang); + mat[2] = x3 = cos(azi2 / rad2ang); + mat[3] = x4 = sin(azi2 / rad2ang); det = (x1 * x4) - ( x3 * x2 ); if(abs(det) <= 0.001, { inv_mat[0] = 0.0; @@ -562,8 +580,8 @@ VBAPSpeakerSet { // triplet or pair var <>inv_mx; *new {|chanOffsets| - ^super.newCopyArgs(chanOffsets); - } + ^super.newCopyArgs(chanOffsets); + } } VBAP : MultiOutUGen {