LCOV - code coverage report
Current view: top level - lib_rend - ivas_vbap.c (source / functions) Hit Total Coverage
Test: Coverage on main -- short test vectors @ 6c9ddc4024a9c0e1ecb8f643f114a84a0e26ec6b Lines: 565 625 90.4 %
Date: 2025-05-23 08:37:30 Functions: 21 21 100.0 %

          Line data    Source code
       1             : /******************************************************************************************************
       2             : 
       3             :    (C) 2022-2025 IVAS codec Public Collaboration with portions copyright Dolby International AB, Ericsson AB,
       4             :    Fraunhofer-Gesellschaft zur Foerderung der angewandten Forschung e.V., Huawei Technologies Co. LTD.,
       5             :    Koninklijke Philips N.V., Nippon Telegraph and Telephone Corporation, Nokia Technologies Oy, Orange,
       6             :    Panasonic Holdings Corporation, Qualcomm Technologies, Inc., VoiceAge Corporation, and other
       7             :    contributors to this repository. All Rights Reserved.
       8             : 
       9             :    This software is protected by copyright law and by international treaties.
      10             :    The IVAS codec Public Collaboration consisting of Dolby International AB, Ericsson AB,
      11             :    Fraunhofer-Gesellschaft zur Foerderung der angewandten Forschung e.V., Huawei Technologies Co. LTD.,
      12             :    Koninklijke Philips N.V., Nippon Telegraph and Telephone Corporation, Nokia Technologies Oy, Orange,
      13             :    Panasonic Holdings Corporation, Qualcomm Technologies, Inc., VoiceAge Corporation, and other
      14             :    contributors to this repository retain full ownership rights in their respective contributions in
      15             :    the software. This notice grants no license of any kind, including but not limited to patent
      16             :    license, nor is any license granted by implication, estoppel or otherwise.
      17             : 
      18             :    Contributors are required to enter into the IVAS codec Public Collaboration agreement before making
      19             :    contributions.
      20             : 
      21             :    This software is provided "AS IS", without any express or implied warranties. The software is in the
      22             :    development stage. It is intended exclusively for experts who have experience with such software and
      23             :    solely for the purpose of inspection. All implied warranties of non-infringement, merchantability
      24             :    and fitness for a particular purpose are hereby disclaimed and excluded.
      25             : 
      26             :    Any dispute, controversy or claim arising under or in relation to providing this software shall be
      27             :    submitted to and settled by the final, binding jurisdiction of the courts of Munich, Germany in
      28             :    accordance with the laws of the Federal Republic of Germany excluding its conflict of law rules and
      29             :    the United Nations Convention on Contracts on the International Sales of Goods.
      30             : 
      31             : *******************************************************************************************************/
      32             : 
      33             : #include <stdint.h>
      34             : #include "options.h"
      35             : #include <assert.h>
      36             : #include <math.h>
      37             : #include "prot.h"
      38             : #include "ivas_prot.h"
      39             : #include "ivas_prot_rend.h"
      40             : #include "ivas_stat_dec.h"
      41             : #ifdef DEBUGGING
      42             : #include "debug.h"
      43             : #endif
      44             : #include "wmc_auto.h"
      45             : 
      46             : /*-----------------------------------------------------------------------*
      47             :  * Local constants
      48             :  *-----------------------------------------------------------------------*/
      49             : 
      50             : /* 128 is maximum num_speaker_nodes number. This relates to memory optimization and maximum of triplets:
      51             :  - triplet indices are unsigned_char (see below structs) --> max triplets is 256
      52             :  - num_speaker_nodes_internal = num_speaker_nodes + 2 (potential virtual node channels, bottom and top)
      53             :  - max_num_triplets = 256 = (max_num_ls_internal - 2) * 2 = (max_num_ls) * 2
      54             :  --> max_num_ls = 128
      55             :  */
      56             : #define VBAP_MAX_NUM_SPEAKER_NODES                  128
      57             : #define VBAP_MAX_NUM_TRIPLETS                       256
      58             : #define VBAP_EPSILON                                0.001f /* The fairly large epsilon is for detecting planes etc and accounts for rounding issues */
      59             : #define VBAP_MAX_PLANES                             50
      60             : #define VBAP_MAX_HORIZONTAL_GAP_FOR_PLANE_DETECTION 140.0f
      61             : /* If a speaker node is found
      62             :  - above VBAP_NO_VIRTUAL_SPEAKER_NODE_ELE_LIMIT, no virtual node is used
      63             :  - above VBAP_DISTRIBUTE_VIRTUAL_SPEAKER_NODE_ELE_LIMIT, energy-spreading virtual node is used
      64             :  - not above VBAP_DISTRIBUTE_VIRTUAL_SPEAKER_NODE_ELE_LIMIT, energy-omitting virtual node is used
      65             :  Same applies for both elevations and inclinations, i.e., the two half-spheres. */
      66             : #define VBAP_NO_VIRTUAL_SPEAKER_NODE_ELE_LIMIT         45.0f
      67             : #define VBAP_DISTRIBUTE_VIRTUAL_SPEAKER_NODE_ELE_LIMIT 20.0f
      68             : #define VBAP_VIRTUAL_BACK_ELE_LIMIT                    45.0f
      69             : #define VBAP_NOT_VALID_CONNECTION                      ( -1 )
      70             : 
      71             : /* Maximum azimuth gap between speaker nodes for detecting zero elevation horizontal plane  */
      72             : #define VBAP_MAX_HORIZONTAL_GAP 170u
      73             : 
      74             : #define VBAP_SEARCH_SECTOR_SIZE ( 360.0f / ( VBAP_NUM_SEARCH_SECTORS ) )
      75             : 
      76             : 
      77             : enum VirtualSpeakerNodeType
      78             : {
      79             :     NO_VIRTUAL_SPEAKER_NODE,
      80             :     VIRTUAL_SPEAKER_NODE_DISCARD_ENERGY,
      81             :     VIRTUAL_SPEAKER_NODE_DISTRIBUTE_ENERGY
      82             : };
      83             : 
      84             : enum ConnectionClass
      85             : {
      86             :     REGULAR_CONNECTION,
      87             :     ELEVATED_PLANE_THIN_TRIANGLE_CONNECTION,
      88             :     CONNECTION_WITH_SPEAKER_NODE_BEHIND
      89             : };
      90             : 
      91             : typedef struct connection_option
      92             : {
      93             :     int16_t chA;
      94             :     int16_t chB;
      95             :     float arc;
      96             :     float arc_weighted;
      97             : } ConnectionOption;
      98             : 
      99             : enum SpeakerNodeGroup
     100             : {
     101             :     SPEAKER_NODE_BOTTOM_HALF,
     102             :     SPEAKER_NODE_HORIZONTAL,
     103             :     SPEAKER_NODE_TOP_HALF,
     104             :     SPEAKER_NODE_BACK,
     105             :     SPEAKER_NODE_ALL
     106             : };
     107             : 
     108             : /* Defines a single speaker node */
     109             : typedef struct vbap_speaker_node_structure
     110             : {
     111             :     float azi_deg;
     112             :     float ele_deg;
     113             :     float unit_vec[3];
     114             :     enum SpeakerNodeGroup group;
     115             : 
     116             : } VBAP_SPEAKER_NODE;
     117             : 
     118             : 
     119             : /*-----------------------------------------------------------------------*
     120             :  * Local function prototypes
     121             :  *-----------------------------------------------------------------------*/
     122             : 
     123             : static uint8_t vector_matrix_multiply_3x3( const float *src_vector, float matrix[3][3], float *result );
     124             : 
     125             : static void init_speaker_node_direction_data( VBAP_SPEAKER_NODE *speaker_node_data, const float *speaker_node_azi_deg, const float *speaker_node_ele_deg, const int16_t num_speaker_nodes );
     126             : 
     127             : static int16_t determine_virtual_surface_triplets( const int16_t num_speaker_nodes, const VBAP_SPEAKER_NODE *speaker_node_data, int16_t connections[][2], const int16_t max_num_connections, VBAP_VS_TRIPLET *triplets, int16_t initial_search_indices[VBAP_NUM_SEARCH_SECTORS], enum SpeakerNodeGroup allowed_group );
     128             : 
     129             : static void determine_initial_search_indices( const int16_t num_triplets, const float triplet_azidegs[VBAP_MAX_NUM_TRIPLETS], int16_t initial_search_indices[VBAP_NUM_SEARCH_SECTORS] );
     130             : 
     131             : static ivas_error determine_connections( const int16_t num_speaker_nodes, const VBAP_SPEAKER_NODE *speaker_node_data, int16_t connections[][2], const int16_t max_num_connections, int16_t *group1_count, int16_t *group2_start, int16_t *group2_count );
     132             : 
     133             : static void formulate_horizontal_connections( const VBAP_SPEAKER_NODE *speaker_node_data, const int16_t num_speaker_nodes, int16_t connections[][2], int16_t *connection_write_index );
     134             : 
     135             : static ivas_error get_half_sphere_connection_options( const VBAP_SPEAKER_NODE *speaker_node_data, const enum SpeakerNodeGroup group, const int16_t num_speaker_nodes, const int16_t num_non_crossing_planes, const float *non_crossing_plane_elevation_deg, ConnectionOption **connection_options_pr, int16_t *num_connection_options );
     136             : 
     137             : static ivas_error formulate_half_sphere_connections( const VBAP_SPEAKER_NODE *speaker_node_data, const int16_t num_speaker_nodes, const enum SpeakerNodeGroup group, int16_t connections[][2], int16_t *connection_write_index, const int16_t max_num_connections, const int16_t num_non_crossing_planes, const float *non_crossing_plane_elevation_deg );
     138             : 
     139             : static int16_t determine_non_crossing_planes( const int16_t num_speaker_nodes, const VBAP_SPEAKER_NODE *node_data, float *non_crossing_plane_elevation_deg );
     140             : 
     141             : static enum VirtualSpeakerNodeType check_need_of_virtual_speaker_node( VBAP_HANDLE hVBAPdata, const float *speaker_node_azi_deg, const float *speaker_node_ele_deg, enum SpeakerNodeGroup group );
     142             : 
     143             : static int16_t determine_best_triplet_and_gains( VBAP_SEARCH_STRUCT *search_struct, const float panning_unit_vec[3], const int16_t azi_deg, float gains[3] );
     144             : 
     145             : static void determine_virtual_speaker_node_division_gains( const int16_t virtual_speaker_node_index, float *virtual_node_division_gains, int16_t connections[][2], const enum VirtualSpeakerNodeType type, const int16_t max_num_connections, const int16_t num_speaker_nodes, const int16_t use_object_mode );
     146             : 
     147             : static void reorder_triplets( VBAP_VS_TRIPLET *triplets, const int16_t *target_order, const int16_t num_triplets );
     148             : 
     149             : 
     150             : /*-------------------------------------------------------------------------*
     151             :  * vbap_init_data()
     152             :  *
     153             :  * Initialize VBAP data structure for the speaker node set
     154             :  *-------------------------------------------------------------------------*/
     155             : 
     156        2668 : ivas_error vbap_init_data(
     157             :     VBAP_HANDLE *hVBAPdata,            /* i/o: handle for VBAP data structure that will be initialized    */
     158             :     const float *speaker_node_azi_deg, /* i  : vector of speaker node azimuths (positive left)            */
     159             :     const float *speaker_node_ele_deg, /* i  : vector of speaker node elevations (positive up)            */
     160             :     const int16_t num_speaker_nodes,   /* i  : number of speaker nodes in the set                         */
     161             :     const IVAS_FORMAT ivas_format      /* i  : IVAS format                                                */
     162             : )
     163             : {
     164             :     /* Variables */
     165             :     int16_t connections[VBAP_MAX_NUM_SPEAKER_NODES][2];
     166             :     int16_t max_num_connections;
     167             :     int16_t is_success;
     168             :     int16_t connection_group1_count;
     169             :     int16_t connection_group2_start;
     170             :     int16_t connection_group2_count;
     171             :     enum VirtualSpeakerNodeType virtual_top_type;
     172             :     enum VirtualSpeakerNodeType virtual_bottom_type;
     173             :     enum VirtualSpeakerNodeType virtual_back_type;
     174             :     float speaker_node_azi_deg_internal[VBAP_MAX_NUM_SPEAKER_NODES];
     175             :     float speaker_node_ele_deg_internal[VBAP_MAX_NUM_SPEAKER_NODES];
     176             :     VBAP_SPEAKER_NODE speaker_node_data[VBAP_MAX_NUM_SPEAKER_NODES];
     177             :     VBAP_DATA *vbap;
     178             :     ivas_error error;
     179             : 
     180        2668 :     push_wmops( "vbap_init" );
     181             : 
     182             :     /* Basic init checks */
     183             :     /* If the requested layout is invalid, hVBAPdata is set to NULL and the signal will
     184             :      * be distributed with an equal gain into all output channels.
     185             :      * The surrounding code needs to handle the NULL pointer properly. */
     186        2668 :     if ( num_speaker_nodes > VBAP_MAX_NUM_SPEAKER_NODES || num_speaker_nodes < 3 )
     187             :     {
     188           0 :         hVBAPdata = NULL;
     189           0 :         pop_wmops();
     190           0 :         return IVAS_ERR_OK;
     191             :     }
     192             : 
     193        2668 :     if ( !speaker_node_azi_deg || !speaker_node_ele_deg )
     194             :     {
     195           0 :         hVBAPdata = NULL;
     196           0 :         return IVAS_ERR_OK;
     197             :     }
     198             : 
     199             :     /* Allocate VBAP structure */
     200        2668 :     if ( ( vbap = (VBAP_HANDLE) malloc( sizeof( VBAP_DATA ) ) ) == NULL )
     201             :     {
     202           0 :         return ( IVAS_ERROR( IVAS_ERR_FAILED_ALLOC, "Can not allocate memory for VBAP data\n" ) );
     203             :     }
     204             : 
     205        2668 :     is_success = 1;
     206        2668 :     vbap->bottom_virtual_speaker_node_index = -1;
     207        2668 :     vbap->top_virtual_speaker_node_index = -1;
     208        2668 :     vbap->back_virtual_speaker_node_index = -1;
     209        2668 :     vbap->bottom_virtual_speaker_node_division_gains = NULL;
     210        2668 :     vbap->top_virtual_speaker_node_division_gains = NULL;
     211        2668 :     vbap->back_virtual_speaker_node_division_gains = NULL;
     212        2668 :     vbap->object_mode_bottom_virtual_speaker_node_division_gains = NULL;
     213        2668 :     vbap->object_mode_top_virtual_speaker_node_division_gains = NULL;
     214        2668 :     vbap->object_mode_back_virtual_speaker_node_division_gains = NULL;
     215        2668 :     vbap->num_speaker_nodes = num_speaker_nodes;
     216        2668 :     vbap->num_speaker_nodes_internal = num_speaker_nodes;
     217             : 
     218             :     /* Check if the speaker node setup needs a virtual top or bottom node
     219             :          (function also increments vbap->num_speaker_nodes_internal when necessary) */
     220        2668 :     virtual_bottom_type = check_need_of_virtual_speaker_node( vbap, speaker_node_azi_deg, speaker_node_ele_deg, SPEAKER_NODE_BOTTOM_HALF );
     221        2668 :     virtual_top_type = check_need_of_virtual_speaker_node( vbap, speaker_node_azi_deg, speaker_node_ele_deg, SPEAKER_NODE_TOP_HALF );
     222        2668 :     virtual_back_type = check_need_of_virtual_speaker_node( vbap, speaker_node_azi_deg, speaker_node_ele_deg, SPEAKER_NODE_BACK );
     223             : 
     224             :     /* Init internal speaker node configuration, which is the original configuration
     225             :          potentially appended with virtual top and/or bottom loudspeakers */
     226        2668 :     mvr2r( speaker_node_azi_deg, speaker_node_azi_deg_internal, num_speaker_nodes );
     227        2668 :     mvr2r( speaker_node_ele_deg, speaker_node_ele_deg_internal, num_speaker_nodes );
     228             : 
     229        2668 :     if ( is_success && virtual_bottom_type != NO_VIRTUAL_SPEAKER_NODE )
     230             :     {
     231        2668 :         if ( ( vbap->bottom_virtual_speaker_node_division_gains = (float *) malloc( num_speaker_nodes * sizeof( float ) ) ) == NULL )
     232             :         {
     233           0 :             return ( IVAS_ERROR( IVAS_ERR_FAILED_ALLOC, "Can not allocate memory for VBAP data\n" ) );
     234             :         }
     235        2668 :         set_zero( vbap->bottom_virtual_speaker_node_division_gains, num_speaker_nodes );
     236        2668 :         is_success &= vbap->bottom_virtual_speaker_node_division_gains != NULL;
     237             : 
     238        2668 :         if ( ivas_format == MASA_ISM_FORMAT )
     239             :         {
     240        1017 :             if ( ( vbap->object_mode_bottom_virtual_speaker_node_division_gains = (float *) malloc( num_speaker_nodes * sizeof( float ) ) ) == NULL )
     241             :             {
     242           0 :                 return ( IVAS_ERROR( IVAS_ERR_FAILED_ALLOC, "Can not allocate memory for VBAP data\n" ) );
     243             :             }
     244        1017 :             set_zero( vbap->object_mode_bottom_virtual_speaker_node_division_gains, num_speaker_nodes );
     245        1017 :             is_success &= vbap->object_mode_bottom_virtual_speaker_node_division_gains != NULL;
     246             :         }
     247             : 
     248        2668 :         speaker_node_azi_deg_internal[vbap->bottom_virtual_speaker_node_index] = 0.0f;
     249        2668 :         speaker_node_ele_deg_internal[vbap->bottom_virtual_speaker_node_index] = -90.0f;
     250             :     }
     251             : 
     252        2668 :     if ( is_success && virtual_top_type != NO_VIRTUAL_SPEAKER_NODE )
     253             :     {
     254        2668 :         if ( ( vbap->top_virtual_speaker_node_division_gains = (float *) malloc( num_speaker_nodes * sizeof( float ) ) ) == NULL )
     255             :         {
     256           0 :             return ( IVAS_ERROR( IVAS_ERR_FAILED_ALLOC, "Can not allocate memory for VBAP data\n" ) );
     257             :         }
     258        2668 :         set_zero( vbap->top_virtual_speaker_node_division_gains, num_speaker_nodes );
     259        2668 :         is_success &= vbap->top_virtual_speaker_node_division_gains != NULL;
     260             : 
     261        2668 :         if ( ivas_format == MASA_ISM_FORMAT )
     262             :         {
     263        1017 :             if ( ( vbap->object_mode_top_virtual_speaker_node_division_gains = (float *) malloc( num_speaker_nodes * sizeof( float ) ) ) == NULL )
     264             :             {
     265           0 :                 return ( IVAS_ERROR( IVAS_ERR_FAILED_ALLOC, "Can not allocate memory for VBAP data\n" ) );
     266             :             }
     267        1017 :             set_zero( vbap->object_mode_top_virtual_speaker_node_division_gains, num_speaker_nodes );
     268        1017 :             is_success &= vbap->object_mode_top_virtual_speaker_node_division_gains != NULL;
     269             :         }
     270             : 
     271        2668 :         speaker_node_azi_deg_internal[vbap->top_virtual_speaker_node_index] = 0.0f;
     272        2668 :         speaker_node_ele_deg_internal[vbap->top_virtual_speaker_node_index] = 90.0f;
     273             :     }
     274             : 
     275        2668 :     if ( is_success && virtual_back_type != NO_VIRTUAL_SPEAKER_NODE )
     276             :     {
     277           0 :         if ( ( vbap->back_virtual_speaker_node_division_gains = (float *) malloc( num_speaker_nodes * sizeof( float ) ) ) == NULL )
     278             :         {
     279           0 :             return ( IVAS_ERROR( IVAS_ERR_FAILED_ALLOC, "Can not allocate memory for VBAP data\n" ) );
     280             :         }
     281           0 :         set_zero( vbap->back_virtual_speaker_node_division_gains, num_speaker_nodes );
     282           0 :         is_success &= vbap->back_virtual_speaker_node_division_gains != NULL;
     283             : 
     284           0 :         if ( ivas_format == MASA_ISM_FORMAT )
     285             :         {
     286           0 :             if ( ( vbap->object_mode_back_virtual_speaker_node_division_gains = (float *) malloc( num_speaker_nodes * sizeof( float ) ) ) == NULL )
     287             :             {
     288           0 :                 return ( IVAS_ERROR( IVAS_ERR_FAILED_ALLOC, "Can not allocate memory for VBAP data\n" ) );
     289             :             }
     290           0 :             set_zero( vbap->object_mode_back_virtual_speaker_node_division_gains, num_speaker_nodes );
     291           0 :             is_success &= vbap->object_mode_back_virtual_speaker_node_division_gains != NULL;
     292             :         }
     293           0 :         speaker_node_azi_deg_internal[vbap->back_virtual_speaker_node_index] = 180.0f;
     294           0 :         speaker_node_ele_deg_internal[vbap->back_virtual_speaker_node_index] = 0.0f;
     295             :     }
     296             : 
     297        2668 :     init_speaker_node_direction_data( speaker_node_data, speaker_node_azi_deg_internal, speaker_node_ele_deg_internal, vbap->num_speaker_nodes_internal );
     298             : 
     299             :     /* Allocate and determine node-node connections */
     300        2668 :     max_num_connections = ( vbap->num_speaker_nodes_internal - 2 ) * 3; /* Theoretical maximum */
     301             : 
     302        2668 :     if ( ( error = determine_connections( vbap->num_speaker_nodes_internal, speaker_node_data, connections, max_num_connections, &connection_group1_count, &connection_group2_start, &connection_group2_count ) ) != IVAS_ERR_OK )
     303             :     {
     304           0 :         return error;
     305             :     }
     306             : 
     307             :     /* Allocate and determine virtual surface speaker node triplets */
     308        2668 :     if ( is_success )
     309             :     {
     310             :         int16_t ch;
     311        2668 :         int16_t speaker_nodes_group1_internal = 0;
     312        2668 :         int16_t speaker_nodes_group2_internal = 0;
     313        2668 :         int16_t speaker_nodes_horiz_internal = 0;
     314        2668 :         uint8_t loop_done = 0;
     315             : 
     316             :         /* Count nodes in different groups to reserve correct memory */
     317       29742 :         for ( ch = 0; ch < vbap->num_speaker_nodes_internal && !loop_done; ch++ )
     318             :         {
     319       27074 :             switch ( speaker_node_data[ch].group )
     320             :             {
     321           0 :                 case SPEAKER_NODE_ALL:
     322             :                     /* If there is even one speaker belonging to "all" group, then all speakers belong to the "all" group.
     323             :                      * We can skip further counts here. */
     324           0 :                     speaker_nodes_group1_internal = vbap->num_speaker_nodes_internal;
     325           0 :                     loop_done = 1;
     326           0 :                     break;
     327        2680 :                 case SPEAKER_NODE_BOTTOM_HALF:
     328        2680 :                     speaker_nodes_group1_internal++;
     329        2680 :                     break;
     330        8940 :                 case SPEAKER_NODE_TOP_HALF:
     331        8940 :                     speaker_nodes_group2_internal++;
     332        8940 :                     break;
     333       15454 :                 case SPEAKER_NODE_HORIZONTAL:
     334             :                 case SPEAKER_NODE_BACK:
     335       15454 :                     speaker_nodes_group1_internal++;
     336       15454 :                     speaker_nodes_group2_internal++;
     337       15454 :                     speaker_nodes_horiz_internal++;
     338       15454 :                     break;
     339             :             }
     340       27074 :         }
     341             : 
     342        2668 :         if ( ( vbap->search_struct[0].triplets = (VBAP_VS_TRIPLET *) malloc( ( ( speaker_nodes_group1_internal - 2 ) * 2 - ( max( 0, ( speaker_nodes_horiz_internal - 2 ) ) ) ) * sizeof( VBAP_VS_TRIPLET ) ) ) == NULL )
     343             :         {
     344           0 :             return ( IVAS_ERROR( IVAS_ERR_FAILED_ALLOC, "Can not allocate memory for VBAP data\n" ) );
     345             :         }
     346        2668 :         is_success &= vbap->search_struct[0].triplets != NULL;
     347             : 
     348        2668 :         if ( speaker_nodes_group2_internal > 0 )
     349             :         {
     350        2668 :             vbap->num_search_structs = 2;
     351        2668 :             if ( ( vbap->search_struct[1].triplets = (VBAP_VS_TRIPLET *) malloc( ( ( speaker_nodes_group2_internal - 2 ) * 2 - ( max( 0, ( speaker_nodes_horiz_internal - 2 ) ) ) ) * sizeof( VBAP_VS_TRIPLET ) ) ) == NULL )
     352             :             {
     353           0 :                 return ( IVAS_ERROR( IVAS_ERR_FAILED_ALLOC, "Can not allocate memory for VBAP data\n" ) );
     354             :             }
     355        2668 :             is_success &= vbap->search_struct[1].triplets != NULL;
     356             :         }
     357             :         else
     358             :         {
     359           0 :             vbap->num_search_structs = 1;
     360           0 :             vbap->search_struct[1].triplets = NULL;
     361             :         }
     362             :     }
     363             : 
     364        2668 :     if ( is_success )
     365             :     {
     366        2668 :         if ( vbap->num_search_structs == 1 )
     367             :         {
     368             :             /* If all speaker nodes belong to ALL set, then we only create one triplet set and search structure */
     369           0 :             vbap->search_struct[0].num_triplets = determine_virtual_surface_triplets( vbap->num_speaker_nodes_internal, speaker_node_data, connections, max_num_connections, vbap->search_struct[0].triplets, vbap->search_struct[0].initial_search_indices, SPEAKER_NODE_ALL );
     370             :         }
     371             :         else
     372             :         {
     373             :             /* Otherwise, we have two sets and can handle them separately for more opmitized processing. */
     374        2668 :             vbap->search_struct[0].num_triplets = determine_virtual_surface_triplets( vbap->num_speaker_nodes_internal, speaker_node_data, connections, connection_group1_count, vbap->search_struct[0].triplets, vbap->search_struct[0].initial_search_indices, SPEAKER_NODE_BOTTOM_HALF );
     375        2668 :             vbap->search_struct[1].num_triplets = determine_virtual_surface_triplets( vbap->num_speaker_nodes_internal, speaker_node_data, connections + connection_group2_start, connection_group2_count, vbap->search_struct[1].triplets, vbap->search_struct[1].initial_search_indices, SPEAKER_NODE_TOP_HALF );
     376             :         }
     377             :     }
     378             : 
     379             :     /* Determine how the virtual node gains should be distributed to real nodes, if necessary (checked within function). */
     380        2668 :     if ( is_success )
     381             :     {
     382        2668 :         determine_virtual_speaker_node_division_gains( vbap->top_virtual_speaker_node_index, vbap->top_virtual_speaker_node_division_gains, connections, virtual_top_type, max_num_connections, num_speaker_nodes, 0 );
     383        2668 :         determine_virtual_speaker_node_division_gains( vbap->bottom_virtual_speaker_node_index, vbap->bottom_virtual_speaker_node_division_gains, connections, virtual_bottom_type, max_num_connections, num_speaker_nodes, 0 );
     384        2668 :         determine_virtual_speaker_node_division_gains( vbap->back_virtual_speaker_node_index, vbap->back_virtual_speaker_node_division_gains, connections, virtual_back_type, max_num_connections, num_speaker_nodes, 0 );
     385        2668 :         if ( ivas_format == MASA_ISM_FORMAT )
     386             :         {
     387        1017 :             determine_virtual_speaker_node_division_gains( vbap->top_virtual_speaker_node_index, vbap->object_mode_top_virtual_speaker_node_division_gains, connections, virtual_top_type == NO_VIRTUAL_SPEAKER_NODE ? NO_VIRTUAL_SPEAKER_NODE : VIRTUAL_SPEAKER_NODE_DISTRIBUTE_ENERGY, max_num_connections, num_speaker_nodes, 1 );
     388        1017 :             determine_virtual_speaker_node_division_gains( vbap->bottom_virtual_speaker_node_index, vbap->object_mode_bottom_virtual_speaker_node_division_gains, connections, virtual_bottom_type == NO_VIRTUAL_SPEAKER_NODE ? NO_VIRTUAL_SPEAKER_NODE : VIRTUAL_SPEAKER_NODE_DISTRIBUTE_ENERGY, max_num_connections, num_speaker_nodes, 1 );
     389        1017 :             determine_virtual_speaker_node_division_gains( vbap->back_virtual_speaker_node_index, vbap->object_mode_back_virtual_speaker_node_division_gains, connections, virtual_back_type == NO_VIRTUAL_SPEAKER_NODE ? NO_VIRTUAL_SPEAKER_NODE : VIRTUAL_SPEAKER_NODE_DISTRIBUTE_ENERGY, max_num_connections, num_speaker_nodes, 1 );
     390             :         }
     391             :     }
     392             : 
     393        2668 :     pop_wmops();
     394             : 
     395        2668 :     if ( is_success )
     396             :     {
     397        2668 :         *hVBAPdata = vbap;
     398             :     }
     399             :     else
     400             :     {
     401           0 :         vbap_free_data( &vbap );
     402             :     }
     403             : 
     404        2668 :     return IVAS_ERR_OK;
     405             : }
     406             : 
     407             : 
     408             : /*-------------------------------------------------------------------------*
     409             :  * vbap_free_data()
     410             :  *
     411             :  * Free VBAP data structure
     412             :  *-------------------------------------------------------------------------*/
     413             : 
     414        6799 : void vbap_free_data(
     415             :     VBAP_HANDLE *hVBAPdata /* i/o: VBAP handle to be freed    */
     416             : )
     417             : {
     418        6799 :     if ( hVBAPdata == NULL || *hVBAPdata == NULL )
     419             :     {
     420        4131 :         return;
     421             :     }
     422             : 
     423        2668 :     if ( ( *hVBAPdata )->bottom_virtual_speaker_node_division_gains != NULL )
     424             :     {
     425        2668 :         free( ( *hVBAPdata )->bottom_virtual_speaker_node_division_gains );
     426             :     }
     427        2668 :     if ( ( *hVBAPdata )->top_virtual_speaker_node_division_gains != NULL )
     428             :     {
     429        2668 :         free( ( *hVBAPdata )->top_virtual_speaker_node_division_gains );
     430             :     }
     431        2668 :     if ( ( *hVBAPdata )->back_virtual_speaker_node_division_gains != NULL )
     432             :     {
     433           0 :         free( ( *hVBAPdata )->back_virtual_speaker_node_division_gains );
     434             :     }
     435        2668 :     if ( ( *hVBAPdata )->object_mode_bottom_virtual_speaker_node_division_gains != NULL )
     436             :     {
     437        1017 :         free( ( *hVBAPdata )->object_mode_bottom_virtual_speaker_node_division_gains );
     438             :     }
     439        2668 :     if ( ( *hVBAPdata )->object_mode_top_virtual_speaker_node_division_gains != NULL )
     440             :     {
     441        1017 :         free( ( *hVBAPdata )->object_mode_top_virtual_speaker_node_division_gains );
     442             :     }
     443        2668 :     if ( ( *hVBAPdata )->object_mode_back_virtual_speaker_node_division_gains != NULL )
     444             :     {
     445           0 :         free( ( *hVBAPdata )->object_mode_back_virtual_speaker_node_division_gains );
     446             :     }
     447        2668 :     if ( ( *hVBAPdata )->search_struct[0].triplets != NULL )
     448             :     {
     449        2668 :         free( ( *hVBAPdata )->search_struct[0].triplets );
     450             :     }
     451        2668 :     if ( ( *hVBAPdata )->num_search_structs == 2 && ( *hVBAPdata )->search_struct[1].triplets != NULL )
     452             :     {
     453        2668 :         free( ( *hVBAPdata )->search_struct[1].triplets );
     454             :     }
     455             : 
     456        2668 :     free( *hVBAPdata );
     457        2668 :     *hVBAPdata = NULL;
     458             : 
     459        2668 :     return;
     460             : }
     461             : 
     462             : 
     463             : /*-------------------------------------------------------------------------*
     464             :  * vbap_determine_gains()
     465             :  *
     466             :  * Obtain panning gains for all speaker nodes based on the given direction
     467             :  *-------------------------------------------------------------------------*/
     468             : 
     469    14347219 : void vbap_determine_gains(
     470             :     const VBAP_HANDLE hVBAPdata,  /* i  : prepared VBAP structure                                       */
     471             :     float *gains,                 /* o  : gain vector for loudspeakers for given direction              */
     472             :     const int16_t azi_deg,        /* i  : azimuth in degrees for panning direction (positive left)      */
     473             :     const int16_t ele_deg,        /* i  : elevation in degrees for panning direction (positive up)      */
     474             :     const int16_t use_object_mode /* i  : select between object mode panning and spatial mode panning   */
     475             : )
     476             : {
     477             :     /* This function formulates gains for the given angle. The triplet-selection has been pre-formulated. */
     478             :     int16_t ch, ch2;
     479             :     int16_t triplet_ch;
     480             :     int16_t triplet_index;
     481             :     float panning_unit_vec[3];
     482             :     float gain_triplet[3];
     483             :     float norm_value;
     484             :     float gain_ene;
     485             :     float azi_rad;
     486             :     float ele_rad;
     487             :     float azi_temp;
     488             :     float ele_temp;
     489             :     int16_t num_speaker_nodes;
     490             :     int16_t bottom_virtual_speaker_node_index;
     491             :     int16_t top_virtual_speaker_node_index;
     492             :     int16_t back_virtual_speaker_node_index;
     493             :     VBAP_VS_TRIPLET *selected_triplet;
     494             :     float *bottom_virtual_speaker_node_division_gains;
     495             :     float *top_virtual_speaker_node_division_gains;
     496             :     float *back_virtual_speaker_node_division_gains;
     497             : 
     498             : #ifdef DEBUGGING
     499             :     assert( hVBAPdata != NULL && "VBAP gain determination requires initialized structure." );
     500             :     assert( gains != NULL && "VBAP gain determination requires reserved memory for gain output." );
     501             : #endif
     502             : 
     503    14347219 :     push_wmops( "vbap_gains" );
     504    14347219 :     num_speaker_nodes = hVBAPdata->num_speaker_nodes;
     505    14347219 :     bottom_virtual_speaker_node_index = hVBAPdata->bottom_virtual_speaker_node_index;
     506    14347219 :     top_virtual_speaker_node_index = hVBAPdata->top_virtual_speaker_node_index;
     507    14347219 :     back_virtual_speaker_node_index = hVBAPdata->back_virtual_speaker_node_index;
     508    14347219 :     if ( use_object_mode )
     509             :     {
     510       98637 :         bottom_virtual_speaker_node_division_gains = hVBAPdata->object_mode_bottom_virtual_speaker_node_division_gains;
     511       98637 :         top_virtual_speaker_node_division_gains = hVBAPdata->object_mode_top_virtual_speaker_node_division_gains;
     512       98637 :         back_virtual_speaker_node_division_gains = hVBAPdata->object_mode_back_virtual_speaker_node_division_gains;
     513             :     }
     514             :     else
     515             :     {
     516    14248582 :         bottom_virtual_speaker_node_division_gains = hVBAPdata->bottom_virtual_speaker_node_division_gains;
     517    14248582 :         top_virtual_speaker_node_division_gains = hVBAPdata->top_virtual_speaker_node_division_gains;
     518    14248582 :         back_virtual_speaker_node_division_gains = hVBAPdata->back_virtual_speaker_node_division_gains;
     519             :     }
     520             : 
     521    14347219 :     panning_wrap_angles( (float) azi_deg, (float) ele_deg, &azi_temp, &ele_temp );
     522    14347219 :     azi_rad = azi_temp * PI_OVER_180;
     523    14347219 :     ele_rad = ele_temp * PI_OVER_180;
     524             : 
     525    14347219 :     panning_unit_vec[0] = cosf( azi_rad ) * cosf( ele_rad );
     526    14347219 :     panning_unit_vec[1] = sinf( azi_rad ) * cosf( ele_rad );
     527    14347219 :     panning_unit_vec[2] = sinf( ele_rad );
     528             : 
     529             :     /* Find the best VS triplet and speaker node gains for the panning direction using the prepared search structures. */
     530    14347219 :     if ( hVBAPdata->num_search_structs == 2 && ele_deg > 0 )
     531             :     {
     532     6112842 :         triplet_index = determine_best_triplet_and_gains( &( hVBAPdata->search_struct[1] ), panning_unit_vec, azi_deg, gain_triplet );
     533     6112842 :         selected_triplet = &hVBAPdata->search_struct[1].triplets[triplet_index];
     534             :     }
     535             :     else
     536             :     {
     537     8234377 :         triplet_index = determine_best_triplet_and_gains( &( hVBAPdata->search_struct[0] ), panning_unit_vec, azi_deg, gain_triplet );
     538     8234377 :         selected_triplet = &hVBAPdata->search_struct[0].triplets[triplet_index];
     539             :     }
     540             : 
     541             :     /* Normalize to unit energy */
     542    14347219 :     gain_ene = 1e-12f; /* Add small value to avoid divide by zero. */
     543    57388876 :     for ( ch = 0; ch < 3; ch++ )
     544             :     {
     545    43041657 :         gain_ene += gain_triplet[ch] * gain_triplet[ch];
     546             :     }
     547             : 
     548    14347219 :     norm_value = inv_sqrt( gain_ene );
     549             : 
     550    57388876 :     for ( ch = 0; ch < 3; ch++ )
     551             :     {
     552    43041657 :         gain_triplet[ch] *= norm_value;
     553             : 
     554             :         /* Sanity check for rounding issues */
     555    43041657 :         if ( gain_triplet[ch] < 0.0f )
     556             :         {
     557        1094 :             gain_triplet[ch] = 0.0f;
     558             :         }
     559             :     }
     560             : 
     561             :     /* Flush gain target */
     562    14347219 :     set_zero( gains, num_speaker_nodes );
     563             : 
     564             :     /* Map gain triplet (internal speaker node configuration) to speaker node output (actual speaker node configuration) */
     565    57388876 :     for ( ch = 0; ch < 3; ch++ )
     566             :     {
     567    43041657 :         triplet_ch = selected_triplet->speaker_node[ch];
     568             : 
     569    43041657 :         if ( triplet_ch == bottom_virtual_speaker_node_index )
     570             :         {
     571    60325609 :             for ( ch2 = 0; ch2 < num_speaker_nodes; ch2++ )
     572             :             {
     573    53334666 :                 gains[ch2] += bottom_virtual_speaker_node_division_gains[ch2] * gain_triplet[ch];
     574             :             }
     575             :         }
     576    36050714 :         else if ( triplet_ch == top_virtual_speaker_node_index )
     577             :         {
     578    33720587 :             for ( ch2 = 0; ch2 < num_speaker_nodes; ch2++ )
     579             :             {
     580    29607447 :                 gains[ch2] += top_virtual_speaker_node_division_gains[ch2] * gain_triplet[ch];
     581             :             }
     582             :         }
     583    31937574 :         else if ( triplet_ch == back_virtual_speaker_node_index )
     584             :         {
     585           0 :             for ( ch2 = 0; ch2 < num_speaker_nodes; ch2++ )
     586             :             {
     587           0 :                 gains[ch2] += back_virtual_speaker_node_division_gains[ch2] * gain_triplet[ch];
     588             :             }
     589             :         }
     590             :         else
     591             :         {
     592    31937574 :             gains[triplet_ch] += gain_triplet[ch];
     593             :         }
     594             :     }
     595             : 
     596    14347219 :     pop_wmops();
     597             : 
     598    14347219 :     return;
     599             : }
     600             : 
     601             : 
     602             : /*-----------------------------------------------------------------------*
     603             :  * Local functions
     604             :  *-----------------------------------------------------------------------*/
     605             : 
     606             : /*-------------------------------------------------------------------------*
     607             :  * vbap_crossp()
     608             :  *
     609             :  * 3-by-3 vector cross product
     610             :  *-------------------------------------------------------------------------*/
     611             : 
     612      353078 : static void vbap_crossp(
     613             :     const float *vec1,  /* i  : input vector 1       */
     614             :     const float *vec2,  /* i  : input vector 2       */
     615             :     float *crossProduct /* o  : cross product output */
     616             : )
     617             : {
     618      353078 :     crossProduct[0] = ( vec1[1] * vec2[2] ) - ( vec1[2] * vec2[1] );
     619      353078 :     crossProduct[1] = ( vec1[2] * vec2[0] ) - ( vec1[0] * vec2[2] );
     620      353078 :     crossProduct[2] = ( vec1[0] * vec2[1] ) - ( vec1[1] * vec2[0] );
     621             : 
     622      353078 :     return;
     623             : }
     624             : 
     625             : 
     626             : /*-------------------------------------------------------------------------*
     627             :  * vector_matrix_multiply_3x3()
     628             :  *
     629             :  * 3-by-3 vector multiply with matrix
     630             :  *-------------------------------------------------------------------------*/
     631             : 
     632             : /*! r: Status result if triplet is usable for panning. Allows early exit. */
     633    51074121 : static uint8_t vector_matrix_multiply_3x3(
     634             :     const float *src_vector, /* i  : input vector  */
     635             :     float matrix[3][3],      /* i  : input matrix  */
     636             :     float *result            /* o  : output vector */
     637             : )
     638             : {
     639    51074121 :     result[0] = src_vector[0] * matrix[0][0];
     640    51074121 :     result[0] += src_vector[1] * matrix[1][0];
     641    51074121 :     result[0] += src_vector[2] * matrix[2][0];
     642             : 
     643    51074121 :     if ( result[0] < -0.01f )
     644             :     {
     645     5439845 :         return 0;
     646             :     }
     647             : 
     648    45634276 :     result[1] = src_vector[0] * matrix[0][1];
     649    45634276 :     result[1] += src_vector[1] * matrix[1][1];
     650    45634276 :     result[1] += src_vector[2] * matrix[2][1];
     651             : 
     652    45634276 :     if ( result[1] < -0.01f )
     653             :     {
     654    20980106 :         return 0;
     655             :     }
     656             : 
     657    24654170 :     result[2] = src_vector[0] * matrix[0][2];
     658    24654170 :     result[2] += src_vector[1] * matrix[1][2];
     659    24654170 :     result[2] += src_vector[2] * matrix[2][2];
     660             : 
     661    24654170 :     if ( result[2] < -0.01f )
     662             :     {
     663    10130979 :         return 0;
     664             :     }
     665             : 
     666    14523191 :     return 1;
     667             : }
     668             : 
     669             : 
     670             : /*----------------------------------------------------------------------------------------------*
     671             :  * determine_best_triplet_and_gains()
     672             :  *
     673             :  * Determine the best speaker node triplet and associated gains for panning to defined direction
     674             :  *----------------------------------------------------------------------------------------------*/
     675             : 
     676             : /*! r: triplet id */
     677    14347219 : static int16_t determine_best_triplet_and_gains(
     678             :     VBAP_SEARCH_STRUCT *search_struct, /* i  : VBAP search struct       */
     679             :     const float panning_unit_vec[3],   /* i  : panning unit vector      */
     680             :     const int16_t azi_deg,             /* i  : panning azimuth          */
     681             :     float gains[3]                     /* o  : panning gains            */
     682             : )
     683             : {
     684             :     int16_t i, tr, k;
     685             :     uint8_t triplet_ok;
     686             :     int16_t best_triplet;
     687             :     float best_min_gain;
     688             :     float min_gain_this;
     689             :     float unnormalized_gains[3];
     690             :     int16_t sector;
     691             :     int16_t first_triplet;
     692             :     int16_t jump;
     693             :     int16_t num_triplets;
     694             : 
     695    14347219 :     num_triplets = search_struct->num_triplets;
     696    14347219 :     best_min_gain = -999.9f;
     697    14347219 :     best_triplet = 0;
     698    14347219 :     set_zero( gains, 3 );
     699             : 
     700             :     /* Determine the correct search sector for that target panning direction using an optimized algorithm for
     701             :      * the chosen four sectors. */
     702    14347219 :     if ( abs( azi_deg ) > 90 )
     703             :     {
     704     6759373 :         sector = azi_deg < 0 ? 2 : 1;
     705             :     }
     706             :     else
     707             :     {
     708     7587846 :         sector = azi_deg < 0 ? 3 : 0;
     709             :     }
     710    14347219 :     first_triplet = search_struct->initial_search_indices[sector];
     711             : 
     712    14347219 :     tr = first_triplet;
     713    14347219 :     jump = 1;
     714    50745031 :     for ( i = 0; i < num_triplets; i++ )
     715             :     {
     716    50743937 :         triplet_ok = vector_matrix_multiply_3x3( panning_unit_vec, search_struct->triplets[tr].inverse_matrix, unnormalized_gains );
     717    50743937 :         if ( triplet_ok )
     718             :         {
     719    14523191 :             min_gain_this = min( ( min( unnormalized_gains[0], unnormalized_gains[1] ) ), unnormalized_gains[2] );
     720             : 
     721    14523191 :             if ( min_gain_this > best_min_gain )
     722             :             {
     723    14522544 :                 best_min_gain = min_gain_this;
     724    14522544 :                 best_triplet = tr;
     725    58090176 :                 for ( k = 0; k < 3; k++ )
     726             :                 {
     727    43567632 :                     gains[k] = unnormalized_gains[k];
     728             :                 }
     729    14522544 :                 if ( !( best_min_gain < 0.00f ) )
     730             :                 {
     731    14346125 :                     return best_triplet;
     732             :                 }
     733             :             }
     734             :         }
     735    36397812 :         tr = first_triplet + jump;
     736    36397812 :         if ( tr < 0 )
     737             :         {
     738      200503 :             tr += num_triplets;
     739             :         }
     740    36197309 :         else if ( tr >= num_triplets )
     741             :         {
     742     9800199 :             tr -= num_triplets;
     743             :         }
     744             : 
     745    36397812 :         jump *= -1;
     746    36397812 :         if ( jump > 0 )
     747             :         {
     748    15341015 :             jump += 1;
     749             :         }
     750             :     }
     751             : 
     752        1094 :     return best_triplet;
     753             : }
     754             : 
     755             : /*-------------------------------------------------------------------------*
     756             :  * determine_virtual_speaker_node_division_gains()
     757             :  *
     758             :  * Determines how the virtual node gains are distributed to real nodes
     759             :  *-------------------------------------------------------------------------*/
     760             : 
     761       11055 : static void determine_virtual_speaker_node_division_gains(
     762             :     const int16_t virtual_speaker_node_index, /* i  : virtual speaker node index          */
     763             :     float *virtual_node_division_gains,       /* o  : virtual speaker node division gains */
     764             :     int16_t connections[][2],                 /* i  : vector of all connections           */
     765             :     const enum VirtualSpeakerNodeType type,   /* i  : virtual speaker node typel          */
     766             :     const int16_t max_num_connections,        /* i  : max number of connections           */
     767             :     const int16_t num_speaker_nodes,          /* i  : max number of speaker nodes         */
     768             :     const int16_t use_object_mode             /* i  : use VBAP in object panning mode vs. spatial panning mode */
     769             : )
     770             : {
     771             :     /* When node type is VIRTUAL_SPEAKER_NODE_DISTRIBUTE_ENERGY, the gains of the virtual node
     772             :      are distributed to all neighboring real speaker nodes. An amplitude-division
     773             :      instead of energy division is utilized just in case to avoid excessive emphasis
     774             :      on the coherent distributed sound. */
     775             :     int16_t c, ch;
     776             :     float sum_val;
     777             : 
     778       11055 :     if ( type == VIRTUAL_SPEAKER_NODE_DISTRIBUTE_ENERGY )
     779             :     {
     780      102114 :         for ( c = 0; c < max_num_connections; c++ )
     781             :         {
     782       98208 :             if ( connections[c][0] != VBAP_NOT_VALID_CONNECTION )
     783             :             {
     784       98208 :                 int16_t connection_node = -1;
     785       98208 :                 if ( connections[c][0] == virtual_speaker_node_index )
     786             :                 {
     787           0 :                     connection_node = connections[c][1];
     788             :                 }
     789       98208 :                 else if ( connections[c][1] == virtual_speaker_node_index )
     790             :                 {
     791       17637 :                     connection_node = connections[c][0];
     792             :                 }
     793             : 
     794             :                 /* The second condition allows division gains only to actual loudspeakers */
     795       98208 :                 if ( connection_node >= 0 && ( connection_node < num_speaker_nodes ) )
     796             :                 {
     797       17637 :                     virtual_node_division_gains[connection_node] = 1.0f;
     798             :                 }
     799             :             }
     800             :         }
     801             : 
     802        3906 :         sum_val = 0.0f;
     803       36642 :         for ( ch = 0; ch < num_speaker_nodes; ch++ )
     804             :         {
     805       32736 :             sum_val += virtual_node_division_gains[ch];
     806             :         }
     807             : 
     808       36642 :         for ( ch = 0; ch < num_speaker_nodes; ch++ )
     809             :         {
     810       32736 :             virtual_node_division_gains[ch] /= sum_val;
     811       32736 :             if ( use_object_mode )
     812             :             {
     813       15930 :                 virtual_node_division_gains[ch] = powf( virtual_node_division_gains[ch], 0.8f );
     814             :             }
     815             :         }
     816             :     }
     817             : 
     818       11055 :     return;
     819             : }
     820             : 
     821             : 
     822             : /*-------------------------------------------------------------------------*
     823             :  * check_need_of_virtual_speaker_node()
     824             :  *
     825             :  * Check if virtual speaker node is required
     826             :  *-------------------------------------------------------------------------*/
     827             : 
     828             : /*! r: virtual speaker node type */
     829        8004 : static enum VirtualSpeakerNodeType check_need_of_virtual_speaker_node(
     830             :     VBAP_HANDLE hVBAPdata,             /* i/o: VBAP structure                            */
     831             :     const float *speaker_node_azi_deg, /* i  : vector of speaker node azimuths           */
     832             :     const float *speaker_node_ele_deg, /* i  : vector of speaker node elevations         */
     833             :     enum SpeakerNodeGroup group        /* i  : group of speaker nodes where this belongs */
     834             : )
     835             : {
     836             :     int16_t ch;
     837        8004 :     float max_elevation = 0.0f;
     838             : 
     839             :     /* The following considers if SPEAKER_NODE_BACK virtual speaker is needed */
     840        8004 :     if ( group == SPEAKER_NODE_BACK )
     841             :     {
     842        2668 :         int16_t virtual_back_needed = 1;
     843        2668 :         const float virtual_back_epsilon = -0.0175f; /* Corresponds to approximately 91 degrees, see code below */
     844       10597 :         for ( ch = 0; ch < hVBAPdata->num_speaker_nodes; ch++ )
     845             :         {
     846       10597 :             if ( fabsf( speaker_node_ele_deg[ch] ) < VBAP_VIRTUAL_BACK_ELE_LIMIT )
     847             :             {
     848       10597 :                 if ( cosf( speaker_node_azi_deg[ch] * PI_OVER_180 ) < virtual_back_epsilon )
     849             :                 {
     850        2668 :                     virtual_back_needed = 0;
     851        2668 :                     break;
     852             :                 }
     853             :             }
     854             :         }
     855             : 
     856        2668 :         if ( virtual_back_needed )
     857             :         {
     858           0 :             hVBAPdata->back_virtual_speaker_node_index = hVBAPdata->num_speaker_nodes_internal;
     859           0 :             hVBAPdata->num_speaker_nodes_internal++;
     860           0 :             return VIRTUAL_SPEAKER_NODE_DISTRIBUTE_ENERGY;
     861             :         }
     862             : 
     863        2668 :         return NO_VIRTUAL_SPEAKER_NODE; /* No virtual back needed */
     864             :     }
     865             : 
     866             :     /* The following considers if TOP or BOTTOM virtual speaker is needed */
     867       48812 :     for ( ch = 0; ch < hVBAPdata->num_speaker_nodes; ch++ )
     868             :     {
     869       43476 :         if ( group == SPEAKER_NODE_TOP_HALF )
     870             :         {
     871       21738 :             if ( speaker_node_ele_deg[ch] > max_elevation )
     872             :             {
     873        1869 :                 max_elevation = speaker_node_ele_deg[ch];
     874             :             }
     875             :         }
     876             :         else
     877             :         {
     878       21738 :             if ( ( -speaker_node_ele_deg[ch] ) > max_elevation )
     879             :             {
     880           3 :                 max_elevation = -speaker_node_ele_deg[ch];
     881             :             }
     882             :         }
     883             :     }
     884             : 
     885        5336 :     if ( max_elevation > VBAP_NO_VIRTUAL_SPEAKER_NODE_ELE_LIMIT - VBAP_EPSILON )
     886             :     {
     887           0 :         return NO_VIRTUAL_SPEAKER_NODE;
     888             :     }
     889             : 
     890             :     /* Use virtual node */
     891        5336 :     if ( group == SPEAKER_NODE_BOTTOM_HALF )
     892             :     {
     893        2668 :         hVBAPdata->bottom_virtual_speaker_node_index = hVBAPdata->num_speaker_nodes_internal;
     894             :     }
     895             :     else
     896             :     {
     897        2668 :         hVBAPdata->top_virtual_speaker_node_index = hVBAPdata->num_speaker_nodes_internal;
     898             :     }
     899             : 
     900        5336 :     hVBAPdata->num_speaker_nodes_internal++;
     901        5336 :     if ( max_elevation > VBAP_DISTRIBUTE_VIRTUAL_SPEAKER_NODE_ELE_LIMIT - VBAP_EPSILON )
     902             :     {
     903        1872 :         return VIRTUAL_SPEAKER_NODE_DISTRIBUTE_ENERGY;
     904             :     }
     905             : 
     906        3464 :     return VIRTUAL_SPEAKER_NODE_DISCARD_ENERGY;
     907             : }
     908             : 
     909             : 
     910             : /*-------------------------------------------------------------------------*
     911             :  * init_speaker_node_direction_data()
     912             :  *
     913             :  * Initialize speaker node data
     914             :  *-------------------------------------------------------------------------*/
     915             : 
     916        2668 : static void init_speaker_node_direction_data(
     917             :     VBAP_SPEAKER_NODE *speaker_node_data, /* o  : storage for speaker node data     */
     918             :     const float *speaker_node_azi_deg,    /* i  : vector of speaker node azimuths   */
     919             :     const float *speaker_node_ele_deg,    /* i  : vector of speaker node elevations */
     920             :     const int16_t num_speaker_nodes       /* i  : number of speaker nodes           */
     921             : )
     922             : {
     923             :     int16_t ch;
     924             :     float azi_rad;
     925             :     float ele_rad;
     926        2668 :     int16_t num_horiz = 0;
     927        2668 :     uint8_t in_all_mode = TRUE;
     928             : 
     929       29742 :     for ( ch = 0; ch < num_speaker_nodes; ch++ )
     930             :     {
     931       27074 :         speaker_node_data[ch].azi_deg = speaker_node_azi_deg[ch];
     932       27074 :         azi_rad = speaker_node_azi_deg[ch] * PI_OVER_180;
     933       27074 :         if ( ( speaker_node_ele_deg[ch] >= -5 ) && ( speaker_node_ele_deg[ch] <= 5 ) )
     934             :         {
     935       15454 :             speaker_node_data[ch].ele_deg = 0.0f;
     936       15454 :             ele_rad = 0.0f;
     937       15454 :             speaker_node_data[ch].group = SPEAKER_NODE_HORIZONTAL;
     938       15454 :             num_horiz++;
     939             :         }
     940             :         else
     941             :         {
     942       11620 :             speaker_node_data[ch].ele_deg = speaker_node_ele_deg[ch];
     943       11620 :             ele_rad = speaker_node_ele_deg[ch] * PI_OVER_180;
     944       11620 :             if ( ele_rad < 0 )
     945             :             {
     946        2680 :                 speaker_node_data[ch].group = SPEAKER_NODE_BOTTOM_HALF;
     947             :             }
     948             :             else
     949             :             {
     950        8940 :                 speaker_node_data[ch].group = SPEAKER_NODE_TOP_HALF;
     951             :             }
     952             :         }
     953             : 
     954       27074 :         speaker_node_data[ch].unit_vec[0] = cosf( azi_rad ) * cosf( ele_rad );
     955       27074 :         speaker_node_data[ch].unit_vec[1] = sinf( azi_rad ) * cosf( ele_rad );
     956       27074 :         speaker_node_data[ch].unit_vec[2] = sinf( ele_rad );
     957             :     }
     958             : 
     959             :     /* Check for largest horizontal gap if there are at least 3 horizontal speaker nodes */
     960        2668 :     if ( num_horiz >= 3 )
     961             :     {
     962             :         int16_t i;
     963             :         uint16_t horiz_azi[VBAP_MAX_NUM_SPEAKER_NODES];
     964             :         uint16_t largest_gap;
     965             :         uint16_t temp;
     966             : 
     967        2668 :         i = 0;
     968       18122 :         for ( ch = 0; ch < num_speaker_nodes && i < num_horiz; ch++ )
     969             :         {
     970       15454 :             if ( speaker_node_data[ch].group == SPEAKER_NODE_HORIZONTAL )
     971             :             {
     972       15454 :                 horiz_azi[i] = speaker_node_azi_deg[ch] < 0.0f ? (uint16_t) floorf( speaker_node_azi_deg[ch] + 360.0f ) : (uint16_t) floorf( speaker_node_azi_deg[ch] );
     973       15454 :                 i++;
     974             :             }
     975             :         }
     976             : 
     977             :         /* Reorder horizontal azi to increasing order */
     978        2668 :         sort( horiz_azi, num_horiz );
     979             : 
     980             :         /* Find largest gap. Initialize with the wrap over gap. */
     981        2668 :         largest_gap = horiz_azi[0] - horiz_azi[num_horiz - 1] + 360;
     982       15454 :         for ( ch = 0; ch < num_horiz - 1; ch++ )
     983             :         {
     984       12786 :             temp = horiz_azi[ch + 1] - horiz_azi[ch];
     985       12786 :             if ( temp > largest_gap )
     986             :             {
     987        5273 :                 largest_gap = temp;
     988             :             }
     989             :         }
     990             : 
     991             :         /* If largest gap is small enough, we have definitive zero elevation plane.
     992             :          * Otherwise, we should assign all speaker nodes to one group. */
     993        2668 :         if ( largest_gap <= VBAP_MAX_HORIZONTAL_GAP )
     994             :         {
     995        2668 :             in_all_mode = FALSE;
     996             :         }
     997             :     }
     998             : 
     999             :     /* Designate all speaker nodes to same group if there was no definitive zero
    1000             :      * elevation plane. */
    1001        2668 :     if ( in_all_mode )
    1002             :     {
    1003           0 :         for ( ch = 0; ch < num_speaker_nodes; ch++ )
    1004             :         {
    1005           0 :             speaker_node_data[ch].group = SPEAKER_NODE_ALL;
    1006             :         }
    1007             :     }
    1008             : 
    1009        2668 :     return;
    1010             : }
    1011             : 
    1012             : /*-------------------------------------------------------------------------*
    1013             :  * matrix_inverse_3x3()
    1014             :  *
    1015             :  * 3-by-3 matrix inverse
    1016             :  *-------------------------------------------------------------------------*/
    1017             : 
    1018       43476 : static void matrix_inverse_3x3(
    1019             :     const float **input_matrix, /* i  : input matrix  */
    1020             :     float inverse_matrix[3][3]  /* o  : output matrix */
    1021             : )
    1022             : {
    1023             :     int16_t k;
    1024             :     float determinant;
    1025             :     float cross_vec[3];
    1026             : 
    1027       43476 :     vbap_crossp( input_matrix[1], input_matrix[2], cross_vec );
    1028             : 
    1029       43476 :     determinant = dotp( input_matrix[0], cross_vec, 3 );
    1030             : 
    1031      173904 :     for ( k = 0; k < 3; k++ )
    1032             :     {
    1033      130428 :         inverse_matrix[k][0] = cross_vec[k] / determinant;
    1034             :     }
    1035             : 
    1036       43476 :     vbap_crossp( input_matrix[2], input_matrix[0], cross_vec );
    1037             : 
    1038      173904 :     for ( k = 0; k < 3; k++ )
    1039             :     {
    1040      130428 :         inverse_matrix[k][1] = cross_vec[k] / determinant;
    1041             :     }
    1042             : 
    1043       43476 :     vbap_crossp( input_matrix[0], input_matrix[1], cross_vec );
    1044             : 
    1045      173904 :     for ( k = 0; k < 3; k++ )
    1046             :     {
    1047      130428 :         inverse_matrix[k][2] = cross_vec[k] / determinant;
    1048             :     }
    1049             : 
    1050       43476 :     return;
    1051             : }
    1052             : 
    1053             : /*-------------------------------------------------------------------------*
    1054             :  * check_and_store_triplet()
    1055             :  *
    1056             :  * Check if the given loudspeaker triplet is a valid one and store data when
    1057             :  * valid triplet is found.
    1058             :  *-------------------------------------------------------------------------*/
    1059             : 
    1060       43476 : static int16_t check_and_store_triplet(
    1061             :     const int16_t chA,                          /* i  : first channel index that forms the loudspeaker triplet  */
    1062             :     const int16_t chB,                          /* i  : second channel index that forms the loudspeaker triplet */
    1063             :     const int16_t chC,                          /* i  : third channel index that forms the loudspeaker triplet  */
    1064             :     const int16_t num_speaker_nodes,            /* i  : number of speaker nodes                                 */
    1065             :     const VBAP_SPEAKER_NODE *speaker_node_data, /* i  : speaker node data structure                             */
    1066             :     VBAP_VS_TRIPLET *triplets,                  /* o  : vector of virtual surface triplets                      */
    1067             :     int16_t *triplet_index,                     /* i/o: index for the next free triplet slot                    */
    1068             :     float *triplet_azidegs,                     /* o  : center azimuths of the found triplets                   */
    1069             :     int16_t *triplet_order                      /* o  : initial order of triplet indices                        */
    1070             : )
    1071             : {
    1072             :     int16_t ch_check;
    1073             :     int16_t k;
    1074             :     int16_t speaker_node_found_inside_triplet;
    1075             :     uint8_t triplet_ok;
    1076             :     float inverse_matrix[3][3], unnormalized_gains[3];
    1077             :     const float *speaker_node_triplet_unit_vec_matrix[3];
    1078             : 
    1079             :     /* Triplet found, determine inverse matrix for VBAP formulation */
    1080       43476 :     speaker_node_triplet_unit_vec_matrix[0] = speaker_node_data[chA].unit_vec;
    1081       43476 :     speaker_node_triplet_unit_vec_matrix[1] = speaker_node_data[chB].unit_vec;
    1082       43476 :     speaker_node_triplet_unit_vec_matrix[2] = speaker_node_data[chC].unit_vec;
    1083       43476 :     matrix_inverse_3x3( speaker_node_triplet_unit_vec_matrix, inverse_matrix );
    1084             : 
    1085             :     /* Check through all speaker nodes that none of them are within the triplet.
    1086             :      * Node within the triplet is identified by that all three panning gains are positive.
    1087             :      * Epsilon-condition is for some small rounding issues.*/
    1088       43476 :     speaker_node_found_inside_triplet = 0;
    1089      504088 :     for ( ch_check = 0; ch_check < num_speaker_nodes; ch_check++ )
    1090             :     {
    1091      460612 :         if ( ( ch_check != chA ) && ( ch_check != chB ) && ( ch_check != chC ) )
    1092             :         {
    1093      330184 :             triplet_ok = vector_matrix_multiply_3x3( speaker_node_data[ch_check].unit_vec, inverse_matrix, unnormalized_gains );
    1094      330184 :             if ( triplet_ok && unnormalized_gains[0] > VBAP_EPSILON && unnormalized_gains[1] > VBAP_EPSILON && unnormalized_gains[2] > VBAP_EPSILON )
    1095             :             {
    1096           0 :                 speaker_node_found_inside_triplet = 1;
    1097           0 :                 break;
    1098             :             }
    1099             :         }
    1100             :     }
    1101             : 
    1102             :     /* No speaker node inside triplet -> appropriate triplet found, save data. */
    1103       43476 :     if ( speaker_node_found_inside_triplet == 0 )
    1104             :     {
    1105       43476 :         triplets[*triplet_index].speaker_node[0] = (uint8_t) chA;
    1106       43476 :         triplets[*triplet_index].speaker_node[1] = (uint8_t) chB;
    1107       43476 :         triplets[*triplet_index].speaker_node[2] = (uint8_t) chC;
    1108      173904 :         for ( k = 0; k < 3; k++ )
    1109             :         {
    1110      130428 :             mvr2r( inverse_matrix[k], triplets[*triplet_index].inverse_matrix[k], 3 );
    1111             :         }
    1112             : 
    1113             :         /* Get center azimuth for fast search use */
    1114       43476 :         triplet_azidegs[*triplet_index] = atan2f( speaker_node_data[chA].unit_vec[1] + speaker_node_data[chB].unit_vec[1] + speaker_node_data[chC].unit_vec[1],
    1115       43476 :                                                   speaker_node_data[chA].unit_vec[0] + speaker_node_data[chB].unit_vec[0] + speaker_node_data[chC].unit_vec[0] ) *
    1116             :                                           _180_OVER_PI;
    1117             : 
    1118             :         /* Store increasing order indices for the later sorting step. */
    1119       43476 :         triplet_order[*triplet_index] = *triplet_index;
    1120             : 
    1121       43476 :         ( *triplet_index )++;
    1122             : 
    1123       43476 :         return 1;
    1124             :     }
    1125             : 
    1126             :     /* Triplet was not good */
    1127           0 :     return 0;
    1128             : }
    1129             : 
    1130             : 
    1131             : /*-------------------------------------------------------------------------*
    1132             :  * determine_virtual_surface_triplets()
    1133             :  *
    1134             :  * Determine virtual surface triples that are used for panning. This
    1135             :  * function is optimized for the use in cases where speaker nodes are in
    1136             :  * one group or divided into two separate groups divided by a horizontal
    1137             :  * layer.
    1138             :  *-------------------------------------------------------------------------*/
    1139             : 
    1140             : /*! r: number of virtual surface triplets */
    1141        5336 : static int16_t determine_virtual_surface_triplets(
    1142             :     const int16_t num_speaker_nodes,                         /* i  : number of speaker nodes                                                            */
    1143             :     const VBAP_SPEAKER_NODE *speaker_node_data,              /* i  : speaker node data structure                                                        */
    1144             :     int16_t connections[][2],                                /* i  : vector of all connections                                                          */
    1145             :     const int16_t max_num_connections,                       /* i  : max number of connections                                                          */
    1146             :     VBAP_VS_TRIPLET *triplets,                               /* o  : vector of virtual surface triplets                                                 */
    1147             :     int16_t initial_search_indices[VBAP_NUM_SEARCH_SECTORS], /* o  : initial search indices for this set of triplets corresponding to the search struct */
    1148             :     enum SpeakerNodeGroup allowed_group                      /* i  : group of allowed speaker nodes for forming the triplets in this call               */
    1149             : )
    1150             : {
    1151             :     int16_t chA, chB, chC, k, l, m;
    1152        5336 :     int16_t num_triplets = 0;
    1153             :     int16_t num_connected_to_chA;
    1154             :     int16_t connected_to_chA[VBAP_MAX_NUM_SPEAKER_NODES];
    1155             :     int16_t connection_uses_left[VBAP_MAX_NUM_SPEAKER_NODES];
    1156             :     float triplet_azidegs[VBAP_MAX_NUM_TRIPLETS];
    1157             :     int16_t triplet_order[VBAP_MAX_NUM_TRIPLETS];
    1158             : 
    1159             :     /* Each connection can be used exactly by two different virtual surface triplets. */
    1160        5336 :     set_s( connection_uses_left, 2, VBAP_MAX_NUM_SPEAKER_NODES );
    1161             : 
    1162       59484 :     for ( chA = 0; chA < num_speaker_nodes; chA++ )
    1163             :     {
    1164             :         /* Early skip if not in correct group. */
    1165       54148 :         if ( speaker_node_data[chA].group != allowed_group )
    1166             :         {
    1167       42528 :             continue;
    1168             :         }
    1169             : 
    1170             :         /* Get all connections connected to current chA that have not been used by
    1171             :          * two triplets yet. */
    1172       11620 :         num_connected_to_chA = 0;
    1173      232416 :         for ( k = 0; k < max_num_connections; k++ )
    1174             :         {
    1175      220796 :             if ( ( connections[k][0] == chA || connections[k][1] == chA ) && connection_uses_left[k] > 0 )
    1176             :             {
    1177       47822 :                 connected_to_chA[num_connected_to_chA] = k;
    1178       47822 :                 num_connected_to_chA++;
    1179             :             }
    1180             :         }
    1181             : 
    1182             :         /* Check that we have enough connections to use. We need at least two available connections to form a triplet.
    1183             :          * This can fail in later stages when all connections are already used. */
    1184       11620 :         if ( num_connected_to_chA < 2 )
    1185             :         {
    1186        1938 :             continue;
    1187             :         }
    1188             : 
    1189             :         /* Try to form triplets from each valid connection. */
    1190       57504 :         for ( k = 0; k < num_connected_to_chA; k++ )
    1191             :         {
    1192       47822 :             int16_t connect_index_k = connected_to_chA[k];
    1193       47822 :             chB = connections[connect_index_k][0] == chA ? connections[connect_index_k][1] : connections[connect_index_k][0];
    1194       83334 :             for ( l = k + 1; l < num_connected_to_chA; l++ )
    1195             :             {
    1196       73652 :                 int16_t connect_index_l = connected_to_chA[l];
    1197       73652 :                 chC = connections[connect_index_l][0] == chA ? connections[connect_index_l][1] : connections[connect_index_l][0];
    1198             : 
    1199             :                 /* With chA, chB, and chC selected, we still need to find connection between chB and chC and verify that the triplet is valid */
    1200      893808 :                 for ( m = 0; m < max_num_connections; m++ )
    1201             :                 {
    1202      863632 :                     if ( ( connections[m][0] == chB && connections[m][1] == chC ) || ( connections[m][1] == chB && connections[m][0] == chC ) )
    1203             :                     {
    1204       43476 :                         if ( check_and_store_triplet( chA, chB, chC, num_speaker_nodes, speaker_node_data, triplets, &num_triplets, triplet_azidegs, triplet_order ) )
    1205             :                         {
    1206       43476 :                             connection_uses_left[connect_index_k]--;
    1207       43476 :                             connection_uses_left[connect_index_l]--;
    1208       43476 :                             connection_uses_left[m]--;
    1209       43476 :                             break;
    1210             :                         }
    1211             :                     }
    1212             :                 }
    1213             : 
    1214             :                 /* Check if chA-chB connection has been used in two triplets already. If yes, then break out of loop
    1215             :                  * as this connection cannot be used for more triplets and we need to continue with another chA-chB
    1216             :                  * connection. */
    1217       73652 :                 if ( connection_uses_left[connect_index_k] < 1 )
    1218             :                 {
    1219       38140 :                     break;
    1220             :                 }
    1221             :             }
    1222             :         }
    1223             :     }
    1224             : 
    1225             :     /* All triplets should be stored now. Sort them for search use and then determine the initial search indices for
    1226             :      * each search sector for this search struct. */
    1227        5336 :     v_sort_ind( triplet_azidegs, triplet_order, num_triplets );
    1228        5336 :     reorder_triplets( triplets, triplet_order, num_triplets );
    1229        5336 :     determine_initial_search_indices( num_triplets, triplet_azidegs, initial_search_indices );
    1230             : 
    1231        5336 :     return num_triplets;
    1232             : }
    1233             : 
    1234             : /*-------------------------------------------------------------------------*
    1235             :  * determine_initial_search_indices()
    1236             :  *
    1237             :  * Determine initial search indices used for fast search of correct triangle
    1238             :  *-------------------------------------------------------------------------*/
    1239             : 
    1240        5336 : static void determine_initial_search_indices(
    1241             :     const int16_t num_triplets,                             /* i  : number of triplets                */
    1242             :     const float triplet_azidegs[VBAP_MAX_NUM_TRIPLETS],     /* i  : azimuths of triplets (in degrees) */
    1243             :     int16_t initial_search_indices[VBAP_NUM_SEARCH_SECTORS] /* o  : initial search indices            */
    1244             : )
    1245             : {
    1246             :     int16_t i, j;
    1247             :     float sector_reference_azideg;
    1248             :     float sector_border_start_azideg;
    1249             :     float sector_border_end_azideg;
    1250             :     int16_t best_index;
    1251             :     float min_azideg_diff;
    1252             :     float azideg_diff;
    1253             : 
    1254       26680 :     for ( i = 0; i < VBAP_NUM_SEARCH_SECTORS; i++ )
    1255             :     {
    1256       21344 :         sector_border_start_azideg = i * VBAP_SEARCH_SECTOR_SIZE;
    1257       21344 :         sector_border_end_azideg = ( i + 1 ) * VBAP_SEARCH_SECTOR_SIZE;
    1258       21344 :         sector_reference_azideg = ( sector_border_start_azideg + sector_border_end_azideg ) / 2.0f;
    1259       21344 :         best_index = 0;
    1260       21344 :         min_azideg_diff = 9999.9f;
    1261             : 
    1262      195248 :         for ( j = 0; j < num_triplets; j++ )
    1263             :         {
    1264      173904 :             azideg_diff = sector_reference_azideg - triplet_azidegs[j];
    1265      173904 :             if ( azideg_diff > 180.0f )
    1266             :             {
    1267       75624 :                 azideg_diff -= 360.0f;
    1268             :             }
    1269       98280 :             else if ( azideg_diff < -180.0f )
    1270             :             {
    1271           0 :                 azideg_diff += 360.0f;
    1272             :             }
    1273      173904 :             azideg_diff = fabsf( azideg_diff );
    1274             : 
    1275      173904 :             if ( azideg_diff < min_azideg_diff )
    1276             :             {
    1277       72982 :                 min_azideg_diff = azideg_diff;
    1278       72982 :                 best_index = j;
    1279             :             }
    1280             :         }
    1281             : 
    1282       21344 :         initial_search_indices[i] = best_index;
    1283             :     }
    1284             : 
    1285        5336 :     return;
    1286             : }
    1287             : 
    1288             : 
    1289             : /*-------------------------------------------------------------------------*
    1290             :  * determine_connections()
    1291             :  *
    1292             :  * Determine all valid connections between all speaker nodes
    1293             :  *-------------------------------------------------------------------------*/
    1294             : 
    1295        2668 : static ivas_error determine_connections(
    1296             :     const int16_t num_speaker_nodes,            /* i  : number of speaker nodes               */
    1297             :     const VBAP_SPEAKER_NODE *speaker_node_data, /* i  : speaker node data                     */
    1298             :     int16_t connections[][2],                   /* o  : vector of connections                 */
    1299             :     const int16_t max_num_connections,          /* i  : max number of connections             */
    1300             :     int16_t *group1_count,                      /* o  : number of connections in first group  */
    1301             :     int16_t *group2_start,                      /* o  : start of second group of connections  */
    1302             :     int16_t *group2_count                       /* o  : number of connections in second group */
    1303             : )
    1304             : {
    1305             :     int16_t num_non_crossing_planes;
    1306             :     int16_t c;
    1307        2668 :     int16_t connection_write_index = 0;
    1308             :     float non_crossing_plane_elevation_deg[VBAP_MAX_PLANES];
    1309             :     ivas_error error;
    1310             : 
    1311        2668 :     set_f( non_crossing_plane_elevation_deg, 0.0f, VBAP_MAX_PLANES );
    1312             : 
    1313             :     /* Reset connection data */
    1314       67882 :     for ( c = 0; c < max_num_connections; c++ )
    1315             :     {
    1316       65214 :         connections[c][0] = VBAP_NOT_VALID_CONNECTION;
    1317             :     }
    1318             : 
    1319             :     /* This function determines some prominent elevated planes, that are favoured in making node-node connections. */
    1320        2668 :     num_non_crossing_planes = determine_non_crossing_planes( num_speaker_nodes, speaker_node_data, non_crossing_plane_elevation_deg );
    1321             : 
    1322             :     /* Process in different mode based on the grouping. It is enough to check for first node. */
    1323        2668 :     if ( speaker_node_data[0].group == SPEAKER_NODE_ALL )
    1324             :     {
    1325           0 :         if ( ( error = formulate_half_sphere_connections( speaker_node_data, num_speaker_nodes, SPEAKER_NODE_ALL, connections, &connection_write_index, max_num_connections, num_non_crossing_planes, non_crossing_plane_elevation_deg ) ) != IVAS_ERR_OK )
    1326             :         {
    1327           0 :             return error;
    1328             :         }
    1329             :     }
    1330             :     else
    1331             :     {
    1332             :         /* The node-node connections are determined in three stages: bottom, horizontal, and top. */
    1333        2668 :         if ( ( error = formulate_half_sphere_connections( speaker_node_data, num_speaker_nodes, SPEAKER_NODE_BOTTOM_HALF, connections, &connection_write_index, max_num_connections, num_non_crossing_planes, non_crossing_plane_elevation_deg ) ) != IVAS_ERR_OK )
    1334             :         {
    1335           0 :             return error;
    1336             :         }
    1337        2668 :         *group2_start = connection_write_index;
    1338             : 
    1339        2668 :         formulate_horizontal_connections( speaker_node_data, num_speaker_nodes, connections, &connection_write_index );
    1340        2668 :         *group1_count = connection_write_index;
    1341             : 
    1342        2668 :         if ( ( error = formulate_half_sphere_connections( speaker_node_data, num_speaker_nodes, SPEAKER_NODE_TOP_HALF, connections, &connection_write_index, max_num_connections, num_non_crossing_planes, non_crossing_plane_elevation_deg ) ) != IVAS_ERR_OK )
    1343             :         {
    1344           0 :             return error;
    1345             :         }
    1346        2668 :         *group2_count = connection_write_index - *group2_start;
    1347             :     }
    1348             : 
    1349        2668 :     return IVAS_ERR_OK;
    1350             : }
    1351             : 
    1352             : 
    1353             : /*-------------------------------------------------------------------------*
    1354             :  * determine_connection_class()
    1355             :  *
    1356             :  * Determine the type of connection
    1357             :  *-------------------------------------------------------------------------*/
    1358             : 
    1359             : /*! r: type of connection */
    1360       81386 : static enum ConnectionClass determine_connection_class(
    1361             :     const VBAP_SPEAKER_NODE *node_data, /* i  : speaker node data       */
    1362             :     const int16_t num_speaker_nodes,    /* i  : number of speaker nodes */
    1363             :     const enum SpeakerNodeGroup group,  /* i  : speaker node group      */
    1364             :     const int16_t chA,                  /* i  : speaker node counter 1  */
    1365             :     const int16_t chB                   /* i  : speaker node counter 2  */
    1366             : )
    1367             : {
    1368             :     int16_t ch, k;
    1369             :     const float *p1, *v2;
    1370             :     float v1v1, v1v2, v2v2, v1p1, v2p1;
    1371             :     float determinant;
    1372             :     float norm_distance_on_v1;
    1373             :     float vec_diff[3];
    1374             :     float v1[3];
    1375             :     float vTarget[3];
    1376             :     float energy_sum;
    1377             :     float eq_value;
    1378             :     float uvecdot;
    1379             : 
    1380             :     /* Check if connection passes through origin. This is not desired.
    1381             :      * When this happens, unit vectors point in opposite directions. */
    1382       81386 :     uvecdot = dotp( node_data[chA].unit_vec, node_data[chB].unit_vec, 3 ) + 1.0f;
    1383       81386 :     if ( uvecdot < VBAP_EPSILON && uvecdot > -VBAP_EPSILON )
    1384             :     {
    1385           0 :         return CONNECTION_WITH_SPEAKER_NODE_BEHIND;
    1386             :     }
    1387             : 
    1388             :     /* This loop checks the chA-chB connection with respect to all loudspeakers:
    1389             :      - in case there is a node behind or nearly behind the connection line. These
    1390             :      connections need to be discarded.
    1391             :      - in case there is a node that is closer to the connection line than 1/5 of the
    1392             :      connection length AND at the same horizontal plane. These connections need to be
    1393             :      weighted with a penalty (a special case, for example avoiding elevated L,R,C triplet)
    1394             :      */
    1395      945410 :     for ( ch = 0; ch < num_speaker_nodes; ch++ )
    1396             :     {
    1397             :         /* Select speaker_node only within TOP or BOTTOM sphere half, not being part of chA-chB pair */
    1398      872668 :         if ( ( group == node_data[ch].group ) && ( ch != chA ) && ( ch != chB ) )
    1399             :         {
    1400             :             /* The following lines formulate the point on the chA-chB-connection that is
    1401             :              nearest to the origo-ch-line */
    1402      194772 :             p1 = node_data[chA].unit_vec;
    1403      779088 :             for ( k = 0; k < 3; k++ )
    1404             :             {
    1405      584316 :                 v1[k] = node_data[chB].unit_vec[k] - node_data[chA].unit_vec[k];
    1406             :             }
    1407      194772 :             v2 = node_data[ch].unit_vec;
    1408      194772 :             v1v1 = dotp( v1, v1, 3 );
    1409      194772 :             v1v2 = dotp( v1, v2, 3 );
    1410      194772 :             v2v2 = 1.0f; /* dotp(v2, v2) is always 1. */
    1411      194772 :             v1p1 = dotp( v1, p1, 3 );
    1412      194772 :             v2p1 = dotp( v2, p1, 3 );
    1413      194772 :             determinant = ( v1v1 ) * ( -v2v2 ) + ( v1v2 * v1v2 );
    1414             :             /* Norm distance = distance parameter on line chA-chB, determines point that is
    1415             :              nearest to origo-ch line. Distance 0 means chA and distance 1 means chB. Can be
    1416             :              outside this region as well.*/
    1417      194772 :             norm_distance_on_v1 = ( 1.0f / determinant ) * ( ( v2v2 ) * ( v1p1 ) + ( -v1v2 ) * ( v2p1 ) );
    1418             : 
    1419             :             /* Continue only if the nearest point is between chA and chB */
    1420      194772 :             if ( norm_distance_on_v1 > 0.0f && norm_distance_on_v1 < 1.0f )
    1421             :             {
    1422             :                 /* Formulate vTarget, that is an unit vector that goes through the determined point on chA-chB connection */
    1423      159602 :                 energy_sum = 0.0f;
    1424      638408 :                 for ( k = 0; k < 3; k++ )
    1425             :                 {
    1426      478806 :                     vTarget[k] = p1[k] + norm_distance_on_v1 * v1[k];
    1427      478806 :                     energy_sum += vTarget[k] * vTarget[k];
    1428      478806 :                     vec_diff[k] = vTarget[k] - v2[k];
    1429             :                 }
    1430      159602 :                 eq_value = sqrtf( 1.0f / energy_sum );
    1431      638408 :                 for ( k = 0; k < 3; k++ )
    1432             :                 {
    1433      478806 :                     vTarget[k] *= eq_value;
    1434             :                 }
    1435             : 
    1436             :                 /* A check if the angle between vTarget and node_data[ch].unit_vec is less than 1 degree.
    1437             :                  Essentially reveals if there is a speaker node too close "behind" the connection. Such
    1438             :                  connections should be rejected.*/
    1439      159602 :                 if ( dotp( vTarget, v2, 3 ) > 0.9998f )
    1440             :                 {
    1441        8644 :                     return CONNECTION_WITH_SPEAKER_NODE_BEHIND;
    1442             :                 }
    1443             : 
    1444             :                 /* A special case, mainly accounting for ELEVATED L,R,C speaker nodes.
    1445             :                  A triplet between these nodes is not desired if there is a top node,
    1446             :                  a penalty is implemented to take care of this. */
    1447      150958 :                 if ( sqrtf( v1v1 ) > 5.0f * sqrtf( dotp( vec_diff, vec_diff, 3 ) ) )
    1448             :                 {
    1449           0 :                     if ( fabsf( node_data[chB].unit_vec[2] - node_data[chA].unit_vec[2] ) < VBAP_EPSILON )
    1450             :                     {
    1451           0 :                         return ELEVATED_PLANE_THIN_TRIANGLE_CONNECTION;
    1452             :                     }
    1453             :                 }
    1454             :             }
    1455             :         }
    1456             :     }
    1457             : 
    1458       72742 :     return REGULAR_CONNECTION;
    1459             : }
    1460             : 
    1461             : 
    1462             : /*-------------------------------------------------------------------------*
    1463             :  * formulate_horizontal_connections()
    1464             :  *
    1465             :  * Formulate connections in the horizontal plane
    1466             :  *-------------------------------------------------------------------------*/
    1467             : 
    1468        2668 : static void formulate_horizontal_connections(
    1469             :     const VBAP_SPEAKER_NODE *speaker_node_data, /* i  : speaker node data         */
    1470             :     const int16_t num_speaker_nodes,            /* i  : number of speaker nodes   */
    1471             :     int16_t connections[][2],                   /* o  : vector of all connections */
    1472             :     int16_t *connection_write_index )
    1473             : {
    1474             :     int16_t ch;
    1475             :     int16_t chCheck;
    1476             :     int16_t next_index;
    1477             :     float min_arc_diff;
    1478             :     float arc_diff;
    1479             : 
    1480       29742 :     for ( ch = 0; ch < num_speaker_nodes; ch++ )
    1481             :     {
    1482             :         /* Find next horizontal speaker node */
    1483       27074 :         if ( speaker_node_data[ch].group == SPEAKER_NODE_HORIZONTAL )
    1484             :         {
    1485       15454 :             next_index = -1;
    1486       15454 :             min_arc_diff = 9999.0f;
    1487      174396 :             for ( chCheck = 0; chCheck < num_speaker_nodes; chCheck++ )
    1488             :             {
    1489      158942 :                 if ( ( ch != chCheck ) && ( speaker_node_data[chCheck].group == SPEAKER_NODE_HORIZONTAL ) )
    1490             :                 {
    1491       76608 :                     arc_diff = speaker_node_data[chCheck].azi_deg - speaker_node_data[ch].azi_deg;
    1492      114912 :                     while ( arc_diff < 0.0f )
    1493             :                     {
    1494       38304 :                         arc_diff += 360.0f;
    1495             :                     }
    1496       76608 :                     if ( arc_diff < min_arc_diff )
    1497             :                     {
    1498       34687 :                         min_arc_diff = arc_diff;
    1499       34687 :                         next_index = chCheck;
    1500             :                     }
    1501             :                 }
    1502             :             }
    1503       15454 :             connections[*connection_write_index][0] = ch;
    1504       15454 :             connections[*connection_write_index][1] = next_index;
    1505       15454 :             ( *connection_write_index )++;
    1506             :         }
    1507             :     }
    1508             : 
    1509        2668 :     return;
    1510             : }
    1511             : 
    1512             : /*-------------------------------------------------------------------------*
    1513             :  * check_plane_crossing()
    1514             :  *
    1515             :  * Check crossing of non-allowed planes
    1516             :  *-------------------------------------------------------------------------*/
    1517             : 
    1518             : /*! r: truth value for crossing */
    1519       72742 : static int16_t check_plane_crossing(
    1520             :     const float ele1_deg,                         /* i  : speaker node 1 elevation            */
    1521             :     const float ele2_deg,                         /* i  : speaker node 2 elevation            */
    1522             :     const int16_t num_non_crossing_planes,        /* i  : number of non-crossing planes       */
    1523             :     const float *non_crossing_plane_elevation_deg /* i  : vector non-crossing plane elevations*/
    1524             : )
    1525             : {
    1526             :     /* Find if the connection crosses a non-crossing plane, with 1-degree threshold. */
    1527             :     int16_t plane;
    1528             : 
    1529      120726 :     for ( plane = 0; plane < num_non_crossing_planes; plane++ )
    1530             :     {
    1531       50416 :         if ( ( ele1_deg > non_crossing_plane_elevation_deg[plane] + 1.0f ) && ( ele2_deg < non_crossing_plane_elevation_deg[plane] - 1.0f ) )
    1532             :         {
    1533          24 :             return 1;
    1534             :         }
    1535             : 
    1536       50392 :         if ( ( ele2_deg > non_crossing_plane_elevation_deg[plane] + 1.0f ) && ( ele1_deg < non_crossing_plane_elevation_deg[plane] - 1.0f ) )
    1537             :         {
    1538        2408 :             return 1;
    1539             :         }
    1540             :     }
    1541             : 
    1542       70310 :     return 0;
    1543             : }
    1544             : 
    1545             : 
    1546             : /*-------------------------------------------------------------------------*
    1547             :  * get_half_sphere_connection_options()
    1548             :  *
    1549             :  * Get list of all potential connections at the half-sphere
    1550             :  *-------------------------------------------------------------------------*/
    1551             : 
    1552        5336 : static ivas_error get_half_sphere_connection_options(
    1553             :     const VBAP_SPEAKER_NODE *speaker_node_data,    /* i  : speaker node data       */
    1554             :     const enum SpeakerNodeGroup group,             /* i  : speaker node group      */
    1555             :     const int16_t num_speaker_nodes,               /* i  : number of speaker nodes */
    1556             :     const int16_t num_non_crossing_planes,         /* i  : number of non-crossing planes           */
    1557             :     const float *non_crossing_plane_elevation_deg, /* i  : vector of non-crossing plane elevations */
    1558             :     ConnectionOption **connection_options_pr,      /* o  : list of connection options */
    1559             :     int16_t *num_connection_options                /* o  : number of connection options */
    1560             : )
    1561             : {
    1562        5336 :     int16_t max_num_connection_options = 0;
    1563        5336 :     int16_t index = 0;
    1564             :     int16_t node, chA, chB, c, c_cmp;
    1565             :     ConnectionOption *c_options, *c_options_reorder;
    1566             : 
    1567             :     /* Count max num connection options at the half sphere */
    1568       59484 :     for ( node = 0; node < num_speaker_nodes; node++ )
    1569             :     {
    1570       54148 :         if ( speaker_node_data[node].group == group || speaker_node_data[node].group == SPEAKER_NODE_HORIZONTAL )
    1571             :         {
    1572       42528 :             max_num_connection_options += index;
    1573       42528 :             index++;
    1574             :         }
    1575             :     }
    1576             : 
    1577             :     /* Init memory for connection options */
    1578        5336 :     if ( ( c_options = (ConnectionOption *) malloc( sizeof( ConnectionOption ) * max_num_connection_options ) ) == NULL )
    1579             :     {
    1580           0 :         return ( IVAS_ERROR( IVAS_ERR_FAILED_ALLOC, "Can not allocate memory for VBAP data\n" ) );
    1581             :     }
    1582             : 
    1583      163330 :     for ( c = 0; c < max_num_connection_options; c++ )
    1584             :     {
    1585      157994 :         c_options[c].chA = -1;
    1586      157994 :         c_options[c].chB = -1;
    1587      157994 :         c_options[c].arc = 1e20f;
    1588      157994 :         c_options[c].arc_weighted = 1e20f;
    1589             :     }
    1590             : 
    1591             :     /* Determine connection options for the half-sphere */
    1592        5336 :     index = 0;
    1593       59484 :     for ( chA = 0; chA < num_speaker_nodes; chA++ )
    1594             :     {
    1595             :         /* First loudspeaker at the connection is at the half sphere */
    1596       54148 :         if ( speaker_node_data[chA].group == group || speaker_node_data[chA].group == SPEAKER_NODE_HORIZONTAL )
    1597             :         {
    1598      276402 :             for ( chB = chA + 1; chB < num_speaker_nodes; chB++ )
    1599             :             {
    1600             :                 /* Second loudspeaker at the connection is at the half sphere, but so that first and second are not both horizontal. */
    1601      233874 :                 if ( ( speaker_node_data[chB].group == group ) || ( speaker_node_data[chB].group == SPEAKER_NODE_HORIZONTAL && speaker_node_data[chA].group == group ) )
    1602             :                 {
    1603       81386 :                     int16_t ConnectionClass = determine_connection_class( speaker_node_data, num_speaker_nodes, group, chA, chB );
    1604             :                     /* Connection is considered only if there is no speaker node behind it */
    1605       81386 :                     if ( ConnectionClass != CONNECTION_WITH_SPEAKER_NODE_BEHIND )
    1606             :                     {
    1607             :                         /* Store connection information */
    1608       72742 :                         c_options[index].chA = chA;
    1609       72742 :                         c_options[index].chB = chB;
    1610       72742 :                         c_options[index].arc = acosf( dotp( speaker_node_data[chA].unit_vec, speaker_node_data[chB].unit_vec, 3 ) );
    1611       72742 :                         c_options[index].arc_weighted = c_options[index].arc;
    1612             : 
    1613             :                         /* A special case, mainly accounting for ELEVATED L,R,C speaker nodes.
    1614             :                         A triplet between these nodes is not desired if there is a top node,
    1615             :                         a penalty is implemented to take care of this. */
    1616       72742 :                         if ( ConnectionClass == ELEVATED_PLANE_THIN_TRIANGLE_CONNECTION )
    1617             :                         {
    1618           0 :                             c_options[index].arc_weighted *= 2.0f;
    1619             :                         }
    1620             : 
    1621             :                         /* If the connection passes a pre-determined plane of speaker nodes, then add further penalty */
    1622       72742 :                         if ( check_plane_crossing( speaker_node_data[chA].ele_deg, speaker_node_data[chB].ele_deg, num_non_crossing_planes, non_crossing_plane_elevation_deg ) )
    1623             :                         {
    1624        2432 :                             c_options[index].arc_weighted *= 2.0f;
    1625             :                         }
    1626       72742 :                         index++;
    1627             :                     }
    1628             :                 }
    1629             :             }
    1630             :         }
    1631             :     }
    1632             :     /* Number of found connection options at the half sphere */
    1633        5336 :     *num_connection_options = index;
    1634             : 
    1635             :     /* Init memory for reordered connection options and order by arc_weighted,
    1636             :      * which informs of the preference order of the connections in case they cross */
    1637        5336 :     if ( ( c_options_reorder = (ConnectionOption *) malloc( sizeof( ConnectionOption ) * ( *num_connection_options ) ) ) == NULL )
    1638             :     {
    1639           0 :         return ( IVAS_ERROR( IVAS_ERR_FAILED_ALLOC, "Can not allocate memory for VBAP data\n" ) );
    1640             :     }
    1641             : 
    1642       78078 :     for ( c = 0; c < *num_connection_options; c++ )
    1643             :     {
    1644             :         float min_arc_weighted;
    1645             :         int16_t min_arc_index;
    1646       72742 :         min_arc_weighted = 1e20f;
    1647       72742 :         min_arc_index = -1;
    1648             : 
    1649     1799280 :         for ( c_cmp = 0; c_cmp < *num_connection_options; c_cmp++ )
    1650             :         {
    1651     1726538 :             if ( c_options[c_cmp].arc_weighted <= min_arc_weighted )
    1652             :             {
    1653      378360 :                 min_arc_weighted = c_options[c_cmp].arc_weighted;
    1654      378360 :                 min_arc_index = c_cmp;
    1655             :             }
    1656             :         }
    1657       72742 :         c_options_reorder[c].chA = c_options[min_arc_index].chA;
    1658       72742 :         c_options_reorder[c].chB = c_options[min_arc_index].chB;
    1659       72742 :         c_options_reorder[c].arc = c_options[min_arc_index].arc;
    1660       72742 :         c_options_reorder[c].arc_weighted = c_options[min_arc_index].arc_weighted;
    1661       72742 :         c_options[min_arc_index].arc_weighted = 1e20f;
    1662             :     }
    1663             : 
    1664             :     /* Set reordered connections as output and free temporary data */
    1665        5336 :     *connection_options_pr = c_options_reorder;
    1666        5336 :     free( c_options );
    1667             : 
    1668        5336 :     return IVAS_ERR_OK;
    1669             : }
    1670             : 
    1671             : 
    1672             : /*-------------------------------------------------------------------------*
    1673             :  * formulate_half_sphere_connections()
    1674             :  *
    1675             :  * Formulate half-sphere connections
    1676             :  *-------------------------------------------------------------------------*/
    1677             : 
    1678        5336 : static ivas_error formulate_half_sphere_connections(
    1679             :     const VBAP_SPEAKER_NODE *speaker_node_data, /* i  : speaker node data                       */
    1680             :     const int16_t num_speaker_nodes,            /* i  : number of speaker nodes                 */
    1681             :     const enum SpeakerNodeGroup group,          /* i  : speaker node group                      */
    1682             :     int16_t connections[][2],                   /* o  : vector of connections                   */
    1683             :     int16_t *connection_write_index,
    1684             :     const int16_t max_num_connections,            /* i  : max number of connections               */
    1685             :     const int16_t num_non_crossing_planes,        /* i  : number of non-crossing planes           */
    1686             :     const float *non_crossing_plane_elevation_deg /* i  : vector of non-crossing plane elevations */
    1687             : )
    1688             : {
    1689             :     /* Variable initializations */
    1690             :     int16_t c, chA, chB, cmp_chA, cmp_chB, k, c_opt;
    1691             :     int16_t new_connection_is_valid;
    1692             :     int16_t within_first_arc;
    1693             :     float planeCrossingVec[3];
    1694             :     float new_arc;
    1695             :     float new_cross[3];
    1696             :     float tmpFloat;
    1697             :     float cmp_arc;
    1698             :     float normVal;
    1699             :     float angleCmp;
    1700             :     float connection_arc[( VBAP_MAX_NUM_SPEAKER_NODES - 2 ) * 3];
    1701             :     float connection_cross[( VBAP_MAX_NUM_SPEAKER_NODES - 2 ) * 3][3];
    1702             :     ConnectionOption *connection_options;
    1703             :     int16_t num_connection_options;
    1704             :     int16_t half_sphere_first_connection;
    1705             :     ivas_error error;
    1706             : 
    1707        5336 :     half_sphere_first_connection = *connection_write_index;
    1708             : 
    1709             :     /* Obtain all connection options (i.e., channel pairs) at the half sphere. The function orders them
    1710             :      * in terms of which connection to keep if two connections would cross each other. */
    1711        5336 :     if ( ( error = get_half_sphere_connection_options(
    1712             :                speaker_node_data,
    1713             :                group,
    1714             :                num_speaker_nodes,
    1715             :                num_non_crossing_planes,
    1716             :                non_crossing_plane_elevation_deg,
    1717             :                &connection_options,
    1718             :                &num_connection_options ) ) != IVAS_ERR_OK )
    1719             :     {
    1720           0 :         return error;
    1721             :     }
    1722             : 
    1723        5336 :     set_f( connection_arc, 0.0f, max_num_connections );
    1724      135764 :     for ( c = 0; c < max_num_connections; c++ )
    1725             :     {
    1726      130428 :         set_f( connection_cross[c], 0.0f, 3 );
    1727             :     }
    1728             : 
    1729             :     /* The following loop goes through all reasonable chA - chB pairs for the half-sphere */
    1730        5336 :     c_opt = 0;
    1731       62854 :     while ( c_opt < num_connection_options && *connection_write_index < max_num_connections )
    1732             :     {
    1733       57518 :         chA = connection_options[c_opt].chA;
    1734       57518 :         chB = connection_options[c_opt].chB;
    1735       57518 :         new_arc = connection_options[c_opt].arc;
    1736             : 
    1737             :         /* Cross-product is needed for later stages */
    1738       57518 :         vbap_crossp( speaker_node_data[chA].unit_vec, speaker_node_data[chB].unit_vec, new_cross );
    1739             : 
    1740             :         /* Determine if new connection between chA and chB is valid */
    1741       57518 :         new_connection_is_valid = 1;
    1742       57518 :         c = half_sphere_first_connection;
    1743      392337 :         while ( ( c < *connection_write_index ) && new_connection_is_valid )
    1744             :         {
    1745      334819 :             cmp_chA = connections[c][0];
    1746      334819 :             cmp_chB = connections[c][1];
    1747             :             /* The connections are compared only if they don't involve same speaker nodes */
    1748      334819 :             if ( ( cmp_chA != chA ) && ( cmp_chA != chB ) && ( cmp_chB != chA ) && ( cmp_chB != chB ) )
    1749             :             {
    1750             :                 /* The following lines determine if the connection chA-chB crosses with the connection cmp_chA-cmp_chB.*/
    1751             :                 /* The connections, i.e., node-pairs determine a plane. The crossing can be determined by
    1752             :                  * studying the intersection of these planes. */
    1753      165132 :                 vbap_crossp( connection_cross[c], new_cross, planeCrossingVec );
    1754      165132 :                 tmpFloat = 1e-12f;
    1755      165132 :                 cmp_arc = connection_arc[c];
    1756      660528 :                 for ( k = 0; k < 3; k++ )
    1757             :                 {
    1758      495396 :                     tmpFloat += planeCrossingVec[k] * planeCrossingVec[k];
    1759             :                 }
    1760      165132 :                 normVal = sqrtf( 1.0f / tmpFloat );
    1761             : 
    1762      660528 :                 for ( k = 0; k < 3; k++ )
    1763             :                 {
    1764      495396 :                     planeCrossingVec[k] *= normVal;
    1765             :                 }
    1766             : 
    1767             :                 /* If the plane intersection is between both connections, then the two connections cross. */
    1768             :                 /* Study first if the crossing is between arc chA-chB */
    1769      165132 :                 angleCmp = acosf( dotp( planeCrossingVec, speaker_node_data[chA].unit_vec, 3 ) );
    1770      165132 :                 angleCmp += acosf( dotp( planeCrossingVec, speaker_node_data[chB].unit_vec, 3 ) );
    1771             : 
    1772      165132 :                 within_first_arc = 0;
    1773      165132 :                 if ( fabsf( new_arc - angleCmp ) < 0.01f )
    1774             :                 {
    1775       35914 :                     within_first_arc = 1;
    1776             :                 }
    1777      129218 :                 else if ( fabsf( new_arc - ( 2.0f * EVS_PI - angleCmp ) ) < 0.01f )
    1778             :                 {
    1779       37687 :                     within_first_arc = 1;
    1780             :                     /* In this case, the plane crossing vector is inverted. The inverse is another
    1781             :                      * plane-crossing vector, and detected to be between chA-chB connection.*/
    1782      150748 :                     for ( k = 0; k < 3; k++ )
    1783             :                     {
    1784      113061 :                         planeCrossingVec[k] *= -1.0f;
    1785             :                     }
    1786             :                 }
    1787             : 
    1788             :                 /* Study if the crossing is also between arc cmp_chA-cmp_chB */
    1789      165132 :                 if ( within_first_arc > 0 )
    1790             :                 {
    1791       73601 :                     angleCmp = acosf( dotp( planeCrossingVec, speaker_node_data[cmp_chA].unit_vec, 3 ) );
    1792       73601 :                     angleCmp += acosf( dotp( planeCrossingVec, speaker_node_data[cmp_chB].unit_vec, 3 ) );
    1793       73601 :                     if ( fabsf( cmp_arc - angleCmp ) < 0.01f )
    1794             :                     {
    1795             :                         /* A crossing is detected. The new connection is not valid, because
    1796             :                          * the connections were ordered in order of preference (arc_weighted) */
    1797        7758 :                         new_connection_is_valid = 0;
    1798             :                     }
    1799             :                 }
    1800             :             }
    1801      334819 :             c++;
    1802             :         }
    1803             : 
    1804             :         /* Store the new connection which has been confirmed valid */
    1805       57518 :         if ( new_connection_is_valid > 0 )
    1806             :         {
    1807       49760 :             connections[*connection_write_index][0] = chA;
    1808       49760 :             connections[*connection_write_index][1] = chB;
    1809       49760 :             connection_arc[*connection_write_index] = new_arc;
    1810       49760 :             connection_cross[*connection_write_index][0] = new_cross[0];
    1811       49760 :             connection_cross[*connection_write_index][1] = new_cross[1];
    1812       49760 :             connection_cross[*connection_write_index][2] = new_cross[2];
    1813       49760 :             ( *connection_write_index )++;
    1814             :         }
    1815       57518 :         c_opt++;
    1816             :     }
    1817             : 
    1818        5336 :     free( connection_options );
    1819             : 
    1820        5336 :     return IVAS_ERR_OK;
    1821             : }
    1822             : 
    1823             : 
    1824             : /*-------------------------------------------------------------------------*
    1825             :  * determine_non_crossing_planes()
    1826             :  *
    1827             :  * Determine non-crossing planes
    1828             :  *-------------------------------------------------------------------------*/
    1829             : 
    1830             : /*! r: number of non-crossing planes */
    1831        2668 : static int16_t determine_non_crossing_planes(
    1832             :     const int16_t num_speaker_nodes,        /* i  : number of speaker nodes                 */
    1833             :     const VBAP_SPEAKER_NODE *node_data,     /* i  : speaker node data                       */
    1834             :     float *non_crossing_plane_elevation_deg /* o  : vector of non-crossing plane elevations */
    1835             : )
    1836             : {
    1837             :     float next_ele_check;
    1838             :     float ele_check;
    1839             :     float max_gap;
    1840             :     float gap_to_next_ls;
    1841             :     int16_t ch, ch_cmp;
    1842             :     int16_t num_planes;
    1843             : 
    1844        2668 :     ele_check = -999.9f;
    1845        2668 :     num_planes = 0;
    1846             : 
    1847             :     /* For each plane, check if a non-crossing plane should be determined */
    1848        9876 :     while ( ele_check < 90.0f )
    1849             :     {
    1850        7208 :         next_ele_check = 999.9f;
    1851             :         /* Find next node elevation that is not in horizontal plane */
    1852       81906 :         for ( ch = 0; ch < num_speaker_nodes; ch++ )
    1853             :         {
    1854       74698 :             if ( ( node_data[ch].group != SPEAKER_NODE_HORIZONTAL ) && ( node_data[ch].ele_deg > ele_check + VBAP_EPSILON ) && ( node_data[ch].ele_deg < next_ele_check - VBAP_EPSILON ) )
    1855             :             {
    1856        9083 :                 next_ele_check = node_data[ch].ele_deg;
    1857             :             }
    1858             :         }
    1859        7208 :         ele_check = next_ele_check;
    1860        7208 :         if ( ele_check > 90.0f )
    1861             :         {
    1862             :             /* When no next node elevation found, break loop */
    1863           0 :             break;
    1864             :         }
    1865             : 
    1866        7208 :         max_gap = -9999.9f;
    1867       81906 :         for ( ch = 0; ch < num_speaker_nodes; ch++ )
    1868             :         {
    1869             :             /* Find gap to the next speaker node at the same plane */
    1870       74698 :             if ( fabsf( node_data[ch].ele_deg - ele_check ) < VBAP_EPSILON )
    1871             :             {
    1872       11620 :                 gap_to_next_ls = 99999.9f;
    1873      137132 :                 for ( ch_cmp = 0; ch_cmp < num_speaker_nodes; ch_cmp++ )
    1874             :                 {
    1875      125512 :                     if ( ch_cmp != ch && fabsf( node_data[ch_cmp].ele_deg - ele_check ) < VBAP_EPSILON )
    1876             :                     {
    1877       16444 :                         float gap = node_data[ch_cmp].azi_deg - node_data[ch].azi_deg;
    1878       24666 :                         while ( gap < 0 )
    1879             :                         {
    1880        8222 :                             gap += 360.0f;
    1881             :                         }
    1882       16444 :                         if ( gap < gap_to_next_ls )
    1883             :                         {
    1884       11364 :                             gap_to_next_ls = gap;
    1885             :                         }
    1886             :                     }
    1887             :                 }
    1888             :                 /* Find maximum gap on that plane */
    1889       11620 :                 if ( gap_to_next_ls > max_gap )
    1890             :                 {
    1891        7888 :                     max_gap = gap_to_next_ls;
    1892             :                 }
    1893             :             }
    1894             :         }
    1895             : 
    1896             :         /* If maximum gap is small enough, then a non-crossing plane is detected */
    1897        7208 :         if ( max_gap < ( VBAP_MAX_HORIZONTAL_GAP_FOR_PLANE_DETECTION + VBAP_EPSILON ) && max_gap > 0.0f )
    1898             :         {
    1899        1270 :             non_crossing_plane_elevation_deg[num_planes] = ele_check;
    1900        1270 :             num_planes++;
    1901        1270 :             if ( num_planes == VBAP_MAX_PLANES )
    1902             :             {
    1903             :                 /* Memory init limit. Does not happen with any real speaker node configuration.
    1904             :                  Triangulation succeeds even if number of non_crossing_planes are limited. */
    1905           0 :                 break;
    1906             :             }
    1907             :         }
    1908             :     }
    1909             : 
    1910        2668 :     return num_planes;
    1911             : }
    1912             : 
    1913             : 
    1914             : /*-------------------------------------------------------------------------*
    1915             :  * reorder_triplets()
    1916             :  *
    1917             :  * Reorder virtual surface triplets into provided target order.
    1918             :  *-------------------------------------------------------------------------*/
    1919             : 
    1920        5336 : static void reorder_triplets(
    1921             :     VBAP_VS_TRIPLET *triplets,   /* i/o: VS triplets to be reordered  */
    1922             :     const int16_t *target_order, /* i  : Target order for VS triplets */
    1923             :     const int16_t num_triplets   /* i  : Number of VS triplets        */
    1924             : )
    1925             : {
    1926             :     VBAP_VS_TRIPLET tempTriplets[VBAP_MAX_NUM_TRIPLETS];
    1927             :     int16_t c;
    1928             : 
    1929             :     /* First copy to temp array */
    1930       48812 :     for ( c = 0; c < num_triplets; c++ )
    1931             :     {
    1932       43476 :         tempTriplets[c] = triplets[c];
    1933             :     }
    1934             : 
    1935             :     /* Then move back in sorted order */
    1936       48812 :     for ( c = 0; c < num_triplets; c++ )
    1937             :     {
    1938       43476 :         triplets[c] = tempTriplets[target_order[c]];
    1939             :     }
    1940             : 
    1941        5336 :     return;
    1942             : }

Generated by: LCOV version 1.14