pickrst.c 44.6 KB
Newer Older
1
/*
2
 * Copyright (c) 2016, Alliance for Open Media. All rights reserved
3
 *
4
5
6
7
8
9
 * This source code is subject to the terms of the BSD 2 Clause License and
 * the Alliance for Open Media Patent License 1.0. If the BSD 2 Clause License
 * was not distributed with this source code in the LICENSE file, you can
 * obtain it at www.aomedia.org/license/software. If the Alliance for Open
 * Media Patent License 1.0 was not distributed with this source code in the
 * PATENTS file, you can obtain it at www.aomedia.org/license/patent.
10
11
12
 */

#include <assert.h>
13
#include <float.h>
14
15
16
#include <limits.h>
#include <math.h>

Yaowu Xu's avatar
Yaowu Xu committed
17
#include "./aom_scale_rtcd.h"
18

19
#include "aom_dsp/psnr.h"
Yaowu Xu's avatar
Yaowu Xu committed
20
21
#include "aom_dsp/aom_dsp_common.h"
#include "aom_mem/aom_mem.h"
22
#include "aom_ports/mem.h"
23

24
25
#include "av1/common/onyxc_int.h"
#include "av1/common/quant_common.h"
26

27
28
29
30
#include "av1/encoder/encoder.h"
#include "av1/encoder/picklpf.h"
#include "av1/encoder/pickrst.h"
#include "av1/encoder/quantize.h"
31

32
33
34
typedef double (*search_restore_type)(const YV12_BUFFER_CONFIG *src,
                                      AV1_COMP *cpi, int filter_level,
                                      int partial_frame, RestorationInfo *info,
35
36
                                      double *best_tile_cost,
                                      YV12_BUFFER_CONFIG *dst_frame);
37

38
const int frame_level_restore_bits[RESTORE_TYPES] = { 2, 2, 3, 3, 2 };
39
40

static int64_t sse_restoration_tile(const YV12_BUFFER_CONFIG *src,
41
42
                                    const YV12_BUFFER_CONFIG *dst,
                                    const AV1_COMMON *cm, int h_start,
43
44
                                    int width, int v_start, int height,
                                    int y_only) {
45
46
47
  int64_t filt_err;
#if CONFIG_AOM_HIGHBITDEPTH
  if (cm->use_highbitdepth) {
48
49
    filt_err =
        aom_highbd_get_y_sse_part(src, dst, h_start, width, v_start, height);
50
51
52
53
54
55
56
57
58
    if (!y_only) {
      filt_err += aom_highbd_get_u_sse_part(
          src, dst, h_start >> cm->subsampling_x, width >> cm->subsampling_x,
          v_start >> cm->subsampling_y, height >> cm->subsampling_y);
      filt_err += aom_highbd_get_v_sse_part(
          src, dst, h_start >> cm->subsampling_x, width >> cm->subsampling_x,
          v_start >> cm->subsampling_y, height >> cm->subsampling_y);
    }
    return filt_err;
59
60
  }
#endif  // CONFIG_AOM_HIGHBITDEPTH
61
62
63
64
65
66
67
68
69
  filt_err = aom_get_y_sse_part(src, dst, h_start, width, v_start, height);
  if (!y_only) {
    filt_err += aom_get_u_sse_part(
        src, dst, h_start >> cm->subsampling_x, width >> cm->subsampling_x,
        v_start >> cm->subsampling_y, height >> cm->subsampling_y);
    filt_err += aom_get_u_sse_part(
        src, dst, h_start >> cm->subsampling_x, width >> cm->subsampling_x,
        v_start >> cm->subsampling_y, height >> cm->subsampling_y);
  }
70
71
72
73
74
  return filt_err;
}

static int64_t try_restoration_tile(const YV12_BUFFER_CONFIG *src,
                                    AV1_COMP *const cpi, RestorationInfo *rsi,
75
                                    int y_only, int partial_frame, int tile_idx,
76
77
                                    int subtile_idx, int subtile_bits,
                                    YV12_BUFFER_CONFIG *dst_frame) {
78
79
80
81
82
83
84
85
  AV1_COMMON *const cm = &cpi->common;
  int64_t filt_err;
  int tile_width, tile_height, nhtiles, nvtiles;
  int h_start, h_end, v_start, v_end;
  const int ntiles = av1_get_rest_ntiles(cm->width, cm->height, &tile_width,
                                         &tile_height, &nhtiles, &nvtiles);
  (void)ntiles;

86
  av1_loop_restoration_frame(cm->frame_to_show, cm, rsi, y_only, partial_frame,
87
                             dst_frame);
88
89
90
91
  av1_get_rest_tile_limits(tile_idx, subtile_idx, subtile_bits, nhtiles,
                           nvtiles, tile_width, tile_height, cm->width,
                           cm->height, 0, 0, &h_start, &h_end, &v_start,
                           &v_end);
92
  filt_err = sse_restoration_tile(src, dst_frame, cm, h_start, h_end - h_start,
93
                                  v_start, v_end - v_start, y_only);
94
95
96
97
98

  return filt_err;
}

static int64_t try_restoration_frame(const YV12_BUFFER_CONFIG *src,
Yaowu Xu's avatar
Yaowu Xu committed
99
                                     AV1_COMP *const cpi, RestorationInfo *rsi,
100
                                     int y_only, int partial_frame,
101
                                     YV12_BUFFER_CONFIG *dst_frame) {
Yaowu Xu's avatar
Yaowu Xu committed
102
  AV1_COMMON *const cm = &cpi->common;
103
  int64_t filt_err;
104
  av1_loop_restoration_frame(cm->frame_to_show, cm, rsi, y_only, partial_frame,
105
                             dst_frame);
Yaowu Xu's avatar
Yaowu Xu committed
106
#if CONFIG_AOM_HIGHBITDEPTH
107
  if (cm->use_highbitdepth) {
108
    filt_err = aom_highbd_get_y_sse(src, dst_frame);
109
110
111
112
113
    if (!y_only) {
      filt_err += aom_highbd_get_u_sse(src, dst_frame);
      filt_err += aom_highbd_get_v_sse(src, dst_frame);
    }
    return filt_err;
114
  }
Yaowu Xu's avatar
Yaowu Xu committed
115
#endif  // CONFIG_AOM_HIGHBITDEPTH
116
117
118
119
120
  filt_err = aom_get_y_sse(src, dst_frame);
  if (!y_only) {
    filt_err += aom_get_u_sse(src, dst_frame);
    filt_err += aom_get_v_sse(src, dst_frame);
  }
121
122
123
  return filt_err;
}

