celt/bands.c

Bug Summary

File:	libs/opus-1.1-p2/celt/bands.c
Location:	line 595, column 12
Description:	Assigned value is garbage or undefined

Annotated Source Code

Written by Jean-Marc Valin and Gregory Maxwell */

Redistribution and use in source and binary forms, with or without

modification, are permitted provided that the following conditions

are met:

- Redistributions of source code must retain the above copyright

notice, this list of conditions and the following disclaimer.

- Redistributions in binary form must reproduce the above copyright

notice, this list of conditions and the following disclaimer in the

documentation and/or other materials provided with the distribution.

THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS

``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT

LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR

A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER

OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,

EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,

PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR

PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF

LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING

NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS

SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.

#ifdef HAVE_CONFIG_H1

#include "config.h"

#endif

#include <math.h>

#include "bands.h"

#include "modes.h"

#include "vq.h"

#include "cwrs.h"

#include "stack_alloc.h"

#include "os_support.h"

#include "mathops.h"

#include "rate.h"

#include "quant_bands.h"

#include "pitch.h"

int hysteresis_decision(opus_val16 val, const opus_val16 *thresholds, const opus_val16 *hysteresis, int N, int prev)

{

int i;

for (i=0;i<N;i++)

{

if (val < thresholds[i])

break;

}

if (i>prev && val < thresholds[prev]+hysteresis[prev])

i=prev;

if (i<prev && val > thresholds[prev-1]-hysteresis[prev-1])

i=prev;

return i;

}

opus_uint32 celt_lcg_rand(opus_uint32 seed)

{

return 1664525 * seed + 1013904223;

}

/* This is a cos() approximation designed to be bit-exact on any platform. Bit exactness

with this approximation is important because it has an impact on the bit allocation */

static opus_int16 bitexact_cos(opus_int16 x)

{

opus_int32 tmp;

opus_int16 x2;

tmp = (4096+((opus_int32)(x)*(x)))>>13;

celt_assert(tmp<=32767);

x2 = tmp;

x2 = (32767-x2) + FRAC_MUL16(x2, (-7651 + FRAC_MUL16(x2, (8277 + FRAC_MUL16(-626, x2)))))((16384+((opus_int32)(opus_int16)(x2)*(opus_int16)((-7651 + (
(16384+((opus_int32)(opus_int16)(x2)*(opus_int16)((8277 + ((16384
+((opus_int32)(opus_int16)(-626)*(opus_int16)(x2)))>>15
)))))>>15)))))>>15);

celt_assert(x2<=32766);

return 1+x2;

}

static int bitexact_log2tan(int isin,int icos)

{

int lc;

int ls;

lc=EC_ILOG(icos)(((int)sizeof(unsigned)*8)-(__builtin_clz(icos)));

ls=EC_ILOG(isin)(((int)sizeof(unsigned)*8)-(__builtin_clz(isin)));

icos<<=15-lc;

isin<<=15-ls;

return (ls-lc)*(1<<11)

+FRAC_MUL16(isin, FRAC_MUL16(isin, -2597) + 7932)((16384+((opus_int32)(opus_int16)(isin)*(opus_int16)(((16384+
((opus_int32)(opus_int16)(isin)*(opus_int16)(-2597)))>>
15) + 7932)))>>15)

-FRAC_MUL16(icos, FRAC_MUL16(icos, -2597) + 7932)((16384+((opus_int32)(opus_int16)(icos)*(opus_int16)(((16384+
((opus_int32)(opus_int16)(icos)*(opus_int16)(-2597)))>>
15) + 7932)))>>15);

}

#ifdef FIXED_POINT

/* Compute the amplitude (sqrt energy) in each of the bands */

void compute_band_energies(const CELTModeOpusCustomMode *m, const celt_sig *X, celt_ener *bandE, int end, int C, int M)

{

int i, c, N;

const opus_int16 *eBands = m->eBands;

N = M*m->shortMdctSize;

100

c=0; do {

101

for (i=0;i<end;i++)

102

{

103

int j;

104

opus_val32 maxval=0;

105

opus_val32 sum = 0;

106

107

j=M*eBands[i]; do {

108

maxval = MAX32(maxval, X[j+c*N])((maxval) > (X[j+c*N]) ? (maxval) : (X[j+c*N]));

109

maxval = MAX32(maxval, -X[j+c*N])((maxval) > (-X[j+c*N]) ? (maxval) : (-X[j+c*N]));

110

} while (++j<M*eBands[i+1]);

111

112

if (maxval > 0)

113

{

114

int shift = celt_ilog2(maxval)-10;

115

j=M*eBands[i]; do {

116

sum = MAC16_16(sum, EXTRACT16(VSHR32(X[j+c*N],shift)),((sum)+(opus_val32)(((X[j+c*N])))*(opus_val32)(((X[j+c*N]))))

117

EXTRACT16(VSHR32(X[j+c*N],shift)))((sum)+(opus_val32)(((X[j+c*N])))*(opus_val32)(((X[j+c*N]))));

118

} while (++j<M*eBands[i+1]);

119

/* We're adding one here to ensure the normalized band isn't larger than unity norm */

120

bandE[i+c*m->nbEBands] = EPSILON1e-15f+VSHR32(EXTEND32(celt_sqrt(sum)),-shift)((((float)sqrt(sum))));

121

} else {

122

bandE[i+c*m->nbEBands] = EPSILON1e-15f;

123

}

124

/*printf ("%f ", bandE[i+c*m->nbEBands]);*/

125

}

126

} while (++c<C);

127

/*printf ("\n");*/

128

}

129

130

/* Normalise each band such that the energy is one. */

131

void normalise_bands(const CELTModeOpusCustomMode *m, const celt_sig * OPUS_RESTRICT__restrict freq, celt_norm * OPUS_RESTRICT__restrict X, const celt_ener *bandE, int end, int C, int M)

132

{

133

int i, c, N;

134

const opus_int16 *eBands = m->eBands;

135

N = M*m->shortMdctSize;

136

c=0; do {

137

i=0; do {

138

opus_val16 g;

139

int j,shift;

140

opus_val16 E;

141

shift = celt_zlog2(bandE[i+c*m->nbEBands])-13;

142

E = VSHR32(bandE[i+c*m->nbEBands], shift)(bandE[i+c*m->nbEBands]);

143

g = EXTRACT16(celt_rcp(SHL32(E,3)))((1.f/((E))));

144

j=M*eBands[i]; do {

145

X[j+c*N] = MULT16_16_Q15(VSHR32(freq[j+c*N],shift-1),g)(((freq[j+c*N]))*(g));

146

} while (++j<M*eBands[i+1]);

147

} while (++i<end);

148

} while (++c<C);

149

}

150

151

#else /* FIXED_POINT */

152

/* Compute the amplitude (sqrt energy) in each of the bands */

153

void compute_band_energies(const CELTModeOpusCustomMode *m, const celt_sig *X, celt_ener *bandE, int end, int C, int M)

154

{

155

int i, c, N;

156

const opus_int16 *eBands = m->eBands;

157

N = M*m->shortMdctSize;

158

c=0; do {

159

for (i=0;i<end;i++)

160

{

161

int j;

162

opus_val32 sum = 1e-27f;

163

for (j=M*eBands[i];j<M*eBands[i+1];j++)

164

sum += X[j+c*N]*X[j+c*N];

165

bandE[i+c*m->nbEBands] = celt_sqrt(sum)((float)sqrt(sum));

166

/*printf ("%f ", bandE[i+c*m->nbEBands]);*/

167

}

168

} while (++c<C);

169

/*printf ("\n");*/

170

}

171

172

/* Normalise each band such that the energy is one. */

173

void normalise_bands(const CELTModeOpusCustomMode *m, const celt_sig * OPUS_RESTRICT__restrict freq, celt_norm * OPUS_RESTRICT__restrict X, const celt_ener *bandE, int end, int C, int M)

174