124
125
126
127
128
static int64_t get_pixel_proj_error(uint8_t *src8, int width, int height,
                                    int src_stride, uint8_t *dat8,
                                    int dat_stride, int bit_depth,
                                    int32_t *flt1, int flt1_stride,
                                    int32_t *flt2, int flt2_stride, int *xqd) {
129
130
131
132
  int i, j;
  int64_t err = 0;
  int xq[2];
  decode_xq(xqd, xq);
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
  if (bit_depth == 8) {
    const uint8_t *src = src8;
    const uint8_t *dat = dat8;
    for (i = 0; i < height; ++i) {
      for (j = 0; j < width; ++j) {
        const int32_t u =
            (int32_t)(dat[i * dat_stride + j] << SGRPROJ_RST_BITS);
        const int32_t f1 = (int32_t)flt1[i * flt1_stride + j] - u;
        const int32_t f2 = (int32_t)flt2[i * flt2_stride + j] - u;
        const int64_t v = xq[0] * f1 + xq[1] * f2 + (u << SGRPROJ_PRJ_BITS);
        const int32_t e =
            ROUND_POWER_OF_TWO(v, SGRPROJ_RST_BITS + SGRPROJ_PRJ_BITS) -
            src[i * src_stride + j];
        err += e * e;
      }
    }
  } else {
    const uint16_t *src = CONVERT_TO_SHORTPTR(src8);
    const uint16_t *dat = CONVERT_TO_SHORTPTR(dat8);
    for (i = 0; i < height; ++i) {
      for (j = 0; j < width; ++j) {
        const int32_t u =
            (int32_t)(dat[i * dat_stride + j] << SGRPROJ_RST_BITS);
        const int32_t f1 = (int32_t)flt1[i * flt1_stride + j] - u;
        const int32_t f2 = (int32_t)flt2[i * flt2_stride + j] - u;
        const int64_t v = xq[0] * f1 + xq[1] * f2 + (u << SGRPROJ_PRJ_BITS);
        const int32_t e =
            ROUND_POWER_OF_TWO(v, SGRPROJ_RST_BITS + SGRPROJ_PRJ_BITS) -
            src[i * src_stride + j];
        err += e * e;
      }
164
165
166
167
168
    }
  }
  return err;
}

169
170
171
172
static void get_proj_subspace(uint8_t *src8, int width, int height,
                              int src_stride, uint8_t *dat8, int dat_stride,
                              int bit_depth, int32_t *flt1, int flt1_stride,
                              int32_t *flt2, int flt2_stride, int *xq) {
173
174
175
176
177
178
179
180
181
  int i, j;
  double H[2][2] = { { 0, 0 }, { 0, 0 } };
  double C[2] = { 0, 0 };
  double Det;
  double x[2];
  const int size = width * height;

  xq[0] = -(1 << SGRPROJ_PRJ_BITS) / 4;
  xq[1] = (1 << SGRPROJ_PRJ_BITS) - xq[0];
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
  if (bit_depth == 8) {
    const uint8_t *src = src8;
    const uint8_t *dat = dat8;
    for (i = 0; i < height; ++i) {
      for (j = 0; j < width; ++j) {
        const double u = (double)(dat[i * dat_stride + j] << SGRPROJ_RST_BITS);
        const double s =
            (double)(src[i * src_stride + j] << SGRPROJ_RST_BITS) - u;
        const double f1 = (double)flt1[i * flt1_stride + j] - u;
        const double f2 = (double)flt2[i * flt2_stride + j] - u;
        H[0][0] += f1 * f1;
        H[1][1] += f2 * f2;
        H[0][1] += f1 * f2;
        C[0] += f1 * s;
        C[1] += f2 * s;
      }
    }
  } else {
    const uint16_t *src = CONVERT_TO_SHORTPTR(src8);
    const uint16_t *dat = CONVERT_TO_SHORTPTR(dat8);
    for (i = 0; i < height; ++i) {
      for (j = 0; j < width; ++j) {
        const double u = (double)(dat[i * dat_stride + j] << SGRPROJ_RST_BITS);
        const double s =
            (double)(src[i * src_stride + j] << SGRPROJ_RST_BITS) - u;
        const double f1 = (double)flt1[i * flt1_stride + j] - u;
        const double f2 = (double)flt2[i * flt2_stride + j] - u;
        H[0][0] += f1 * f1;
        H[1][1] += f2 * f2;
        H[0][1] += f1 * f2;
        C[0] += f1 * s;
        C[1] += f2 * s;
      }
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
    }
  }
  H[0][0] /= size;
  H[0][1] /= size;
  H[1][1] /= size;
  H[1][0] = H[0][1];
  C[0] /= size;
  C[1] /= size;
  Det = (H[0][0] * H[1][1] - H[0][1] * H[1][0]);
  if (Det < 1e-8) return;  // ill-posed, return default values
  x[0] = (H[1][1] * C[0] - H[0][1] * C[1]) / Det;
  x[1] = (H[0][0] * C[1] - H[1][0] * C[0]) / Det;
  xq[0] = (int)rint(x[0] * (1 << SGRPROJ_PRJ_BITS));
  xq[1] = (int)rint(x[1] * (1 << SGRPROJ_PRJ_BITS));
}

void encode_xq(int *xq, int *xqd) {
  xqd[0] = -xq[0];
  xqd[0] = clamp(xqd[0], SGRPROJ_PRJ_MIN0, SGRPROJ_PRJ_MAX0);
  xqd[1] = (1 << SGRPROJ_PRJ_BITS) + xqd[0] - xq[1];
  xqd[1] = clamp(xqd[1], SGRPROJ_PRJ_MIN1, SGRPROJ_PRJ_MAX1);
}

static void search_selfguided_restoration(uint8_t *dat8, int width, int height,
                                          int dat_stride, uint8_t *src8,
                                          int src_stride, int bit_depth,
241
242
                                          int *eps, int *xqd, int32_t *rstbuf) {
  int32_t *flt1 = rstbuf;
243
  int32_t *flt2 = flt1 + RESTORATION_TILEPELS_MAX;
244
  int32_t *tmpbuf2 = flt2 + RESTORATION_TILEPELS_MAX;
245
246
247
  int i, j, ep, bestep = 0;
  int64_t err, besterr = -1;
  int exqd[2], bestxqd[2] = { 0, 0 };
248

249
250
251
252
253
254
  for (ep = 0; ep < SGRPROJ_PARAMS; ep++) {
    int exq[2];
    if (bit_depth > 8) {
      uint16_t *dat = CONVERT_TO_SHORTPTR(dat8);
      for (i = 0; i < height; ++i) {
        for (j = 0; j < width; ++j) {
255
256
          flt1[i * width + j] = (int32_t)dat[i * dat_stride + j];
          flt2[i * width + j] = (int32_t)dat[i * dat_stride + j];
257
258
259
260
261
262
263
264
        }
      }
    } else {
      uint8_t *dat = dat8;
      for (i = 0; i < height; ++i) {
        for (j = 0; j < width; ++j) {
          const int k = i * width + j;
          const int l = i * dat_stride + j;
265
266
          flt1[k] = (int32_t)dat[l];
          flt2[k] = (int32_t)dat[l];
267
268
269
270
271
272
273
        }
      }
    }
    av1_selfguided_restoration(flt1, width, height, width, bit_depth,
                               sgr_params[ep].r1, sgr_params[ep].e1, tmpbuf2);
    av1_selfguided_restoration(flt2, width, height, width, bit_depth,
                               sgr_params[ep].r2, sgr_params[ep].e2, tmpbuf2);
274
275
    get_proj_subspace(src8, width, height, src_stride, dat8, dat_stride,
                      bit_depth, flt1, width, flt2, width, exq);
276
    encode_xq(exq, exqd);
277
278
279
    err =
        get_pixel_proj_error(src8, width, height, src_stride, dat8, dat_stride,
                             bit_depth, flt1, width, flt2, width, exqd);
280
281
282
283
284
285
286
287
288
289
290
291
292
293
    if (besterr == -1 || err < besterr) {
      bestep = ep;
      besterr = err;
      bestxqd[0] = exqd[0];
      bestxqd[1] = exqd[1];
    }
  }
  *eps = bestep;
  xqd[0] = bestxqd[0];
  xqd[1] = bestxqd[1];
}