{

175

int i, c, N;

176

const opus_int16 *eBands = m->eBands;

177

N = M*m->shortMdctSize;

178

c=0; do {

179

for (i=0;i<end;i++)

180

{

181

int j;

182

opus_val16 g = 1.f/(1e-27f+bandE[i+c*m->nbEBands]);

183

for (j=M*eBands[i];j<M*eBands[i+1];j++)

184

X[j+c*N] = freq[j+c*N]*g;

185

}

186

} while (++c<C);

187

}

188

189

#endif /* FIXED_POINT */

190

191

/* De-normalise the energy to produce the synthesis from the unit-energy bands */

192

void denormalise_bands(const CELTModeOpusCustomMode *m, const celt_norm * OPUS_RESTRICT__restrict X,

193

celt_sig * OPUS_RESTRICT__restrict freq, const opus_val16 *bandLogE, int start, int end, int C, int M)

194

{

195

int i, c, N;

196

const opus_int16 *eBands = m->eBands;

197

N = M*m->shortMdctSize;

198

celt_assert2(C<=2, "denormalise_bands() not implemented for >2 channels");

199

c=0; do {

200

celt_sig * OPUS_RESTRICT__restrict f;

201

const celt_norm * OPUS_RESTRICT__restrict x;

202

f = freq+c*N;

203

x = X+c*N+M*eBands[start];

204

for (i=0;i<M*eBands[start];i++)

205

*f++ = 0;

206

for (i=start;i<end;i++)

207

{

208

int j, band_end;

209

opus_val16 g;

210

opus_val16 lg;

211

#ifdef FIXED_POINT

212

int shift;

213

#endif

214

j=M*eBands[i];

215

band_end = M*eBands[i+1];

216

lg = ADD16(bandLogE[i+c*m->nbEBands], SHL16((opus_val16)eMeans[i],6))((bandLogE[i+c*m->nbEBands])+(((opus_val16)eMeans[i])));

217

#ifndef FIXED_POINT

218

g = celt_exp2(lg)((float)exp(0.6931471805599453094*(lg)));

219

#else

220

/* Handle the integer part of the log energy */

221

shift = 16-(lg>>DB_SHIFT);

222

if (shift>31)

223

{

224

shift=0;

225

g=0;

226

} else {

227

/* Handle the fractional part. */

228

g = celt_exp2_frac(lg&((1<<DB_SHIFT)-1));

229

}

230

/* Handle extreme gains with negative shift. */

231

if (shift<0)

232

{

233

/* For shift < -2 we'd be likely to overflow, so we're capping

234

the gain here. This shouldn't happen unless the bitstream is

235

already corrupted. */

236

if (shift < -2)

237

{

238

g = 32767;

239

shift = -2;

240

}

241

do {

242

*f++ = SHL32(MULT16_16(*x++, g), -shift)(((opus_val32)(*x++)*(opus_val32)(g)));

243

} while (++j<band_end);

244

} else

245

#endif

246

/* Be careful of the fixed-point "else" just above when changing this code */

247

do {

248

*f++ = SHR32(MULT16_16(*x++, g), shift)(((opus_val32)(*x++)*(opus_val32)(g)));

249

} while (++j<band_end);

250

}

251

celt_assert(start <= end);

252

for (i=M*eBands[end];i<N;i++)

253

*f++ = 0;

254

} while (++c<C);

255

}

256

257

/* This prevents energy collapse for transients with multiple short MDCTs */

258

void anti_collapse(const CELTModeOpusCustomMode *m, celt_norm *X_, unsigned char *collapse_masks, int LM, int C, int size,

259

int start, int end, opus_val16 *logE, opus_val16 *prev1logE,

260

opus_val16 *prev2logE, int *pulses, opus_uint32 seed)

261

{

262

int c, i, j, k;

263

for (i=start;i<end;i++)

264

{

265

int N0;

266

opus_val16 thresh, sqrt_1;

267

int depth;

268

#ifdef FIXED_POINT

269

int shift;

270

opus_val32 thresh32;

271

#endif

272

273

N0 = m->eBands[i+1]-m->eBands[i];

274

/* depth in 1/8 bits */

275

depth = (1+pulses[i])/((m->eBands[i+1]-m->eBands[i])<<LM);

276

277

#ifdef FIXED_POINT

278

thresh32 = SHR32(celt_exp2(-SHL16(depth, 10-BITRES)),1)(((float)exp(0.6931471805599453094*(-(depth)))));

279

thresh = MULT16_32_Q15(QCONST16(0.5f, 15), MIN32(32767,thresh32))(((0.5f))*(((32767) < (thresh32) ? (32767) : (thresh32))));

280

{

281

opus_val32 t;

282

t = N0<<LM;

283

shift = celt_ilog2(t)>>1;

284

t = SHL32(t, (7-shift)<<1)(t);

285

sqrt_1 = celt_rsqrt_norm(t)((1.f/((float)sqrt(t))));

286

}

287

#else

288

thresh = .5f*celt_exp2(-.125f*depth)((float)exp(0.6931471805599453094*(-.125f*depth)));

289

sqrt_1 = celt_rsqrt(N0<<LM)(1.f/((float)sqrt(N0<<LM)));

290

#endif

291

292

c=0; do

293

{

294

celt_norm *X;

295

opus_val16 prev1;

296

opus_val16 prev2;

297

opus_val32 Ediff;

298

opus_val16 r;

299

int renormalize=0;

300

prev1 = prev1logE[c*m->nbEBands+i];

301

prev2 = prev2logE[c*m->nbEBands+i];

302

if (C==1)

303

{

304

prev1 = MAX16(prev1,prev1logE[m->nbEBands+i])((prev1) > (prev1logE[m->nbEBands+i]) ? (prev1) : (prev1logE
[m->nbEBands+i]));

305

prev2 = MAX16(prev2,prev2logE[m->nbEBands+i])((prev2) > (prev2logE[m->nbEBands+i]) ? (prev2) : (prev2logE
[m->nbEBands+i]));

306

}

307

Ediff = EXTEND32(logE[c*m->nbEBands+i])(logE[c*m->nbEBands+i])-EXTEND32(MIN16(prev1,prev2))(((prev1) < (prev2) ? (prev1) : (prev2)));

308

Ediff = MAX32(0, Ediff)((0) > (Ediff) ? (0) : (Ediff));

309

310

#ifdef FIXED_POINT

311

if (Ediff < 16384)

312

{

313

opus_val32 r32 = SHR32(celt_exp2(-EXTRACT16(Ediff)),1)(((float)exp(0.6931471805599453094*(-(Ediff)))));

314

r = 2*MIN16(16383,r32)((16383) < (r32) ? (16383) : (r32));

315

} else {

316

r = 0;

317

}

318

if (LM==3)

319

r = MULT16_16_Q14(23170, MIN32(23169, r))((23170)*(((23169) < (r) ? (23169) : (r))));

320

r = SHR16(MIN16(thresh, r),1)(((thresh) < (r) ? (thresh) : (r)));

321

r = SHR32(MULT16_16_Q15(sqrt_1, r),shift)(((sqrt_1)*(r)));

322

#else

323

/* r needs to be multiplied by 2 or 2*sqrt(2) depending on LM because

324

short blocks don't have the same energy as long */

325

r = 2.f*celt_exp2(-Ediff)((float)exp(0.6931471805599453094*(-Ediff)));

326

if (LM==3)

327

r *= 1.41421356f;

328

r = MIN16(thresh, r)((thresh) < (r) ? (thresh) : (r));

329

r = r*sqrt_1;

330

#endif

331

X = X_+c*size+(m->eBands[i]<<LM);

332

for (k=0;k<1<<LM;k++)

333

{

334

/* Detect collapse */

335

if (!(collapse_masks[i*C+c]&1<<k))

336

{

337

/* Fill with noise */

338

for (j=0;j<N0;j++)

339

{

340

seed = celt_lcg_rand(seed);

341

X[(j<<LM)+k] = (seed&0x8000 ? r : -r);

342

}

343

renormalize = 1;

344

}

345

}

346

/* We just added some energy, so we need to renormalise */

347

if (renormalize)

348

renormalise_vector(X, N0<<LM, Q15ONE1.0f);

349

} while (++c<C);

350

}

351

}

352

353

static void intensity_stereo(const CELTModeOpusCustomMode *m, celt_norm *X, celt_norm *Y, const celt_ener *bandE, int bandID, int N)