static double search_sgrproj(const YV12_BUFFER_CONFIG *src, AV1_COMP *cpi,
                             int filter_level, int partial_frame,
294
295
                             RestorationInfo *info, double *best_tile_cost,
                             YV12_BUFFER_CONFIG *dst_frame) {
296
297
298
299
300
301
  SgrprojInfo *sgrproj_info = info->sgrproj_info;
  double err, cost_norestore, cost_sgrproj;
  int bits;
  MACROBLOCK *x = &cpi->td.mb;
  AV1_COMMON *const cm = &cpi->common;
  const YV12_BUFFER_CONFIG *dgd = cm->frame_to_show;
302
  RestorationInfo *rsi = &cpi->rst_search;
303
304
  int tile_idx, tile_width, tile_height, nhtiles, nvtiles;
  int h_start, h_end, v_start, v_end;
305
  // Allocate for the src buffer at high precision
306
307
308
309
310
311
312
313
  const int ntiles = av1_get_rest_ntiles(cm->width, cm->height, &tile_width,
                                         &tile_height, &nhtiles, &nvtiles);
  //  Make a copy of the unfiltered / processed recon buffer
  aom_yv12_copy_y(cm->frame_to_show, &cpi->last_frame_uf);
  av1_loop_filter_frame(cm->frame_to_show, cm, &cpi->td.mb.e_mbd, filter_level,
                        1, partial_frame);
  aom_yv12_copy_y(cm->frame_to_show, &cpi->last_frame_db);

314
  rsi->frame_restoration_type = RESTORE_SGRPROJ;
315
316

  for (tile_idx = 0; tile_idx < ntiles; ++tile_idx)
317
    rsi->sgrproj_info[tile_idx].level = 0;
318
319
320
321
322
  // Compute best Sgrproj filters for each tile
  for (tile_idx = 0; tile_idx < ntiles; ++tile_idx) {
    av1_get_rest_tile_limits(tile_idx, 0, 0, nhtiles, nvtiles, tile_width,
                             tile_height, cm->width, cm->height, 0, 0, &h_start,
                             &h_end, &v_start, &v_end);
323
    err = sse_restoration_tile(src, cm->frame_to_show, cm, h_start,
324
                               h_end - h_start, v_start, v_end - v_start, 1);
325
326
327
328
329
330
331
332
333
334
335
336
337
    // #bits when a tile is not restored
    bits = av1_cost_bit(RESTORE_NONE_SGRPROJ_PROB, 0);
    cost_norestore = RDCOST_DBL(x->rdmult, x->rddiv, (bits >> 4), err);
    best_tile_cost[tile_idx] = DBL_MAX;
    search_selfguided_restoration(
        dgd->y_buffer + v_start * dgd->y_stride + h_start, h_end - h_start,
        v_end - v_start, dgd->y_stride,
        src->y_buffer + v_start * src->y_stride + h_start, src->y_stride,
#if CONFIG_AOM_HIGHBITDEPTH
        cm->bit_depth,
#else
        8,
#endif  // CONFIG_AOM_HIGHBITDEPTH
338
        &rsi->sgrproj_info[tile_idx].ep, rsi->sgrproj_info[tile_idx].xqd,
339
        cm->rst_internal.tmpbuf);
340
341
    rsi->sgrproj_info[tile_idx].level = 1;
    err = try_restoration_tile(src, cpi, rsi, 1, partial_frame, tile_idx, 0, 0,
342
                               dst_frame);
343
344
345
346
347
348
    bits = SGRPROJ_BITS << AV1_PROB_COST_SHIFT;
    bits += av1_cost_bit(RESTORE_NONE_SGRPROJ_PROB, 1);
    cost_sgrproj = RDCOST_DBL(x->rdmult, x->rddiv, (bits >> 4), err);
    if (cost_sgrproj >= cost_norestore) {
      sgrproj_info[tile_idx].level = 0;
    } else {
349
      memcpy(&sgrproj_info[tile_idx], &rsi->sgrproj_info[tile_idx],
350
351
352
353
354
355
             sizeof(sgrproj_info[tile_idx]));
      bits = SGRPROJ_BITS << AV1_PROB_COST_SHIFT;
      best_tile_cost[tile_idx] = RDCOST_DBL(
          x->rdmult, x->rddiv,
          (bits + cpi->switchable_restore_cost[RESTORE_SGRPROJ]) >> 4, err);
    }
356
    rsi->sgrproj_info[tile_idx].level = 0;
357
358
  }
  // Cost for Sgrproj filtering
359
  bits = frame_level_restore_bits[rsi->frame_restoration_type]
360
361
362
363
         << AV1_PROB_COST_SHIFT;
  for (tile_idx = 0; tile_idx < ntiles; ++tile_idx) {
    bits +=
        av1_cost_bit(RESTORE_NONE_SGRPROJ_PROB, sgrproj_info[tile_idx].level);
364
    memcpy(&rsi->sgrproj_info[tile_idx], &sgrproj_info[tile_idx],
365
366
367
368
369
           sizeof(sgrproj_info[tile_idx]));
    if (sgrproj_info[tile_idx].level) {
      bits += (SGRPROJ_BITS << AV1_PROB_COST_SHIFT);
    }
  }
370
  err = try_restoration_frame(src, cpi, rsi, 1, partial_frame, dst_frame);
371
372
373
374
375
376
  cost_sgrproj = RDCOST_DBL(x->rdmult, x->rddiv, (bits >> 4), err);

  aom_yv12_copy_y(&cpi->last_frame_uf, cm->frame_to_show);
  return cost_sgrproj;
}

377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
static int64_t compute_sse(uint8_t *dgd, int width, int height, int dgd_stride,
                           uint8_t *src, int src_stride) {
  int64_t sse = 0;
  int i, j;
  for (i = 0; i < height; ++i) {
    for (j = 0; j < width; ++j) {
      const int diff =
          (int)dgd[i * dgd_stride + j] - (int)src[i * src_stride + j];
      sse += diff * diff;
    }
  }
  return sse;
}

#if CONFIG_AOM_HIGHBITDEPTH
static int64_t compute_sse_highbd(uint16_t *dgd, int width, int height,
                                  int dgd_stride, uint16_t *src,
                                  int src_stride) {
  int64_t sse = 0;
  int i, j;
  for (i = 0; i < height; ++i) {
    for (j = 0; j < width; ++j) {
      const int diff =
          (int)dgd[i * dgd_stride + j] - (int)src[i * src_stride + j];
      sse += diff * diff;
    }
  }
  return sse;
}
#endif  // CONFIG_AOM_HIGHBITDEPTH

static void search_domaintxfmrf_restoration(uint8_t *dgd8, int width,
                                            int height, int dgd_stride,
                                            uint8_t *src8, int src_stride,
411
                                            int bit_depth, int *sigma_r,
412
                                            uint8_t *fltbuf, int32_t *tmpbuf) {
413
414
415
416
417
  const int first_p_step = 8;
  const int second_p_range = first_p_step >> 1;
  const int second_p_step = 2;
  const int third_p_range = second_p_step >> 1;
  const int third_p_step = 1;
418
  int p, best_p0, best_p = -1;
419
420
  int64_t best_sse = INT64_MAX, sse;
  if (bit_depth == 8) {
421
    uint8_t *flt = fltbuf;
422
423
424
425
    uint8_t *dgd = dgd8;
    uint8_t *src = src8;
    // First phase
    for (p = first_p_step / 2; p < DOMAINTXFMRF_PARAMS; p += first_p_step) {
426
      av1_domaintxfmrf_restoration(dgd, width, height, dgd_stride, p, flt,
427
                                   width, tmpbuf);
428
      sse = compute_sse(flt, width, height, width, src, src_stride);
429
430
431
432
433
434
435
436
437
438
      if (sse < best_sse || best_p == -1) {
        best_p = p;
        best_sse = sse;
      }
    }
    // Second Phase
    best_p0 = best_p;
    for (p = best_p0 - second_p_range; p <= best_p0 + second_p_range;
         p += second_p_step) {
      if (p < 0 || p == best_p || p >= DOMAINTXFMRF_PARAMS) continue;
439
      av1_domaintxfmrf_restoration(dgd, width, height, dgd_stride, p, flt,
440
                                   width, tmpbuf);
441
      sse = compute_sse(flt, width, height, width, src, src_stride);
442
443
444
445
446
447
448
449
450
451
      if (sse < best_sse) {
        best_p = p;
        best_sse = sse;
      }
    }
    // Third Phase
    best_p0 = best_p;
    for (p = best_p0 - third_p_range; p <= best_p0 + third_p_range;
         p += third_p_step) {
      if (p < 0 || p == best_p || p >= DOMAINTXFMRF_PARAMS) continue;
452
      av1_domaintxfmrf_restoration(dgd, width, height, dgd_stride, p, flt,
453
                                   width, tmpbuf);
454
      sse = compute_sse(flt, width, height, width, src, src_stride);
455
456
457
458
459
460
461
      if (sse < best_sse) {
        best_p = p;
        best_sse = sse;
      }
    }
  } else {
#if CONFIG_AOM_HIGHBITDEPTH
462
    uint16_t *flt = (uint16_t *)fltbuf;
463
464
465
466
    uint16_t *dgd = CONVERT_TO_SHORTPTR(dgd8);
    uint16_t *src = CONVERT_TO_SHORTPTR(src8);
    // First phase
    for (p = first_p_step / 2; p < DOMAINTXFMRF_PARAMS; p += first_p_step) {
467
      av1_domaintxfmrf_restoration_highbd(dgd, width, height, dgd_stride, p,
468
469
                                          bit_depth, flt, width, tmpbuf);
      sse = compute_sse_highbd(flt, width, height, width, src, src_stride);
470
471
472
473
474
475
476
477
478
479
      if (sse < best_sse || best_p == -1) {
        best_p = p;
        best_sse = sse;
      }
    }
    // Second Phase
    best_p0 = best_p;
    for (p = best_p0 - second_p_range; p <= best_p0 + second_p_range;
         p += second_p_step) {
      if (p < 0 || p == best_p || p >= DOMAINTXFMRF_PARAMS) continue;
480
      av1_domaintxfmrf_restoration_highbd(dgd, width, height, dgd_stride, p,
481
482
                                          bit_depth, flt, width, tmpbuf);
      sse = compute_sse_highbd(flt, width, height, width, src, src_stride);
483
484
485
486
487
488
489
490
491
492
      if (sse < best_sse) {
        best_p = p;
        best_sse = sse;
      }
    }
    // Third Phase
    best_p0 = best_p;
    for (p = best_p0 - third_p_range; p <= best_p0 + third_p_range;
         p += third_p_step) {
      if (p < 0 || p == best_p || p >= DOMAINTXFMRF_PARAMS) continue;
493
      av1_domaintxfmrf_restoration_highbd(dgd, width, height, dgd_stride, p,
494
495
                                          bit_depth, flt, width, tmpbuf);
      sse = compute_sse_highbd(flt, width, height, width, src, src_stride);
496
497
498
499
500
501
502
503
504
505
506
507
508
509
      if (sse < best_sse) {
        best_p = p;
        best_sse = sse;
      }
    }
#else
    assert(0);
#endif  // CONFIG_AOM_HIGHBITDEPTH
  }
  *sigma_r = best_p;
}