354

{

355

int i = bandID;

356

int j;

357

opus_val16 a1, a2;

358

opus_val16 left, right;

359

opus_val16 norm;

360

#ifdef FIXED_POINT

361

int shift = celt_zlog2(MAX32(bandE[i], bandE[i+m->nbEBands])((bandE[i]) > (bandE[i+m->nbEBands]) ? (bandE[i]) : (bandE
[i+m->nbEBands])))-13;

362

#endif

363

left = VSHR32(bandE[i],shift)(bandE[i]);

364

right = VSHR32(bandE[i+m->nbEBands],shift)(bandE[i+m->nbEBands]);

365

norm = EPSILON1e-15f + celt_sqrt(EPSILON+MULT16_16(left,left)+MULT16_16(right,right))((float)sqrt(1e-15f +((opus_val32)(left)*(opus_val32)(left))+
((opus_val32)(right)*(opus_val32)(right))));

366

a1 = DIV32_16(SHL32(EXTEND32(left),14),norm)(((opus_val32)(((left))))/(opus_val16)(norm));

367

a2 = DIV32_16(SHL32(EXTEND32(right),14),norm)(((opus_val32)(((right))))/(opus_val16)(norm));

368

for (j=0;j<N;j++)

369

{

370

celt_norm r, l;

371

l = X[j];

372

r = Y[j];

373

X[j] = MULT16_16_Q14(a1,l)((a1)*(l)) + MULT16_16_Q14(a2,r)((a2)*(r));

374

/* Side is not encoded, no need to calculate */

375

}

376

}

377

378

static void stereo_split(celt_norm *X, celt_norm *Y, int N)

379

{

380

int j;

381

for (j=0;j<N;j++)

382

{

383

celt_norm r, l;

384

l = MULT16_16_Q15(QCONST16(.70710678f,15), X[j])(((.70710678f))*(X[j]));

385

r = MULT16_16_Q15(QCONST16(.70710678f,15), Y[j])(((.70710678f))*(Y[j]));

386

X[j] = l+r;

387

Y[j] = r-l;

388

}

389

}

390

391

static void stereo_merge(celt_norm *X, celt_norm *Y, opus_val16 mid, int N)

392

{

393

int j;

394

opus_val32 xp=0, side=0;

395

opus_val32 El, Er;

396

opus_val16 mid2;

397

#ifdef FIXED_POINT

398

int kl, kr;

399

#endif

400

opus_val32 t, lgain, rgain;

401

402

/* Compute the norm of X+Y and X-Y as |X|^2 + |Y|^2 +/- sum(xy) */

403

dual_inner_prod(Y, X, Y, N, &xp, &side);

404

/* Compensating for the mid normalization */

405

xp = MULT16_32_Q15(mid, xp)((mid)*(xp));

406

/* mid and side are in Q15, not Q14 like X and Y */

407

mid2 = SHR32(mid, 1)(mid);

408

El = MULT16_16(mid2, mid2)((opus_val32)(mid2)*(opus_val32)(mid2)) + side - 2*xp;

409

Er = MULT16_16(mid2, mid2)((opus_val32)(mid2)*(opus_val32)(mid2)) + side + 2*xp;

410

if (Er < QCONST32(6e-4f, 28)(6e-4f) || El < QCONST32(6e-4f, 28)(6e-4f))

411

{

412

for (j=0;j<N;j++)

413

Y[j] = X[j];

414

return;

415

}

416

417

#ifdef FIXED_POINT

418

kl = celt_ilog2(El)>>1;

419

kr = celt_ilog2(Er)>>1;

420

#endif

421

t = VSHR32(El, (kl-7)<<1)(El);

422

lgain = celt_rsqrt_norm(t)((1.f/((float)sqrt(t))));

423

t = VSHR32(Er, (kr-7)<<1)(Er);

424

rgain = celt_rsqrt_norm(t)((1.f/((float)sqrt(t))));

425

426

#ifdef FIXED_POINT

427

if (kl < 7)

428

kl = 7;

429

if (kr < 7)

430

kr = 7;

431

#endif

432

433

for (j=0;j<N;j++)

434

{

435

celt_norm r, l;

436

/* Apply mid scaling (side is already scaled) */

437

l = MULT16_16_Q15(mid, X[j])((mid)*(X[j]));

438

r = Y[j];

439

X[j] = EXTRACT16(PSHR32(MULT16_16(lgain, SUB16(l,r)), kl+1))((((opus_val32)(lgain)*(opus_val32)(((l)-(r))))));

440

Y[j] = EXTRACT16(PSHR32(MULT16_16(rgain, ADD16(l,r)), kr+1))((((opus_val32)(rgain)*(opus_val32)(((l)+(r))))));

441

}

442

}

443

444

/* Decide whether we should spread the pulses in the current frame */

445

int spreading_decision(const CELTModeOpusCustomMode *m, celt_norm *X, int *average,

446

int last_decision, int *hf_average, int *tapset_decision, int update_hf,

447

int end, int C, int M)

448

{

449

int i, c, N0;

450

int sum = 0, nbBands=0;

451

const opus_int16 * OPUS_RESTRICT__restrict eBands = m->eBands;

452

int decision;

453

int hf_sum=0;

454

455

celt_assert(end>0);

456

457

N0 = M*m->shortMdctSize;

458

459

if (M*(eBands[end]-eBands[end-1]) <= 8)

460

return SPREAD_NONE(0);

461

c=0; do {

462

for (i=0;i<end;i++)

463

{

464

int j, N, tmp=0;

465

int tcount[3] = {0,0,0};

466

celt_norm * OPUS_RESTRICT__restrict x = X+M*eBands[i]+c*N0;

467

N = M*(eBands[i+1]-eBands[i]);

468

if (N<=8)

469

continue;

470

/* Compute rough CDF of |x[j]| */

471

for (j=0;j<N;j++)

472

{

473

opus_val32 x2N; /* Q13 */

474

475

x2N = MULT16_16(MULT16_16_Q15(x[j], x[j]), N)((opus_val32)(((x[j])*(x[j])))*(opus_val32)(N));

476

if (x2N < QCONST16(0.25f,13)(0.25f))

477

tcount[0]++;

478

if (x2N < QCONST16(0.0625f,13)(0.0625f))

479

tcount[1]++;

480

if (x2N < QCONST16(0.015625f,13)(0.015625f))

481

tcount[2]++;

482

}

483

484

/* Only include four last bands (8 kHz and up) */

485

if (i>m->nbEBands-4)

486

hf_sum += 32*(tcount[1]+tcount[0])/N;

487

tmp = (2*tcount[2] >= N) + (2*tcount[1] >= N) + (2*tcount[0] >= N);

488

sum += tmp*256;

489

nbBands++;

490

}

491

} while (++c<C);

492

493

if (update_hf)

494

{

495

if (hf_sum)

496

hf_sum /= C*(4-m->nbEBands+end);

497

*hf_average = (*hf_average+hf_sum)>>1;

498

hf_sum = *hf_average;

499

if (*tapset_decision==2)

500

hf_sum += 4;

501

else if (*tapset_decision==0)

502

hf_sum -= 4;

503

if (hf_sum > 22)

504

*tapset_decision=2;

505

else if (hf_sum > 18)

506

*tapset_decision=1;

507

else

508

*tapset_decision=0;

509

}

510

/*printf("%d %d %d\n", hf_sum, *hf_average, *tapset_decision);*/

511

celt_assert(nbBands>0); /* end has to be non-zero */

512

sum /= nbBands;

513

/* Recursive averaging */

514

sum = (sum+*average)>>1;

515

*average = sum;

516

/* Hysteresis */

517

sum = (3*sum + (((3-last_decision)<<7) + 64) + 2)>>2;

518

if (sum < 80)

519

{

520

decision = SPREAD_AGGRESSIVE(3);

521

} else if (sum < 256)

522

{

523

decision = SPREAD_NORMAL(2);

524

} else if (sum < 384)

525

{

526

decision = SPREAD_LIGHT(1);

527

} else {

528

decision = SPREAD_NONE(0);

529

}

530

#ifdef FUZZING

531

decision = rand()&0x3;

532

*tapset_decision=rand()%3;

533

#endif

534

return decision;

535

}