static double search_domaintxfmrf(const YV12_BUFFER_CONFIG *src, AV1_COMP *cpi,
                                  int filter_level, int partial_frame,
510
511
                                  RestorationInfo *info, double *best_tile_cost,
                                  YV12_BUFFER_CONFIG *dst_frame) {
512
  DomaintxfmrfInfo *domaintxfmrf_info = info->domaintxfmrf_info;
513
514
  double cost_norestore, cost_domaintxfmrf;
  int64_t err;
515
516
517
518
  int bits;
  MACROBLOCK *x = &cpi->td.mb;
  AV1_COMMON *const cm = &cpi->common;
  const YV12_BUFFER_CONFIG *dgd = cm->frame_to_show;
519
  RestorationInfo *rsi = &cpi->rst_search;
520
521
522
523
524
525
526
527
528
529
  int tile_idx, tile_width, tile_height, nhtiles, nvtiles;
  int h_start, h_end, v_start, v_end;
  const int ntiles = av1_get_rest_ntiles(cm->width, cm->height, &tile_width,
                                         &tile_height, &nhtiles, &nvtiles);
  //  Make a copy of the unfiltered / processed recon buffer
  aom_yv12_copy_y(cm->frame_to_show, &cpi->last_frame_uf);
  av1_loop_filter_frame(cm->frame_to_show, cm, &cpi->td.mb.e_mbd, filter_level,
                        1, partial_frame);
  aom_yv12_copy_y(cm->frame_to_show, &cpi->last_frame_db);

530
531
  rsi->frame_restoration_type = RESTORE_DOMAINTXFMRF;
  rsi->domaintxfmrf_info = cpi->rst_search.domaintxfmrf_info;
532
533

  for (tile_idx = 0; tile_idx < ntiles; ++tile_idx)
534
    rsi->domaintxfmrf_info[tile_idx].level = 0;
535
536
537
538
539
  // Compute best Domaintxfm filters for each tile
  for (tile_idx = 0; tile_idx < ntiles; ++tile_idx) {
    av1_get_rest_tile_limits(tile_idx, 0, 0, nhtiles, nvtiles, tile_width,
                             tile_height, cm->width, cm->height, 0, 0, &h_start,
                             &h_end, &v_start, &v_end);
540
    err = sse_restoration_tile(src, cm->frame_to_show, cm, h_start,
541
                               h_end - h_start, v_start, v_end - v_start, 1);
542
543
544
545
546
547
548
549
550
551
552
553
554
555
    // #bits when a tile is not restored
    bits = av1_cost_bit(RESTORE_NONE_DOMAINTXFMRF_PROB, 0);
    cost_norestore = RDCOST_DBL(x->rdmult, x->rddiv, (bits >> 4), err);
    best_tile_cost[tile_idx] = DBL_MAX;

    search_domaintxfmrf_restoration(
        dgd->y_buffer + v_start * dgd->y_stride + h_start, h_end - h_start,
        v_end - v_start, dgd->y_stride,
        src->y_buffer + v_start * src->y_stride + h_start, src->y_stride,
#if CONFIG_AOM_HIGHBITDEPTH
        cm->bit_depth,
#else
        8,
#endif  // CONFIG_AOM_HIGHBITDEPTH
556
557
        &rsi->domaintxfmrf_info[tile_idx].sigma_r, cpi->extra_rstbuf,
        cm->rst_internal.tmpbuf);
558

559
560
    rsi->domaintxfmrf_info[tile_idx].level = 1;
    err = try_restoration_tile(src, cpi, rsi, 1, partial_frame, tile_idx, 0, 0,
561
                               dst_frame);
562
563
564
565
566
567
    bits = DOMAINTXFMRF_PARAMS_BITS << AV1_PROB_COST_SHIFT;
    bits += av1_cost_bit(RESTORE_NONE_DOMAINTXFMRF_PROB, 1);
    cost_domaintxfmrf = RDCOST_DBL(x->rdmult, x->rddiv, (bits >> 4), err);
    if (cost_domaintxfmrf >= cost_norestore) {
      domaintxfmrf_info[tile_idx].level = 0;
    } else {
568
      memcpy(&domaintxfmrf_info[tile_idx], &rsi->domaintxfmrf_info[tile_idx],
569
570
571
572
573
574
575
             sizeof(domaintxfmrf_info[tile_idx]));
      bits = DOMAINTXFMRF_PARAMS_BITS << AV1_PROB_COST_SHIFT;
      best_tile_cost[tile_idx] = RDCOST_DBL(
          x->rdmult, x->rddiv,
          (bits + cpi->switchable_restore_cost[RESTORE_DOMAINTXFMRF]) >> 4,
          err);
    }
576
    rsi->domaintxfmrf_info[tile_idx].level = 0;
577
578
  }
  // Cost for Domaintxfmrf filtering
579
  bits = frame_level_restore_bits[rsi->frame_restoration_type]
580
581
582
583
         << AV1_PROB_COST_SHIFT;
  for (tile_idx = 0; tile_idx < ntiles; ++tile_idx) {
    bits += av1_cost_bit(RESTORE_NONE_DOMAINTXFMRF_PROB,
                         domaintxfmrf_info[tile_idx].level);
584
    memcpy(&rsi->domaintxfmrf_info[tile_idx], &domaintxfmrf_info[tile_idx],
585
586
587
588
589
           sizeof(domaintxfmrf_info[tile_idx]));
    if (domaintxfmrf_info[tile_idx].level) {
      bits += (DOMAINTXFMRF_PARAMS_BITS << AV1_PROB_COST_SHIFT);
    }
  }
590
  err = try_restoration_frame(src, cpi, rsi, 1, partial_frame, dst_frame);
591
592
593
594
595
596
  cost_domaintxfmrf = RDCOST_DBL(x->rdmult, x->rddiv, (bits >> 4), err);

  aom_yv12_copy_y(&cpi->last_frame_uf, cm->frame_to_show);
  return cost_domaintxfmrf;
}

597
598
static double find_average(uint8_t *src, int h_start, int h_end, int v_start,
                           int v_end, int stride) {
599
600
601
  uint64_t sum = 0;
  double avg = 0;
  int i, j;
602
603
604
  for (i = v_start; i < v_end; i++)
    for (j = h_start; j < h_end; j++) sum += src[i * stride + j];
  avg = (double)sum / ((v_end - v_start) * (h_end - h_start));
605
606
607
  return avg;
}

608
609
610
static void compute_stats(uint8_t *dgd, uint8_t *src, int h_start, int h_end,
                          int v_start, int v_end, int dgd_stride,
                          int src_stride, double *M, double *H) {
611
  int i, j, k, l;
612
  double Y[WIENER_WIN2];
613
614
  const double avg =
      find_average(dgd, h_start, h_end, v_start, v_end, dgd_stride);
615

616
617
  memset(M, 0, sizeof(*M) * WIENER_WIN2);
  memset(H, 0, sizeof(*H) * WIENER_WIN2 * WIENER_WIN2);
618
619
  for (i = v_start; i < v_end; i++) {
    for (j = h_start; j < h_end; j++) {
620
621
      const double X = (double)src[i * src_stride + j] - avg;
      int idx = 0;
622
623
      for (k = -WIENER_HALFWIN; k <= WIENER_HALFWIN; k++) {
        for (l = -WIENER_HALFWIN; l <= WIENER_HALFWIN; l++) {
624
625
626
627
          Y[idx] = (double)dgd[(i + l) * dgd_stride + (j + k)] - avg;
          idx++;
        }
      }
628
      for (k = 0; k < WIENER_WIN2; ++k) {
629
        M[k] += Y[k] * X;
630
631
        H[k * WIENER_WIN2 + k] += Y[k] * Y[k];
        for (l = k + 1; l < WIENER_WIN2; ++l) {
632
          double value = Y[k] * Y[l];
633
634
          H[k * WIENER_WIN2 + l] += value;
          H[l * WIENER_WIN2 + k] += value;
635
636
637
638
639
640
        }
      }
    }
  }
}

Yaowu Xu's avatar
Yaowu Xu committed
641
#if CONFIG_AOM_HIGHBITDEPTH
642
643
static double find_average_highbd(uint16_t *src, int h_start, int h_end,
                                  int v_start, int v_end, int stride) {
644
645
646
  uint64_t sum = 0;
  double avg = 0;
  int i, j;
647
648
649
  for (i = v_start; i < v_end; i++)
    for (j = h_start; j < h_end; j++) sum += src[i * stride + j];
  avg = (double)sum / ((v_end - v_start) * (h_end - h_start));
650
651
652
  return avg;
}

653
654
655
656
static void compute_stats_highbd(uint8_t *dgd8, uint8_t *src8, int h_start,
                                 int h_end, int v_start, int v_end,
                                 int dgd_stride, int src_stride, double *M,
                                 double *H) {
657
  int i, j, k, l;
658
  double Y[WIENER_WIN2];
659
660
  uint16_t *src = CONVERT_TO_SHORTPTR(src8);
  uint16_t *dgd = CONVERT_TO_SHORTPTR(dgd8);
661
662
  const double avg =
      find_average_highbd(dgd, h_start, h_end, v_start, v_end, dgd_stride);
663

664
665
  memset(M, 0, sizeof(*M) * WIENER_WIN2);
  memset(H, 0, sizeof(*H) * WIENER_WIN2 * WIENER_WIN2);
666
667
  for (i = v_start; i < v_end; i++) {
    for (j = h_start; j < h_end; j++) {
668
669
      const double X = (double)src[i * src_stride + j] - avg;
      int idx = 0;
670
671
      for (k = -WIENER_HALFWIN; k <= WIENER_HALFWIN; k++) {
        for (l = -WIENER_HALFWIN; l <= WIENER_HALFWIN; l++) {
672
673
674
675
          Y[idx] = (double)dgd[(i + l) * dgd_stride + (j + k)] - avg;
          idx++;
        }
      }
676
      for (k = 0; k < WIENER_WIN2; ++k) {
677
        M[k] += Y[k] * X;
678
679
        H[k * WIENER_WIN2 + k] += Y[k] * Y[k];
        for (l = k + 1; l < WIENER_WIN2; ++l) {
680
          double value = Y[k] * Y[l];
681
682
          H[k * WIENER_WIN2 + l] += value;
          H[l * WIENER_WIN2 + k] += value;
683
684
685
686
687
        }
      }
    }
  }
}
Yaowu Xu's avatar
Yaowu Xu committed
688
#endif  // CONFIG_AOM_HIGHBITDEPTH
689
690
691
692
693
694
695
696
697
698
699
700
701
702
703
704
705
706
707
708
709
710

// Solves Ax = b, where x and b are column vectors
static int linsolve(int n, double *A, int stride, double *b, double *x) {
  int i, j, k;
  double c;
  // Partial pivoting
  for (i = n - 1; i > 0; i--) {
    if (A[(i - 1) * stride] < A[i * stride]) {
      for (j = 0; j < n; j++) {
        c = A[i * stride + j];
        A[i * stride + j] = A[(i - 1) * stride + j];
        A[(i - 1) * stride + j] = c;
      }
      c = b[i];
      b[i] = b[i - 1];
      b[i - 1] = c;
    }
  }
  // Forward elimination
  for (k = 0; k < n - 1; k++) {
    for (i = k; i < n - 1; i++) {
      c = A[(i + 1) * stride + k] / A[k * stride + k];
711
      for (j = 0; j < n; j++) A[(i + 1) * stride + j] -= c * A[k * stride + j];
712
713
714
715
716
      b[i + 1] -= c * b[k];
    }
  }
  // Backward substitution
  for (i = n - 1; i >= 0; i--) {
717
    if (fabs(A[i * stride + i]) < 1e-10) return 0;
718
    c = 0;
719
    for (j = i + 1; j <= n - 1; j++) c += A[i * stride + j] * x[j];
720
721
722
723
724
725
    x[i] = (b[i] - c) / A[i * stride + i];
  }
  return 1;
}

static INLINE int wrap_index(int i) {
726
  return (i >= WIENER_HALFWIN1 ? WIENER_WIN - 1 - i : i);
727
728
729
730
731
}

// Fix vector b, update vector a
static void update_a_sep_sym(double **Mc, double **Hc, double *a, double *b) {
  int i, j;
732
733
  double S[WIENER_WIN];
  double A[WIENER_WIN], B[WIENER_WIN2];
Aamir Anis's avatar
Aamir Anis committed
734
  int w, w2;
735
736
  memset(A, 0, sizeof(A));
  memset(B, 0, sizeof(B));
737
738
  for (i = 0; i < WIENER_WIN; i++) {
    for (j = 0; j < WIENER_WIN; ++j) {
739
740
741
742
      const int jj = wrap_index(j);
      A[jj] += Mc[i][j] * b[i];
    }
  }
743
744
  for (i = 0; i < WIENER_WIN; i++) {
    for (j = 0; j < WIENER_WIN; j++) {
745
      int k, l;
746
747
      for (k = 0; k < WIENER_WIN; ++k)
        for (l = 0; l < WIENER_WIN; ++l) {
748
749
          const int kk = wrap_index(k);
          const int ll = wrap_index(l);
750
751
          B[ll * WIENER_HALFWIN1 + kk] +=
              Hc[j * WIENER_WIN + i][k * WIENER_WIN2 + l] * b[i] * b[j];
752
753
754
        }
    }
  }
Aamir Anis's avatar
Aamir Anis committed
755
  // Normalization enforcement in the system of equations itself
756
  w = WIENER_WIN;
Aamir Anis's avatar
Aamir Anis committed
757
758
  w2 = (w >> 1) + 1;
  for (i = 0; i < w2 - 1; ++i)
759
760
    A[i] -=
        A[w2 - 1] * 2 + B[i * w2 + w2 - 1] - 2 * B[(w2 - 1) * w2 + (w2 - 1)];
Aamir Anis's avatar
Aamir Anis committed
761
762
763
764
765
766
767
768
769
  for (i = 0; i < w2 - 1; ++i)
    for (j = 0; j < w2 - 1; ++j)
      B[i * w2 + j] -= 2 * (B[i * w2 + (w2 - 1)] + B[(w2 - 1) * w2 + j] -
                            2 * B[(w2 - 1) * w2 + (w2 - 1)]);
  if (linsolve(w2 - 1, B, w2, A, S)) {
    S[w2 - 1] = 1.0;
    for (i = w2; i < w; ++i) {
      S[i] = S[w - 1 - i];
      S[w2 - 1] -= 2 * S[i];
770
    }
Aamir Anis's avatar
Aamir Anis committed
771
    memcpy(a, S, w * sizeof(*a));
772
773
774
775
776
777
  }
}

// Fix vector a, update vector b
static void update_b_sep_sym(double **Mc, double **Hc, double *a, double *b) {
  int i, j;
778
779
  double S[WIENER_WIN];
  double A[WIENER_WIN], B[WIENER_WIN2];
Aamir Anis's avatar
Aamir Anis committed
780
  int w, w2;
781
782
  memset(A, 0, sizeof(A));
  memset(B, 0, sizeof(B));
783
  for (i = 0; i < WIENER_WIN; i++) {
784
    const int ii = wrap_index(i);
785
    for (j = 0; j < WIENER_WIN; j++) A[ii] += Mc[i][j] * a[j];
786
787
  }

788
789
  for (i = 0; i < WIENER_WIN; i++) {
    for (j = 0; j < WIENER_WIN; j++) {
790
791
792
      const int ii = wrap_index(i);
      const int jj = wrap_index(j);
      int k, l;
793
794
795
796
      for (k = 0; k < WIENER_WIN; ++k)
        for (l = 0; l < WIENER_WIN; ++l)
          B[jj * WIENER_HALFWIN1 + ii] +=
              Hc[i * WIENER_WIN + j][k * WIENER_WIN2 + l] * a[k] * a[l];
797
798
    }
  }
Aamir Anis's avatar
Aamir Anis committed
799
  // Normalization enforcement in the system of equations itself
800
801
  w = WIENER_WIN;
  w2 = WIENER_HALFWIN1;
Aamir Anis's avatar
Aamir Anis committed
802
  for (i = 0; i < w2 - 1; ++i)
803
804
    A[i] -=
        A[w2 - 1] * 2 + B[i * w2 + w2 - 1] - 2 * B[(w2 - 1) * w2 + (w2 - 1)];
Aamir Anis's avatar
Aamir Anis committed
805
806
807
808
809
810
811
812
813
  for (i = 0; i < w2 - 1; ++i)
    for (j = 0; j < w2 - 1; ++j)
      B[i * w2 + j] -= 2 * (B[i * w2 + (w2 - 1)] + B[(w2 - 1) * w2 + j] -
                            2 * B[(w2 - 1) * w2 + (w2 - 1)]);
  if (linsolve(w2 - 1, B, w2, A, S)) {
    S[w2 - 1] = 1.0;
    for (i = w2; i < w; ++i) {
      S[i] = S[w - 1 - i];
      S[w2 - 1] -= 2 * S[i];
814
    }
Aamir Anis's avatar
Aamir Anis committed
815
    memcpy(b, S, w * sizeof(*b));
816
817
818
  }
}