536

537

/* Indexing table for converting from natural Hadamard to ordery Hadamard

538

This is essentially a bit-reversed Gray, on top of which we've added

539

an inversion of the order because we want the DC at the end rather than

540

the beginning. The lines are for N=2, 4, 8, 16 */

541

static const int ordery_table[] = {

542

1, 0,

543

3, 0, 2, 1,

544

7, 0, 4, 3, 6, 1, 5, 2,

545

15, 0, 8, 7, 12, 3, 11, 4, 14, 1, 9, 6, 13, 2, 10, 5,

546

};

547

548

static void deinterleave_hadamard(celt_norm *X, int N0, int stride, int hadamard)

549

{

550

int i,j;

551

VARDECL(celt_norm, tmp);

552

int N;

553

SAVE_STACK;

554

N = N0*stride;

555

ALLOC(tmp, N, celt_norm)celt_norm tmp[N];

556

celt_assert(stride>0);

557

if (hadamard)

558

{

559

const int *ordery = ordery_table+stride-2;

560

for (i=0;i<stride;i++)

561

{

562

for (j=0;j<N0;j++)

563

tmp[ordery[i]*N0+j] = X[j*stride+i];

564

}

565

} else {

566

for (i=0;i<stride;i++)

567

for (j=0;j<N0;j++)

568

tmp[i*N0+j] = X[j*stride+i];

569

}

570

for (j=0;j<N;j++)

571

X[j] = tmp[j];

572

RESTORE_STACK;

573

}

574

575

static void interleave_hadamard(celt_norm *X, int N0, int stride, int hadamard)

576

{

577

int i,j;

578

VARDECL(celt_norm, tmp);

579

int N;

580

SAVE_STACK;

581

N = N0*stride;

582

ALLOC(tmp, N, celt_norm)celt_norm tmp[N];

583

if (hadamard)

←

Taking false branch

→

584

{

585

const int *ordery = ordery_table+stride-2;

586

for (i=0;i<stride;i++)

587

for (j=0;j<N0;j++)

588

tmp[j*stride+i] = X[ordery[i]*N0+j];

589

} else {

590

for (i=0;i<stride;i++)

←

Loop condition is true. Entering loop body

→

←

Loop condition is true. Entering loop body

→

←

Assuming 'i' is >= 'stride'

→

←

Loop condition is false. Execution continues on line 594

→

591

for (j=0;j<N0;j++)

←

Assuming 'j' is >= 'N0'

→

←

Loop condition is false. Execution continues on line 590

→

←

Loop condition is false. Execution continues on line 590

→

592

tmp[j*stride+i] = X[i*N0+j];

593

}

594

for (j=0;j<N;j++)

←

Assuming 'j' is < 'N'

→

←

Loop condition is true. Entering loop body

→

595

X[j] = tmp[j];

←

Assigned value is garbage or undefined

596

RESTORE_STACK;

597

}

598

599

void haar1(celt_norm *X, int N0, int stride)

600

{

601

int i, j;

602

N0 >>= 1;

603

for (i=0;i<stride;i++)

604

for (j=0;j<N0;j++)

605

{

606

celt_norm tmp1, tmp2;

607

tmp1 = MULT16_16_Q15(QCONST16(.70710678f,15), X[stride*2*j+i])(((.70710678f))*(X[stride*2*j+i]));

608

tmp2 = MULT16_16_Q15(QCONST16(.70710678f,15), X[stride*(2*j+1)+i])(((.70710678f))*(X[stride*(2*j+1)+i]));

609

X[stride*2*j+i] = tmp1 + tmp2;

610

X[stride*(2*j+1)+i] = tmp1 - tmp2;

611

}

612

}

613

614

static int compute_qn(int N, int b, int offset, int pulse_cap, int stereo)

615

{

616

static const opus_int16 exp2_table8[8] =

617

{16384, 17866, 19483, 21247, 23170, 25267, 27554, 30048};

618

int qn, qb;

619

int N2 = 2*N-1;

620

if (stereo && N==2)

621

N2--;

622

/* The upper limit ensures that in a stereo split with itheta==16384, we'll

623

always have enough bits left over to code at least one pulse in the

624

side; otherwise it would collapse, since it doesn't get folded. */

625

qb = IMIN(b-pulse_cap-(4<<BITRES), (b+N2*offset)/N2)((b-pulse_cap-(4<<3)) < ((b+N2*offset)/N2) ? (b-pulse_cap
-(4<<3)) : ((b+N2*offset)/N2));

626

627

qb = IMIN(8<<BITRES, qb)((8<<3) < (qb) ? (8<<3) : (qb));

628

629

if (qb<(1<<BITRES3>>1)) {

630

qn = 1;

631

} else {

632

qn = exp2_table8[qb&0x7]>>(14-(qb>>BITRES3));

633

qn = (qn+1)>>1<<1;

634

}

635

celt_assert(qn <= 256);

636

return qn;

637

}

638

639

struct band_ctx {

640

int encode;

641

const CELTModeOpusCustomMode *m;

642

int i;

643

int intensity;

644

int spread;

645

int tf_change;

646

ec_ctx *ec;

647

opus_int32 remaining_bits;

648

const celt_ener *bandE;

649

opus_uint32 seed;

650

};

651

652

struct split_ctx {

653

int inv;

654

int imid;

655

int iside;

656

int delta;

657

int itheta;

658

int qalloc;

659

};

660

661

static void compute_theta(struct band_ctx *ctx, struct split_ctx *sctx,

662

celt_norm *X, celt_norm *Y, int N, int *b, int B, int B0,

663

int LM,

664

int stereo, int *fill)

665

{

666

int qn;

667

int itheta=0;

668

int delta;

669

int imid, iside;

670

int qalloc;

671

int pulse_cap;

672

int offset;

673

opus_int32 tell;

674

int inv=0;

675

int encode;

676

const CELTModeOpusCustomMode *m;

677

int i;

678

int intensity;

679

ec_ctx *ec;

680

const celt_ener *bandE;

681

682

encode = ctx->encode;

683

m = ctx->m;

684

i = ctx->i;

685

intensity = ctx->intensity;

686

ec = ctx->ec;

687

bandE = ctx->bandE;

688

689

/* Decide on the resolution to give to the split parameter theta */

690

pulse_cap = m->logN[i]+LM*(1<<BITRES3);

691

offset = (pulse_cap>>1) - (stereo&&N==2 ? QTHETA_OFFSET_TWOPHASE16 : QTHETA_OFFSET4);

692

qn = compute_qn(N, *b, offset, pulse_cap, stereo);

693

if (stereo && i>=intensity)

694

qn = 1;

695

if (encode)

696

{

697

/* theta is the atan() of the ratio between the (normalized)

698

side and mid. With just that parameter, we can re-scale both

699

mid and side because we know that 1) they have unit norm and

700

2) they are orthogonal. */

701

itheta = stereo_itheta(X, Y, stereo, N);

702

}

703

tell = ec_tell_frac(ec);

704

if (qn!=1)

705

{

706

if (encode)

707

itheta = (itheta*qn+8192)>>14;

708

709

/* Entropy coding of the angle. We use a uniform pdf for the

710

time split, a step for stereo, and a triangular one for the rest. */

711

if (stereo && N>2)

712

{

713

int p0 = 3;

714

int x = itheta;

715

int x0 = qn/2;

716

int ft = p0*(x0+1) + x0;

717

/* Use a probability of p0 up to itheta=8192 and then use 1 after */

718

if (encode)

719

{

720

ec_encode(ec,x<=x0?p0*x:(x-1-x0)+(x0+1)*p0,x<=x0?p0*(x+1):(x-x0)+(x0+1)*p0,ft);

721

} else {

722

int fs;

723

fs=ec_decode(ec,ft);

724

if (fs<(x0+1)*p0)

725

x=fs/p0;

726

else

727

x=x0+1+(fs-(x0+1)*p0);

728

ec_dec_update(ec,x<=x0?p0*x:(x-1-x0)+(x0+1)*p0,x<=x0?p0*(x+1):(x-x0)+(x0+1)*p0,ft);

729

itheta = x;

730

}

731

} else if (B0>1 || stereo) {

732

/* Uniform pdf */

733

if (encode)

734

ec_enc_uint(ec, itheta, qn+1);

735

else

736

itheta = ec_dec_uint(ec, qn+1);

737

} else {

738

int fs=1, ft;

739

ft = ((qn>>1)+1)*((qn>>1)+1);

740

if (encode)

741

{

742

int fl;

743

744

fs = itheta <= (qn>>1) ? itheta + 1 : qn + 1 - itheta;

745

fl = itheta <= (qn>>1) ? itheta*(itheta + 1)>>1 :

746

ft - ((qn + 1 - itheta)*(qn + 2 - itheta)>>1);

747

748

ec_encode(ec, fl, fl+fs, ft);

749

} else {

750

/* Triangular pdf */

751

int fl=0;

752

int fm;

753

fm = ec_decode(ec, ft);

754

755

if (fm < ((qn>>1)*((qn>>1) + 1)>>1))

756

{

757

itheta = (isqrt32(8*(opus_uint32)fm + 1) - 1)>>1;

758

fs = itheta + 1;

759

fl = itheta*(itheta + 1)>>1;

760

}

761

else

762

{

763

itheta = (2*(qn + 1)

764

- isqrt32(8*(opus_uint32)(ft - fm - 1) + 1))>>1;

765

fs = qn + 1 - itheta;

766

fl = ft - ((qn + 1 - itheta)*(qn + 2 - itheta)>>1);

767

}

768

769

ec_dec_update(ec, fl, fl+fs, ft);

770

}

771

}

772

itheta = (opus_int32)itheta*16384/qn;

773

if (encode && stereo)

774

{

775

if (itheta==0)

776

intensity_stereo(m, X, Y, bandE, i, N);

777

else

778

stereo_split(X, Y, N);

779

}

780

/* NOTE: Renormalising X and Y *may* help fixed-point a bit at very high rate.

781

Let's do that at higher complexity */

782

} else if (stereo) {

783

if (encode)

784

{

785

inv = itheta > 8192;

786

if (inv)

787

{

788

int j;

789

for (j=0;j<N;j++)

790

Y[j] = -Y[j];

791

}

792

intensity_stereo(m, X, Y, bandE, i, N);

793

}

794

if (*b>2<<BITRES3 && ctx->remaining_bits > 2<<BITRES3)

795

{

796

if (encode)

797

ec_enc_bit_logp(ec, inv, 2);

798

else

799

inv = ec_dec_bit_logp(ec, 2);

800

} else

801

inv = 0;

802

itheta = 0;

803

}

804

qalloc = ec_tell_frac(ec) - tell;

805

*b -= qalloc;

806

807

if (itheta == 0)

808

{

809

imid = 32767;

810

iside = 0;

811

*fill &= (1<<B)-1;

812

delta = -16384;

813

} else if (itheta == 16384)

814

{

815

imid = 0;

816

iside = 32767;

817

*fill &= ((1<<B)-1)<<B;

818

delta = 16384;

819

} else {

820

imid = bitexact_cos((opus_int16)itheta);

821

iside = bitexact_cos((opus_int16)(16384-itheta));

822

/* This is the mid vs side allocation that minimizes squared error

823

in that band. */

824

delta = FRAC_MUL16((N-1)<<7,bitexact_log2tan(iside,imid))((16384+((opus_int32)(opus_int16)((N-1)<<7)*(opus_int16
)(bitexact_log2tan(iside,imid))))>>15);

825

}

826

827

sctx->inv = inv;

828

sctx->imid = imid;

829

sctx->iside = iside;

830

sctx->delta = delta;

831

sctx->itheta = itheta;

832

sctx->qalloc = qalloc;

833

}

834

static unsigned quant_band_n1(struct band_ctx *ctx, celt_norm *X, celt_norm *Y, int b,

835

celt_norm *lowband_out)

836

{

837

#ifdef RESYNTH

838

int resynth = 1;

839

#else

840

int resynth = !ctx->encode;

841

#endif

842

int c;

843

int stereo;

844

celt_norm *x = X;

845

int encode;

846

ec_ctx *ec;

847

848

encode = ctx->encode;

849

ec = ctx->ec;

850

851

stereo = Y != NULL((void*)0);

852

c=0; do {

853

int sign=0;

854

if (ctx->remaining_bits>=1<<BITRES3)

855

{

856

if (encode)

857

{

858

sign = x[0]<0;

859

ec_enc_bits(ec, sign, 1);

860

} else {

861

sign = ec_dec_bits(ec, 1);

862

}

863

ctx->remaining_bits -= 1<<BITRES3;

864

b-=1<<BITRES3;

865

}

866

if (resynth)

867

x[0] = sign ? -NORM_SCALING1.f : NORM_SCALING1.f;

868

x = Y;

869

} while (++c<1+stereo);

870

if (lowband_out)

871

lowband_out[0] = SHR16(X[0],4)(X[0]);

872

return 1;

873

}

874

875

/* This function is responsible for encoding and decoding a mono partition.

876

It can split the band in two and transmit the energy difference with

877

the two half-bands. It can be called recursively so bands can end up being

878

split in 8 parts. */

879

static unsigned quant_partition(struct band_ctx *ctx, celt_norm *X,

880

int N, int b, int B, celt_norm *lowband,

881

int LM,

882

opus_val16 gain, int fill)

883