819
820
static int wiener_decompose_sep_sym(double *M, double *H, double *a,
                                    double *b) {
821
  static const double init_filt[WIENER_WIN] = {
822
    0.035623, -0.127154, 0.211436, 0.760190, 0.211436, -0.127154, 0.035623,
823
824
  };
  int i, j, iter;
825
826
827
828
829
830
831
  double *Hc[WIENER_WIN2];
  double *Mc[WIENER_WIN];
  for (i = 0; i < WIENER_WIN; i++) {
    Mc[i] = M + i * WIENER_WIN;
    for (j = 0; j < WIENER_WIN; j++) {
      Hc[i * WIENER_WIN + j] =
          H + i * WIENER_WIN * WIENER_WIN2 + j * WIENER_WIN;
832
833
    }
  }
834
835
  memcpy(a, init_filt, sizeof(*a) * WIENER_WIN);
  memcpy(b, init_filt, sizeof(*b) * WIENER_WIN);
836
837
838
839
840
841
842

  iter = 1;
  while (iter < 10) {
    update_a_sep_sym(Mc, Hc, a, b);
    update_b_sep_sym(Mc, Hc, a, b);
    iter++;
  }
843
  return 1;
844
845
}

Aamir Anis's avatar
Aamir Anis committed
846
847
848
// Computes the function x'*A*x - x'*b for the learned filters, and compares
// against identity filters; Final score is defined as the difference between
// the function values
849
static double compute_score(double *M, double *H, int *vfilt, int *hfilt) {
850
  double ab[WIENER_WIN * WIENER_WIN];
Aamir Anis's avatar
Aamir Anis committed
851
852
853
854
855
  int i, k, l;
  double P = 0, Q = 0;
  double iP = 0, iQ = 0;
  double Score, iScore;
  int w;
856
857
858
859
860
861
862
863
  double a[WIENER_WIN], b[WIENER_WIN];
  w = WIENER_WIN;
  a[WIENER_HALFWIN] = b[WIENER_HALFWIN] = 1.0;
  for (i = 0; i < WIENER_HALFWIN; ++i) {
    a[i] = a[WIENER_WIN - i - 1] = (double)vfilt[i] / WIENER_FILT_STEP;
    b[i] = b[WIENER_WIN - i - 1] = (double)hfilt[i] / WIENER_FILT_STEP;
    a[WIENER_HALFWIN] -= 2 * a[i];
    b[WIENER_HALFWIN] -= 2 * b[i];
Aamir Anis's avatar
Aamir Anis committed
864
865
  }
  for (k = 0; k < w; ++k) {
866
    for (l = 0; l < w; ++l) ab[k * w + l] = a[l] * b[k];
Aamir Anis's avatar
Aamir Anis committed
867
868
869
  }
  for (k = 0; k < w * w; ++k) {
    P += ab[k] * M[k];
870
    for (l = 0; l < w * w; ++l) Q += ab[k] * H[k * w * w + l] * ab[l];
Aamir Anis's avatar
Aamir Anis committed
871
872
873
874
875
876
877
878
879
880
  }
  Score = Q - 2 * P;

  iP = M[(w * w) >> 1];
  iQ = H[((w * w) >> 1) * w * w + ((w * w) >> 1)];
  iScore = iQ - 2 * iP;

  return Score - iScore;
}

881
882
static void quantize_sym_filter(double *f, int *fi) {
  int i;
883
884
  for (i = 0; i < WIENER_HALFWIN; ++i) {
    fi[i] = RINT(f[i] * WIENER_FILT_STEP);
885
886
887
888
889
890
891
  }
  // Specialize for 7-tap filter
  fi[0] = CLIP(fi[0], WIENER_FILT_TAP0_MINV, WIENER_FILT_TAP0_MAXV);
  fi[1] = CLIP(fi[1], WIENER_FILT_TAP1_MINV, WIENER_FILT_TAP1_MAXV);
  fi[2] = CLIP(fi[2], WIENER_FILT_TAP2_MINV, WIENER_FILT_TAP2_MAXV);
}