{

884

const unsigned char *cache;

885

int q;

886

int curr_bits;

887

int imid=0, iside=0;

888

int B0=B;

889

opus_val16 mid=0, side=0;

890

unsigned cm=0;

891

#ifdef RESYNTH

892

int resynth = 1;

893

#else

894

int resynth = !ctx->encode;

895

#endif

896

celt_norm *Y=NULL((void*)0);

897

int encode;

898

const CELTModeOpusCustomMode *m;

899

int i;

900

int spread;

901

ec_ctx *ec;

902

903

encode = ctx->encode;

904

m = ctx->m;

905

i = ctx->i;

906

spread = ctx->spread;

907

ec = ctx->ec;

908

909

/* If we need 1.5 more bit than we can produce, split the band in two. */

910

cache = m->cache.bits + m->cache.index[(LM+1)*m->nbEBands+i];

911

if (LM != -1 && b > cache[cache[0]]+12 && N>2)

912

{

913

int mbits, sbits, delta;

914

int itheta;

915

int qalloc;

916

struct split_ctx sctx;

917

celt_norm *next_lowband2=NULL((void*)0);

918

opus_int32 rebalance;

919

920

N >>= 1;

921

Y = X+N;

922

LM -= 1;

923

if (B==1)

924

fill = (fill&1)|(fill<<1);

925

B = (B+1)>>1;

926

927

compute_theta(ctx, &sctx, X, Y, N, &b, B, B0,

928

LM, 0, &fill);

929

imid = sctx.imid;

930

iside = sctx.iside;

931

delta = sctx.delta;

932

itheta = sctx.itheta;

933

qalloc = sctx.qalloc;

934

#ifdef FIXED_POINT

935

mid = imid;

936

side = iside;

937

#else

938

mid = (1.f/32768)*imid;

939

side = (1.f/32768)*iside;

940

#endif

941

942

/* Give more bits to low-energy MDCTs than they would otherwise deserve */

943

if (B0>1 && (itheta&0x3fff))

944

{

945

if (itheta > 8192)

946

/* Rough approximation for pre-echo masking */

947

delta -= delta>>(4-LM);

948

else

949

/* Corresponds to a forward-masking slope of 1.5 dB per 10 ms */

950

delta = IMIN(0, delta + (N<<BITRES>>(5-LM)))((0) < (delta + (N<<3>>(5-LM))) ? (0) : (delta
+ (N<<3>>(5-LM))));

951

}

952

mbits = IMAX(0, IMIN(b, (b-delta)/2))((0) > (((b) < ((b-delta)/2) ? (b) : ((b-delta)/2))) ? (
0) : (((b) < ((b-delta)/2) ? (b) : ((b-delta)/2))));

953

sbits = b-mbits;

954

ctx->remaining_bits -= qalloc;

955

956

if (lowband)

957

next_lowband2 = lowband+N; /* >32-bit split case */

958

959

rebalance = ctx->remaining_bits;

960

if (mbits >= sbits)

961

{

962

cm = quant_partition(ctx, X, N, mbits, B,

963

lowband, LM,

964

MULT16_16_P15(gain,mid)((gain)*(mid)), fill);

965

rebalance = mbits - (rebalance-ctx->remaining_bits);

966

if (rebalance > 3<<BITRES3 && itheta!=0)

967

sbits += rebalance - (3<<BITRES3);

968

cm |= quant_partition(ctx, Y, N, sbits, B,

969

next_lowband2, LM,

970

MULT16_16_P15(gain,side)((gain)*(side)), fill>>B)<<(B0>>1);

971

} else {

972

cm = quant_partition(ctx, Y, N, sbits, B,

973

next_lowband2, LM,

974

MULT16_16_P15(gain,side)((gain)*(side)), fill>>B)<<(B0>>1);

975

rebalance = sbits - (rebalance-ctx->remaining_bits);

976

if (rebalance > 3<<BITRES3 && itheta!=16384)

977

mbits += rebalance - (3<<BITRES3);

978

cm |= quant_partition(ctx, X, N, mbits, B,

979

lowband, LM,

980

MULT16_16_P15(gain,mid)((gain)*(mid)), fill);

981

}

982

} else {

983

/* This is the basic no-split case */

984

q = bits2pulses(m, i, LM, b);

985

curr_bits = pulses2bits(m, i, LM, q);

986

ctx->remaining_bits -= curr_bits;

987

988

/* Ensures we can never bust the budget */

989

while (ctx->remaining_bits < 0 && q > 0)

990

{

991

ctx->remaining_bits += curr_bits;

992

q--;

993

curr_bits = pulses2bits(m, i, LM, q);

994

ctx->remaining_bits -= curr_bits;

995

}

996

997

if (q!=0)

998

{

999

int K = get_pulses(q);

1000

1001

/* Finally do the actual quantization */

1002

if (encode)

1003

{

1004

cm = alg_quant(X, N, K, spread, B, ec

1005

#ifdef RESYNTH

1006

, gain

1007

#endif

1008

);

1009

} else {

1010

cm = alg_unquant(X, N, K, spread, B, ec, gain);

1011

}

1012

} else {

1013

/* If there's no pulse, fill the band anyway */

1014

int j;

1015

if (resynth)

1016

{

1017

unsigned cm_mask;

1018

/* B can be as large as 16, so this shift might overflow an int on a

1019

16-bit platform; use a long to get defined behavior.*/

1020

cm_mask = (unsigned)(1UL<<B)-1;

1021

fill &= cm_mask;

1022

if (!fill)

1023

{

1024

for (j=0;j<N;j++)

1025

X[j] = 0;

1026

} else {

1027

if (lowband == NULL((void*)0))

1028

{

1029

/* Noise */

1030

for (j=0;j<N;j++)

1031

{

1032

ctx->seed = celt_lcg_rand(ctx->seed);

1033

X[j] = (celt_norm)((opus_int32)ctx->seed>>20);

1034

}

1035

cm = cm_mask;

1036

} else {

1037

/* Folded spectrum */

1038

for (j=0;j<N;j++)

1039

{

1040

opus_val16 tmp;

1041

ctx->seed = celt_lcg_rand(ctx->seed);

1042

/* About 48 dB below the "normal" folding level */

1043

tmp = QCONST16(1.0f/256, 10)(1.0f/256);

1044

tmp = (ctx->seed)&0x8000 ? tmp : -tmp;

1045

X[j] = lowband[j]+tmp;

1046

}

1047

cm = fill;

1048

}

1049

renormalise_vector(X, N, gain);

1050

}

1051

}

1052

}

1053

}

1054

1055

return cm;

1056

}

1057

1058

1059

/* This function is responsible for encoding and decoding a band for the mono case. */

1060

static unsigned quant_band(struct band_ctx *ctx, celt_norm *X,

1061

int N, int b, int B, celt_norm *lowband,

1062

int LM, celt_norm *lowband_out,

1063

opus_val16 gain, celt_norm *lowband_scratch, int fill)

1064

{

1065

int N0=N;

1066

int N_B=N;

1067

int N_B0;

1068

int B0=B;

1069

int time_divide=0;

1070

int recombine=0;

1071

int longBlocks;

1072

unsigned cm=0;

1073

#ifdef RESYNTH

1074

int resynth = 1;

1075

#else

1076

int resynth = !ctx->encode;

1077

#endif

1078

int k;

1079

int encode;

1080

int tf_change;

1081

1082

encode = ctx->encode;

1083

tf_change = ctx->tf_change;

1084

1085

longBlocks = B0==1;

Assuming 'B0' is not equal to 1

→

1086

1087

N_B /= B;

1088

1089

/* Special case for one sample */

1090

if (N==1)

←

Assuming 'N' is not equal to 1

→

←

Taking false branch

→

1091

{

1092

return quant_band_n1(ctx, X, NULL((void*)0), b, lowband_out);

1093

}

1094

1095

if (tf_change>0)

←

Assuming 'tf_change' is <= 0

→

←

Taking false branch

→

1096

recombine = tf_change;

1097

/* Band recombining to increase frequency resolution */

1098

1099

if (lowband_scratch && lowband && (recombine || ((N_B&1) == 0 && tf_change<0) || B0>1))

1100

{

1101

int j;

1102

for (j=0;j<N;j++)

1103

lowband_scratch[j] = lowband[j];

1104

lowband = lowband_scratch;

1105

}

1106

1107

for (k=0;k<recombine;k++)

←

Loop condition is false. Execution continues on line 1118

→

1108

{

1109

static const unsigned char bit_interleave_table[16]={

1110

0,1,1,1,2,3,3,3,2,3,3,3,2,3,3,3

1111

};

1112

if (encode)

1113

haar1(X, N>>k, 1<<k);

1114

if (lowband)

1115

haar1(lowband, N>>k, 1<<k);

1116

fill = bit_interleave_table[fill&0xF]|bit_interleave_table[fill>>4]<<2;

1117

}

1118

B>>=recombine;

1119

N_B<<=recombine;

1120

1121

/* Increasing the time resolution */

1122

while ((N_B&1) == 0 && tf_change<0)

1123

{

1124

if (encode)

1125

haar1(X, N_B, B);

1126

if (lowband)

1127

haar1(lowband, N_B, B);

1128

fill |= fill<<B;

1129

B <<= 1;

1130

N_B >>= 1;

1131

time_divide++;

1132

tf_change++;

1133

}

1134

B0=B;

1135

N_B0 = N_B;

1136

1137

/* Reorganize the samples in time order instead of frequency order */

1138

if (B0>1)

Assuming 'B0' is > 1

Taking true branch

1139

{

1140

if (encode)

←

Taking false branch

→

1141

deinterleave_hadamard(X, N_B>>recombine, B0<<recombine, longBlocks);

1142

if (lowband)

←

Assuming 'lowband' is null

→

←

Taking false branch

→

1143

deinterleave_hadamard(lowband, N_B>>recombine, B0<<recombine, longBlocks);

1144

}

1145

1146

cm = quant_partition(ctx, X, N, b, B, lowband,

1147

LM, gain, fill);

1148

1149

/* This code is used by the decoder and by the resynthesis-enabled encoder */

1150

if (resynth)

←

Taking true branch

→

1151

{

1152

/* Undo the sample reorganization going from time order to frequency order */

1153

if (B0>1)

←

Taking true branch

→

1154

interleave_hadamard(X, N_B>>recombine, B0<<recombine, longBlocks);

←

Calling 'interleave_hadamard'

→

1155

1156

/* Undo time-freq changes that we did earlier */

1157

N_B = N_B0;

1158

B = B0;

1159

for (k=0;k<time_divide;k++)

1160

{

1161

B >>= 1;

1162

N_B <<= 1;

1163

cm |= cm>>B;

1164

haar1(X, N_B, B);

1165

}

1166

1167

for (k=0;k<recombine;k++)

1168

{

1169

static const unsigned char bit_deinterleave_table[16]={

1170

0x00,0x03,0x0C,0x0F,0x30,0x33,0x3C,0x3F,

1171

0xC0,0xC3,0xCC,0xCF,0xF0,0xF3,0xFC,0xFF

1172

};

1173

cm = bit_deinterleave_table[cm];

1174

haar1(X, N0>>k, 1<<k);

1175

}

1176

B<<=recombine;

1177

1178

/* Scale output for later folding */

1179

if (lowband_out)

1180

{

1181

int j;

1182

opus_val16 n;

1183

n = celt_sqrt(SHL32(EXTEND32(N0),22))((float)sqrt(((N0))));

1184

for (j=0;j<N0;j++)

1185

lowband_out[j] = MULT16_16_Q15(n,X[j])((n)*(X[j]));

1186

}

1187

cm &= (1<<B)-1;

1188

}

1189

return cm;

1190

}

1191

1192

1193

/* This function is responsible for encoding and decoding a band for the stereo case. */

1194

static unsigned quant_band_stereo(struct band_ctx *ctx, celt_norm *X, celt_norm *Y,

1195

int N, int b, int B, celt_norm *lowband,

1196

int LM, celt_norm *lowband_out,

1197

celt_norm *lowband_scratch, int fill)

1198

{

1199

int imid=0, iside=0;

1200

int inv = 0;

1201

opus_val16 mid=0, side=0;

1202

unsigned cm=0;

1203

#ifdef RESYNTH

1204

int resynth = 1;

1205

#else

1206

int resynth = !ctx->encode;

1207

#endif

1208

int mbits, sbits, delta;

1209

int itheta;

1210

int qalloc;

1211

struct split_ctx sctx;

1212

int orig_fill;

1213

int encode;

1214

ec_ctx *ec;

1215

1216

encode = ctx->encode;

1217

ec = ctx->ec;

1218

1219

/* Special case for one sample */

1220

if (N==1)

1221

{

1222

return quant_band_n1(ctx, X, Y, b, lowband_out);

1223

}

1224

1225

orig_fill = fill;

1226

1227

compute_theta(ctx, &sctx, X, Y, N, &b, B, B,

1228

LM, 1, &fill);

1229

inv = sctx.inv;

1230

imid = sctx.imid;

1231

iside = sctx.iside;

1232

delta = sctx.delta;

1233

itheta = sctx.itheta;

1234

qalloc = sctx.qalloc;

1235

#ifdef FIXED_POINT

1236

mid = imid;

1237

side = iside;

1238

#else

1239

mid = (1.f/32768)*imid;

1240

side = (1.f/32768)*iside;

1241

#endif

1242

1243

/* This is a special case for N=2 that only works for stereo and takes

1244

advantage of the fact that mid and side are orthogonal to encode

1245

the side with just one bit. */

1246

if (N==2)

1247

{

1248

int c;

1249

int sign=0;

1250

celt_norm *x2, *y2;

1251

mbits = b;

1252

sbits = 0;

1253

/* Only need one bit for the side. */

1254

if (itheta != 0 && itheta != 16384)

1255

sbits = 1<<BITRES3;

1256

mbits -= sbits;

1257

c = itheta > 8192;

1258

ctx->remaining_bits -= qalloc+sbits;

1259

1260

x2 = c ? Y : X;

1261

y2 = c ? X : Y;

1262

if (sbits)

1263

{

1264

if (encode)

1265

{

1266

/* Here we only need to encode a sign for the side. */

1267

sign = x2[0]*y2[1] - x2[1]*y2[0] < 0;

1268

ec_enc_bits(ec, sign, 1);

1269

} else {

1270

sign = ec_dec_bits(ec, 1);

1271

}

1272

}

1273

sign = 1-2*sign;

1274

/* We use orig_fill here because we want to fold the side, but if

1275

itheta==16384, we'll have cleared the low bits of fill. */

1276

cm = quant_band(ctx, x2, N, mbits, B, lowband,

1277

LM, lowband_out, Q15ONE1.0f, lowband_scratch, orig_fill);

1278

/* We don't split N=2 bands, so cm is either 1 or 0 (for a fold-collapse),

1279

and there's no need to worry about mixing with the other channel. */

1280

y2[0] = -sign*x2[1];

1281

y2[1] = sign*x2[0];

1282

if (resynth)

1283

{

1284

celt_norm tmp;

1285

X[0] = MULT16_16_Q15(mid, X[0])((mid)*(X[0]));

1286

X[1] = MULT16_16_Q15(mid, X[1])((mid)*(X[1]));

1287

Y[0] = MULT16_16_Q15(side, Y[0])((side)*(Y[0]));

1288

Y[1] = MULT16_16_Q15(side, Y[1])((side)*(Y[1]));

1289

tmp = X[0];

1290

X[0] = SUB16(tmp,Y[0])((tmp)-(Y[0]));

1291

Y[0] = ADD16(tmp,Y[0])((tmp)+(Y[0]));

1292

tmp = X[1];

1293

X[1] = SUB16(tmp,Y[1])((tmp)-(Y[1]));

1294

Y[1] = ADD16(tmp,Y[1])((tmp)+(Y[1]));

1295

}

1296

} else {

1297

/* "Normal" split code */

1298

opus_int32 rebalance;

1299

1300

mbits = IMAX(0, IMIN(b, (b-delta)/2))((0) > (((b) < ((b-delta)/2) ? (b) : ((b-delta)/2))) ? (
0) : (((b) < ((b-delta)/2) ? (b) : ((b-delta)/2))));

1301

sbits = b-mbits;

1302

ctx->remaining_bits -= qalloc;

1303

1304

rebalance = ctx->remaining_bits;

1305

if (mbits >= sbits)

1306

{

1307

/* In stereo mode, we do not apply a scaling to the mid because we need the normalized

1308

mid for folding later. */

1309

cm = quant_band(ctx, X, N, mbits, B,

1310

lowband, LM, lowband_out,

1311

Q15ONE1.0f, lowband_scratch, fill);

1312

rebalance = mbits - (rebalance-ctx->remaining_bits);

1313

if (rebalance > 3<<BITRES3 && itheta!=0)

1314

sbits += rebalance - (3<<BITRES3);

1315

1316

/* For a stereo split, the high bits of fill are always zero, so no

1317

folding will be done to the side. */

1318

cm |= quant_band(ctx, Y, N, sbits, B,

1319

NULL((void*)0), LM, NULL((void*)0),

1320

side, NULL((void*)0), fill>>B);

1321

} else {

1322

/* For a stereo split, the high bits of fill are always zero, so no

1323

folding will be done to the side. */

1324

cm = quant_band(ctx, Y, N, sbits, B,

1325

NULL((void*)0), LM, NULL((void*)0),

1326

side, NULL((void*)0), fill>>B);

1327

rebalance = sbits - (rebalance-ctx->remaining_bits);

1328

if (rebalance > 3<<BITRES3 && itheta!=16384)

1329

mbits += rebalance - (3<<BITRES3);

1330

/* In stereo mode, we do not apply a scaling to the mid because we need the normalized

1331

mid for folding later. */

1332

cm |= quant_band(ctx, X, N, mbits, B,

1333

lowband, LM, lowband_out,

1334

Q15ONE1.0f, lowband_scratch, fill);

1335

}

1336

}

1337

1338

1339

/* This code is used by the decoder and by the resynthesis-enabled encoder */

1340

if (resynth)

1341

{

1342

if (N!=2)

1343

stereo_merge(X, Y, mid, N);

1344

if (inv)

1345

{

1346

int j;

1347

for (j=0;j<N;j++)

1348

Y[j] = -Y[j];

1349

}

1350

}

1351

return cm;

1352

}

1353

1354

1355

void quant_all_bands(int encode, const CELTModeOpusCustomMode *m, int start, int end,

1356

celt_norm *X_, celt_norm *Y_, unsigned char *collapse_masks, const celt_ener *bandE, int *pulses,

1357

int shortBlocks, int spread, int dual_stereo, int intensity, int *tf_res,

1358

opus_int32 total_bits, opus_int32 balance, ec_ctx *ec, int LM, int codedBands, opus_uint32 *seed)