892
893
static double search_wiener(const YV12_BUFFER_CONFIG *src, AV1_COMP *cpi,
                            int filter_level, int partial_frame,
894
895
                            RestorationInfo *info, double *best_tile_cost,
                            YV12_BUFFER_CONFIG *dst_frame) {
896
  WienerInfo *wiener_info = info->wiener_info;
Yaowu Xu's avatar
Yaowu Xu committed
897
  AV1_COMMON *const cm = &cpi->common;
898
  RestorationInfo *rsi = &cpi->rst_search;
899
900
  int64_t err;
  int bits;
901
  double cost_wiener, cost_norestore;
902
  MACROBLOCK *x = &cpi->td.mb;
903
904
905
  double M[WIENER_WIN2];
  double H[WIENER_WIN2 * WIENER_WIN2];
  double vfilterd[WIENER_WIN], hfilterd[WIENER_WIN];
906
907
908
909
910
  const YV12_BUFFER_CONFIG *dgd = cm->frame_to_show;
  const int width = cm->width;
  const int height = cm->height;
  const int src_stride = src->y_stride;
  const int dgd_stride = dgd->y_stride;
Aamir Anis's avatar
Aamir Anis committed
911
  double score;
912
  int tile_idx, tile_width, tile_height, nhtiles, nvtiles;
913
  int h_start, h_end, v_start, v_end;
914
  int i;
915

916
917
  const int ntiles = av1_get_rest_ntiles(width, height, &tile_width,
                                         &tile_height, &nhtiles, &nvtiles);
918
919
920
921
922
923
  assert(width == dgd->y_crop_width);
  assert(height == dgd->y_crop_height);
  assert(width == src->y_crop_width);
  assert(height == src->y_crop_height);

  //  Make a copy of the unfiltered / processed recon buffer
Yaowu Xu's avatar
Yaowu Xu committed
924
925
926
927
  aom_yv12_copy_y(cm->frame_to_show, &cpi->last_frame_uf);
  av1_loop_filter_frame(cm->frame_to_show, cm, &cpi->td.mb.e_mbd, filter_level,
                        1, partial_frame);
  aom_yv12_copy_y(cm->frame_to_show, &cpi->last_frame_db);
928

929
  rsi->frame_restoration_type = RESTORE_WIENER;
930

931
  for (tile_idx = 0; tile_idx < ntiles; ++tile_idx)
932
    rsi->wiener_info[tile_idx].level = 0;
933
934
935

  // Compute best Wiener filters for each tile
  for (tile_idx = 0; tile_idx < ntiles; ++tile_idx) {
936
937
938
    av1_get_rest_tile_limits(tile_idx, 0, 0, nhtiles, nvtiles, tile_width,
                             tile_height, width, height, 0, 0, &h_start, &h_end,
                             &v_start, &v_end);
939
    err = sse_restoration_tile(src, cm->frame_to_show, cm, h_start,
940
                               h_end - h_start, v_start, v_end - v_start, 1);
941
942
    // #bits when a tile is not restored
    bits = av1_cost_bit(RESTORE_NONE_WIENER_PROB, 0);
943
    cost_norestore = RDCOST_DBL(x->rdmult, x->rddiv, (bits >> 4), err);
944
    best_tile_cost[tile_idx] = DBL_MAX;
945
946

    av1_get_rest_tile_limits(tile_idx, 0, 0, nhtiles, nvtiles, tile_width,
947
948
949
                             tile_height, width, height, WIENER_HALFWIN,
                             WIENER_HALFWIN, &h_start, &h_end, &v_start,
                             &v_end);
Yaowu Xu's avatar
Yaowu Xu committed
950
#if CONFIG_AOM_HIGHBITDEPTH
951
952
953
954
    if (cm->use_highbitdepth)
      compute_stats_highbd(dgd->y_buffer, src->y_buffer, h_start, h_end,
                           v_start, v_end, dgd_stride, src_stride, M, H);
    else
Yaowu Xu's avatar
Yaowu Xu committed
955
#endif  // CONFIG_AOM_HIGHBITDEPTH
956
957
958
      compute_stats(dgd->y_buffer, src->y_buffer, h_start, h_end, v_start,
                    v_end, dgd_stride, src_stride, M, H);

959
    wiener_info[tile_idx].level = 1;
960
    if (!wiener_decompose_sep_sym(M, H, vfilterd, hfilterd)) {
961
      wiener_info[tile_idx].level = 0;
962
963
      continue;
    }
964
965
    quantize_sym_filter(vfilterd, rsi->wiener_info[tile_idx].vfilter);
    quantize_sym_filter(hfilterd, rsi->wiener_info[tile_idx].hfilter);
966
967
968
969

    // Filter score computes the value of the function x'*A*x - x'*b for the
    // learned filter and compares it against identity filer. If there is no
    // reduction in the function, the filter is reverted back to identity
970
971
    score = compute_score(M, H, rsi->wiener_info[tile_idx].vfilter,
                          rsi->wiener_info[tile_idx].hfilter);
972
    if (score > 0.0) {
973
      wiener_info[tile_idx].level = 0;
974
975
      continue;
    }
976

977
978
    rsi->wiener_info[tile_idx].level = 1;
    err = try_restoration_tile(src, cpi, rsi, 1, partial_frame, tile_idx, 0, 0,
979
                               dst_frame);
980
981
982
    bits = WIENER_FILT_BITS << AV1_PROB_COST_SHIFT;
    bits += av1_cost_bit(RESTORE_NONE_WIENER_PROB, 1);
    cost_wiener = RDCOST_DBL(x->rdmult, x->rddiv, (bits >> 4), err);
983
    if (cost_wiener >= cost_norestore) {
984
985
986
      wiener_info[tile_idx].level = 0;
    } else {
      wiener_info[tile_idx].level = 1;
987
988
989
990
991
      for (i = 0; i < WIENER_HALFWIN; ++i) {
        wiener_info[tile_idx].vfilter[i] =
            rsi->wiener_info[tile_idx].vfilter[i];
        wiener_info[tile_idx].hfilter[i] =
            rsi->wiener_info[tile_idx].hfilter[i];
992
      }
993
994
995
996
997
      bits = WIENER_FILT_BITS << AV1_PROB_COST_SHIFT;
      best_tile_cost[tile_idx] = RDCOST_DBL(
          x->rdmult, x->rddiv,
          (bits + cpi->switchable_restore_cost[RESTORE_WIENER]) >> 4, err);
    }
998
    rsi->wiener_info[tile_idx].level = 0;
999
  }
1000
  // Cost for Wiener filtering
1001
  bits = frame_level_restore_bits[rsi->frame_restoration_type]
1002
         << AV1_PROB_COST_SHIFT;
1003
  for (tile_idx = 0; tile_idx < ntiles; ++tile_idx) {
1004
    bits += av1_cost_bit(RESTORE_NONE_WIENER_PROB, wiener_info[tile_idx].level);
1005
    rsi->wiener_info[tile_idx].level = wiener_info[tile_idx].level;
1006
    if (wiener_info[tile_idx].level) {
1007
      bits += (WIENER_FILT_BITS << AV1_PROB_COST_SHIFT);
1008
1009
1010
1011
1012
      for (i = 0; i < WIENER_HALFWIN; ++i) {
        rsi->wiener_info[tile_idx].vfilter[i] =
            wiener_info[tile_idx].vfilter[i];
        rsi->wiener_info[tile_idx].hfilter[i] =
            wiener_info[tile_idx].hfilter[i];
1013
1014
      }
    }
Aamir Anis's avatar
Aamir Anis committed
1015
  }
1016
  err = try_restoration_frame(src, cpi, rsi, 1, partial_frame, dst_frame);
1017
  cost_wiener = RDCOST_DBL(x->rdmult, x->rddiv, (bits >> 4), err);
1018

Yaowu Xu's avatar
Yaowu Xu committed
1019
  aom_yv12_copy_y(&cpi->last_frame_uf, cm->frame_to_show);
1020
  return cost_wiener;
1021
1022
}

1023
1024
static double search_norestore(const YV12_BUFFER_CONFIG *src, AV1_COMP *cpi,
                               int filter_level, int partial_frame,
1025
1026
                               RestorationInfo *info, double *best_tile_cost,
                               YV12_BUFFER_CONFIG *dst_frame) {
1027
1028
  double err, cost_norestore;
  int bits;
1029
  MACROBLOCK *x = &cpi->td.mb;
1030
1031
  AV1_COMMON *const cm = &cpi->common;
  int tile_idx, tile_width, tile_height, nhtiles, nvtiles;
1032
1033
1034
  int h_start, h_end, v_start, v_end;
  const int ntiles = av1_get_rest_ntiles(cm->width, cm->height, &tile_width,
                                         &tile_height, &nhtiles, &nvtiles);
1035
  (void)info;
1036
  (void)dst_frame;
1037
1038
1039
1040
1041
1042
1043
1044
1045
1046
1047

  //  Make a copy of the unfiltered / processed recon buffer
  aom_yv12_copy_y(cm->frame_to_show, &cpi->last_frame_uf);
  av1_loop_filter_frame(cm->frame_to_show, cm, &cpi->td.mb.e_mbd, filter_level,
                        1, partial_frame);
  aom_yv12_copy_y(cm->frame_to_show, &cpi->last_frame_db);

  for (tile_idx = 0; tile_idx < ntiles; ++tile_idx) {
    av1_get_rest_tile_limits(tile_idx, 0, 0, nhtiles, nvtiles, tile_width,
                             tile_height, cm->width, cm->height, 0, 0, &h_start,
                             &h_end, &v_start, &v_end);
1048
    err = sse_restoration_tile(src, cm->frame_to_show, cm, h_start,
1049
                               h_end - h_start, v_start, v_end - v_start, 1);
1050
    best_tile_cost[tile_idx] =