1359

{

1360

int i;

1361

opus_int32 remaining_bits;

1362

const opus_int16 * OPUS_RESTRICT__restrict eBands = m->eBands;

1363

celt_norm * OPUS_RESTRICT__restrict norm, * OPUS_RESTRICT__restrict norm2;

1364

VARDECL(celt_norm, _norm);

1365

celt_norm *lowband_scratch;

1366

int B;

1367

int M;

1368

int lowband_offset;

1369

int update_lowband = 1;

1370

int C = Y_ != NULL((void*)0) ? 2 : 1;

1371

int norm_offset;

1372

#ifdef RESYNTH

1373

int resynth = 1;

1374

#else

1375

int resynth = !encode;

1376

#endif

1377

struct band_ctx ctx;

1378

SAVE_STACK;

1379

1380

M = 1<<LM;

1381

B = shortBlocks ? M : 1;

1382

norm_offset = M*eBands[start];

1383

/* No need to allocate norm for the last band because we don't need an

1384

output in that band. */

1385

ALLOC(_norm, C*(M*eBands[m->nbEBands-1]-norm_offset), celt_norm)celt_norm _norm[C*(M*eBands[m->nbEBands-1]-norm_offset)];

1386

norm = _norm;

1387

norm2 = norm + M*eBands[m->nbEBands-1]-norm_offset;

1388

/* We can use the last band as scratch space because we don't need that

1389

scratch space for the last band. */

1390

lowband_scratch = X_+M*eBands[m->nbEBands-1];

1391

1392

lowband_offset = 0;

1393

ctx.bandE = bandE;

1394

ctx.ec = ec;

1395

ctx.encode = encode;

1396

ctx.intensity = intensity;

1397

ctx.m = m;

1398

ctx.seed = *seed;

1399

ctx.spread = spread;

1400

for (i=start;i<end;i++)

1401

{

1402

opus_int32 tell;

1403

int b;

1404

int N;

1405

opus_int32 curr_balance;

1406

int effective_lowband=-1;

1407

celt_norm * OPUS_RESTRICT__restrict X, * OPUS_RESTRICT__restrict Y;

1408

int tf_change=0;

1409

unsigned x_cm;

1410

unsigned y_cm;

1411

int last;

1412

1413

ctx.i = i;

1414

last = (i==end-1);

1415

1416

X = X_+M*eBands[i];

1417

if (Y_!=NULL((void*)0))

1418

Y = Y_+M*eBands[i];

1419

else

1420

Y = NULL((void*)0);

1421

N = M*eBands[i+1]-M*eBands[i];

1422

tell = ec_tell_frac(ec);

1423

1424

/* Compute how many bits we want to allocate to this band */

1425

if (i != start)

1426

balance -= tell;

1427

remaining_bits = total_bits-tell-1;

1428

ctx.remaining_bits = remaining_bits;

1429

if (i <= codedBands-1)

1430

{

1431

curr_balance = balance / IMIN(3, codedBands-i)((3) < (codedBands-i) ? (3) : (codedBands-i));

1432

b = IMAX(0, IMIN(16383, IMIN(remaining_bits+1,pulses[i]+curr_balance)))((0) > (((16383) < (((remaining_bits+1) < (pulses[i]
+curr_balance) ? (remaining_bits+1) : (pulses[i]+curr_balance
))) ? (16383) : (((remaining_bits+1) < (pulses[i]+curr_balance
) ? (remaining_bits+1) : (pulses[i]+curr_balance))))) ? (0) :
(((16383) < (((remaining_bits+1) < (pulses[i]+curr_balance
) ? (remaining_bits+1) : (pulses[i]+curr_balance))) ? (16383)
: (((remaining_bits+1) < (pulses[i]+curr_balance) ? (remaining_bits
+1) : (pulses[i]+curr_balance))))));

1433

} else {

1434

b = 0;

1435

}

1436

1437

if (resynth && M*eBands[i]-N >= M*eBands[start] && (update_lowband || lowband_offset==0))

1438

lowband_offset = i;

1439

1440

tf_change = tf_res[i];

1441

ctx.tf_change = tf_change;

1442

if (i>=m->effEBands)

1443

{

1444

X=norm;

1445

if (Y_!=NULL((void*)0))

1446

Y = norm;

1447

lowband_scratch = NULL((void*)0);

1448

}

1449

if (i==end-1)

1450

lowband_scratch = NULL((void*)0);

1451

1452

/* Get a conservative estimate of the collapse_mask's for the bands we're

1453

going to be folding from. */

1454

if (lowband_offset != 0 && (spread!=SPREAD_AGGRESSIVE(3) || B>1 || tf_change<0))

1455

{

1456

int fold_start;

1457

int fold_end;

1458

int fold_i;

1459

/* This ensures we never repeat spectral content within one band */

1460

effective_lowband = IMAX(0, M*eBands[lowband_offset]-norm_offset-N)((0) > (M*eBands[lowband_offset]-norm_offset-N) ? (0) : (M
*eBands[lowband_offset]-norm_offset-N));

1461

fold_start = lowband_offset;

1462

while(M*eBands[--fold_start] > effective_lowband+norm_offset);

1463

fold_end = lowband_offset-1;

1464

while(M*eBands[++fold_end] < effective_lowband+norm_offset+N);

1465

x_cm = y_cm = 0;

1466

fold_i = fold_start; do {

1467

x_cm |= collapse_masks[fold_i*C+0];

1468

y_cm |= collapse_masks[fold_i*C+C-1];

1469

} while (++fold_i<fold_end);

1470

}

1471

/* Otherwise, we'll be using the LCG to fold, so all blocks will (almost

1472

always) be non-zero. */

1473

else

1474

x_cm = y_cm = (1<<B)-1;

1475

1476

if (dual_stereo && i==intensity)

1477

{

1478

int j;

1479

1480

/* Switch off dual stereo to do intensity. */

1481

dual_stereo = 0;

1482

if (resynth)

1483

for (j=0;j<M*eBands[i]-norm_offset;j++)

1484

norm[j] = HALF32(norm[j]+norm2[j])(.5f*(norm[j]+norm2[j]));

1485

}

1486

if (dual_stereo)

1487

{

1488

x_cm = quant_band(&ctx, X, N, b/2, B,

1489

effective_lowband != -1 ? norm+effective_lowband : NULL((void*)0), LM,

1490

last?NULL((void*)0):norm+M*eBands[i]-norm_offset, Q15ONE1.0f, lowband_scratch, x_cm);

1491

y_cm = quant_band(&ctx, Y, N, b/2, B,

1492

effective_lowband != -1 ? norm2+effective_lowband : NULL((void*)0), LM,

1493

last?NULL((void*)0):norm2+M*eBands[i]-norm_offset, Q15ONE1.0f, lowband_scratch, y_cm);

1494

} else {

1495

if (Y!=NULL((void*)0))

1496

{

1497

x_cm = quant_band_stereo(&ctx, X, Y, N, b, B,

1498

effective_lowband != -1 ? norm+effective_lowband : NULL((void*)0), LM,

1499

last?NULL((void*)0):norm+M*eBands[i]-norm_offset, lowband_scratch, x_cm|y_cm);

1500

} else {

1501

x_cm = quant_band(&ctx, X, N, b, B,

1502

effective_lowband != -1 ? norm+effective_lowband : NULL((void*)0), LM,

1503

last?NULL((void*)0):norm+M*eBands[i]-norm_offset, Q15ONE1.0f, lowband_scratch, x_cm|y_cm);

1504

}

1505

y_cm = x_cm;

1506

}

1507

collapse_masks[i*C+0] = (unsigned char)x_cm;

1508

collapse_masks[i*C+C-1] = (unsigned char)y_cm;

1509

balance += pulses[i] + tell;

1510

1511

/* Update the folding position only as long as we have 1 bit/sample depth. */

1512

update_lowband = b>(N<<BITRES3);

1513

}

1514

*seed = ctx.seed;

1515

1516

RESTORE_STACK;

1517

}

1518