IDSDDownsampling.cpp 14.4 KB
Newer Older
sarthou's avatar
sarthou committed
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
//
// Created by Jonathan on 5/17/2016.
//

#include "IDSDDownsampling.h"

IDSDDownsampling::IDSDDownsampling(InpaintingDataSource * const source, int const downsamplingLevels, int const downsamplingRate)
: IDSDecorator(source), idsLevels(new InpaintingDataSource*[downsamplingLevels]),
  downsamplingLevels(downsamplingLevels), downsamplingRate(downsamplingRate), curLevel(downsamplingLevels - 1)
{
    idsLevels[0] = decoratedIDS; // The base downsampling is the same as the source.

    PaddedMatrixDim * const oriDim = idsLevels[0]->GetDim();
    int x = oriDim->GetSizeX();
    int y = oriDim->GetSizeY();
    int z = oriDim->GetSizeZ(); // No Z downsampling. Could be added.
    for(int i = 1; i < downsamplingLevels; ++i){
        x /= downsamplingRate;
        y /= downsamplingRate;
        // z /= downsamplingRate;

        PaddedMatrixDim * const nextDim = new PaddedMatrixDim(x, y, z, oriDim->GetBuffX(), oriDim->GetBuffY(), oriDim->GetBuffZ());
        idsLevels[i] = new IDSConcrete(nextDim); // All the matrices are empty.

        DownsampleColorsAndOcclusion(i); // Fill the RGB and OCC channels of the "i" IDS.
        DownsampleCasellesDescriptors(i);
    }
}

IDSDDownsampling::~IDSDDownsampling(){
    for(int i = 1; i < downsamplingLevels; ++i)
        delete idsLevels[i]; // Index 0 deletes the ids contained in the base class. (Which is bad, since it doesn't own it.)

    delete[] idsLevels;
}

void IDSDDownsampling::SetCurLevel(int const curLevel){
    this->curLevel = curLevel;
}

void IDSDDownsampling::UpsampleShiftMatrices(){
    PaddedMatrixDim * fromDim = idsLevels[curLevel]->GetDim();
    int* fromXDispData = idsLevels[curLevel]->GetDispX()->GetData();
    int* fromYDispData = idsLevels[curLevel]->GetDispY()->GetData();
    int* fromZDispData = idsLevels[curLevel]->GetDispZ()->GetData();
    float* fromADispData = idsLevels[curLevel]->GetDispA()->GetData();

    PaddedMatrixDim * toDim = idsLevels[curLevel - 1]->GetDim();
    int* toXDispData = idsLevels[curLevel - 1]->GetDispX()->GetData();
    int* toYDispData = idsLevels[curLevel - 1]->GetDispY()->GetData();
    int* toZDispData = idsLevels[curLevel - 1]->GetDispZ()->GetData();
    float* toADispData = idsLevels[curLevel - 1]->GetDispA()->GetData();

    PaddedMatrixIt nextVolIt = PaddedMatrixIt(toDim);
    PaddedMatrixIt sourceVolIt = PaddedMatrixIt(fromDim);
    do{
        int curPos = nextVolIt.GetPosition(); // The same position
        int nexPos = nextVolIt.GetPosition() + 1; // The pixel to the right of the current.
        int lowPos = nextVolIt.GetPosition() + toDim->GetBuffSizeX(); // The pixel below the current.
        int diaPos = nextVolIt.GetPosition() + toDim->GetBuffSizeX() + 1; // The pixel below and to to right - in diagonal.

        int nextXDisp = downsamplingRate * fromXDispData[sourceVolIt.GetPosition()]; // There are more X pixels in the next level, so we multiply the current value.
        toXDispData[curPos] = nextXDisp;
        toXDispData[nexPos] = nextXDisp;
        toXDispData[lowPos] = nextXDisp;
        toXDispData[diaPos] = nextXDisp;

        int nextYDisp = downsamplingRate * fromYDispData[sourceVolIt.GetPosition()]; // There are more Y pixels in the next level, so we multiply the current value.
        toYDispData[curPos] = nextYDisp;
        toYDispData[nexPos] = nextYDisp;
        toYDispData[lowPos] = nextYDisp;
        toYDispData[diaPos] = nextYDisp;

        int nextZDisp = fromZDispData[sourceVolIt.GetPosition()]; // There is no Z downsampling... So we keep the same value.
        toZDispData[curPos] = nextZDisp;
        toZDispData[nexPos] = nextZDisp;
        toZDispData[lowPos] = nextZDisp;
        toZDispData[diaPos] = nextZDisp;

        float nextADisp = fromADispData[sourceVolIt.GetPosition()]; // We can't evaluate the true value of the displacement until it gets reconstructed. (0 means uninitialized!)
        toADispData[curPos] = nextADisp;
        toADispData[nexPos] = nextADisp;
        toADispData[lowPos] = nextADisp;
        toADispData[diaPos] = nextADisp;
    }while(nextVolIt.Next(downsamplingRate, downsamplingRate) && sourceVolIt.Next());
}

88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
void IDSDDownsampling::MatlabUpsampleShiftMatrix(){
  /*This method is for interpolate the displacements fields, with an algorithm
  more similar than what is done in the Matlab example, to compare efficiency */
  PaddedMatrixDim * fromDim = idsLevels[curLevel]->GetDim();
  int* fromXDispData = idsLevels[curLevel]->GetDispX()->GetData();
  int* fromYDispData = idsLevels[curLevel]->GetDispY()->GetData();
  int* fromZDispData = idsLevels[curLevel]->GetDispZ()->GetData();
  float* fromADispData = idsLevels[curLevel]->GetDispA()->GetData();

  PaddedMatrixDim * toDim = idsLevels[curLevel - 1]->GetDim();
  int* toXDispData = idsLevels[curLevel - 1]->GetDispX()->GetData();
  int* toYDispData = idsLevels[curLevel - 1]->GetDispY()->GetData();
  int* toZDispData = idsLevels[curLevel - 1]->GetDispZ()->GetData();
  float* toADispData = idsLevels[curLevel - 1]->GetDispA()->GetData();

  PaddedMatrix<int> * GridX = new PaddedMatrix<int>(toDim);
  PaddedMatrix<int> * GridY = new PaddedMatrix<int>(toDim);
  PaddedMatrix<int> * GridZ = new PaddedMatrix<int>(toDim);

107
108
109
110
  PaddedMatrix<int> * CoordX = new PaddedMatrix<int>(toDim);
  PaddedMatrix<int> * CoordY = new PaddedMatrix<int>(toDim);
  PaddedMatrix<int> * CoordZ = new PaddedMatrix<int>(toDim);

111
112
113
114
115
116
117
118
119
120
121
122
123
  // Creating Meshgrid
  for(int i = 0; i < toDim->GetSizeX(); i++){
    for(int j = 0; j < toDim->GetSizeY(); j++){
      for(int k = 0; k < toDim->GetSizeZ(); k++){
        int temp_x = std::max(std::min((int)roundf(i/2),fromDim->GetSizeX()-1),0);
        int temp_y = std::max(std::min((int)roundf(j/2),fromDim->GetSizeY()-1),0);
        int temp_z = std::max(std::min(k,fromDim->GetSizeZ()-1),0);
        GridX->Set(i,j,k,temp_x);
        GridY->Set(i,j,k,temp_y);
        GridZ->Set(i,j,k,temp_z);
      }
    }
  }
124

125
126
127
128
129
130
  PaddedMatrixIt nextVolIt = PaddedMatrixIt(toDim);
  do {
    toXDispData[nextVolIt.GetPosition()] = downsamplingRate * \
     fromXDispData[fromDim->DeterminePosition(GridX->GetData()[nextVolIt.GetPosition()],\
     GridY->GetData()[nextVolIt.GetPosition()],GridZ->GetData()[nextVolIt.GetPosition()])];

131
132
133
134
    // CoordX->GetData()[nextVolIt.GetPosition()] = GridX->GetData()[nextVolIt.GetPosition()];
    // CoordY->GetData()[nextVolIt.GetPosition()] = GridY->GetData()[nextVolIt.GetPosition()];
    // CoordZ->GetData()[nextVolIt.GetPosition()] = GridZ->GetData()[nextVolIt.GetPosition()];

135
136
137
138
    toYDispData[nextVolIt.GetPosition()] = downsamplingRate * \
     fromYDispData[fromDim->DeterminePosition(GridX->GetData()[nextVolIt.GetPosition()],\
     GridY->GetData()[nextVolIt.GetPosition()],GridZ->GetData()[nextVolIt.GetPosition()])];

139
140
141
142
143
144
145
146

    // On garde les mêmes valeurs car il n'y a pas de downsampling sur Z et on "clamp" les valeurs
    int tempDispZ = fromZDispData[fromDim->DeterminePosition(GridX->GetData()[nextVolIt.GetPosition()],\
    GridY->GetData()[nextVolIt.GetPosition()],GridZ->GetData()[nextVolIt.GetPosition()])];

    // TODO Here the 2 (half patch size is hardcoded, have to change it later)
    int minBarrier = 2;
    int maxBarrier = toDim->GetSizeZ() - 2;
147
    toZDispData[nextVolIt.GetPosition()] =\
148
149
150
     std::max(std::min(tempDispZ,maxBarrier-GridZ->GetData()[nextVolIt.GetPosition()]),minBarrier-GridZ->GetData()[nextVolIt.GetPosition()]);
    // toZDispData[nextVolIt.GetPosition()] = tempDispZ;

151
152
153
154
155
156
157

    toADispData[nextVolIt.GetPosition()] =\
     fromADispData[fromDim->DeterminePosition(GridX->GetData()[nextVolIt.GetPosition()],\
     GridY->GetData()[nextVolIt.GetPosition()],GridZ->GetData()[nextVolIt.GetPosition()])];

  } while(nextVolIt.Next());

158
159
160
161
162
163
164

  // for(int nbFrame=0;nbFrame < toDim->GetSizeZ();nbFrame++)
  // {
  //   CoordX->GetFrame(nbFrame)->MprintToCSV("Data/upsample_ant/coord/Cpp_UA_CoordX_F"+std::to_string(nbFrame));
  //   CoordY->GetFrame(nbFrame)->MprintToCSV("Data/upsample_ant/coord/Cpp_UA_CoordY_F"+std::to_string(nbFrame));
  //   CoordZ->GetFrame(nbFrame)->MprintToCSV("Data/upsample_ant/coord/Cpp_UA_CoordZ_F"+std::to_string(nbFrame));
  // }
165
166
}

sarthou's avatar
sarthou committed
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
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
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
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
PaddedMatrixDim * const IDSDDownsampling::GetDim(){
    return idsLevels[curLevel]->GetDim();
}

PaddedMatrix<float> * const IDSDDownsampling::GetRed(){
    return idsLevels[curLevel]->GetRed();
}

PaddedMatrix<float> * const IDSDDownsampling::GetBlue(){
    return idsLevels[curLevel]->GetBlue();
}

PaddedMatrix<float> * const IDSDDownsampling::GetGreen(){
    return idsLevels[curLevel]->GetGreen();
}

PaddedMatrix<float> * const IDSDDownsampling::GetOcc(){
    return idsLevels[curLevel]->GetOcc();
}

PaddedMatrix<float> * const IDSDDownsampling::GetGradX(){
    return idsLevels[curLevel]->GetGradX();
}

PaddedMatrix<float> * const IDSDDownsampling::GetGradY(){
    return idsLevels[curLevel]->GetGradY();
}

PaddedMatrix<int> * const IDSDDownsampling::GetDispX(){
    return idsLevels[curLevel]->GetDispX();
}

PaddedMatrix<int> * const IDSDDownsampling::GetDispY(){
    return idsLevels[curLevel]->GetDispY();
}

PaddedMatrix<int> * const IDSDDownsampling::GetDispZ(){
    return idsLevels[curLevel]->GetDispZ();
}

PaddedMatrix<float> * const IDSDDownsampling::GetDispA(){
    return idsLevels[curLevel]->GetDispA();
}

// TODO : This function assumes a minimum of "downsamplingRate" padding on the X and Y axis.
void IDSDDownsampling::DownsampleColorsAndOcclusion(int const targetLevel){
    float gaussianFilter[] = {0.307801329123470f, 0.384397341753060f, 0.307801329123470f}; // ImageManipulator::GaussianFilter(0.5, downsamplingRate + 1, 1); // TODO : Soft code the 0.5

    PaddedMatrixDim xFilterDim = PaddedMatrixDim(3,1,1);
    PaddedMatrix<float> gaussFilterX = PaddedMatrix<float>(&xFilterDim);

    PaddedMatrixDim yFilterDim = PaddedMatrixDim(1,3,1);
    PaddedMatrix<float> gaussFilterY = PaddedMatrix<float>(&yFilterDim);

    for(int i = 0; i < 3; i++){
        gaussFilterX.GetData()[i] = gaussianFilter[i];
        gaussFilterY.GetData()[i] = gaussianFilter[i];
    }

    InpaintingDataSource * const sourceIds = idsLevels[targetLevel - 1];
    InpaintingDataSource * const targetIds = idsLevels[targetLevel];

    PaddedMatrix<float> xFilteredRed = ImageManipulator::ImageFilter(sourceIds->GetRed(), &gaussFilterX);
    PaddedMatrix<float> xFilteredGreen = ImageManipulator::ImageFilter(sourceIds->GetGreen(), &gaussFilterX);
    PaddedMatrix<float> xFilteredBlue = ImageManipulator::ImageFilter(sourceIds->GetBlue(), &gaussFilterX);
    PaddedMatrix<float> xFilteredOcc = ImageManipulator::ImageFilter(sourceIds->GetOcc(), &gaussFilterX);

    PaddedMatrix<float> fullFilteredRed = ImageManipulator::ImageFilter(&xFilteredRed, &gaussFilterY);
    PaddedMatrix<float> fullFilteredGreen = ImageManipulator::ImageFilter(&xFilteredGreen, &gaussFilterY);
    PaddedMatrix<float> fullFilteredBlue = ImageManipulator::ImageFilter(&xFilteredBlue, &gaussFilterY);
    PaddedMatrix<float> fullFilteredOcc = ImageManipulator::ImageFilter(&xFilteredOcc, &gaussFilterY);

    PaddedMatrixIt sourceLevelIt = PaddedMatrixIt(sourceIds->GetDim());
    PaddedMatrixIt targetLevelIt = PaddedMatrixIt(targetIds->GetDim());
    do{
        targetIds->GetRed()->GetData()[targetLevelIt.GetPosition()] = fullFilteredRed.GetData()[sourceLevelIt.GetPosition()];
        targetIds->GetGreen()->GetData()[targetLevelIt.GetPosition()] = fullFilteredGreen.GetData()[sourceLevelIt.GetPosition()];
        targetIds->GetBlue()->GetData()[targetLevelIt.GetPosition()] = fullFilteredBlue.GetData()[sourceLevelIt.GetPosition()];
        targetIds->GetOcc()->GetData()[targetLevelIt.GetPosition()] = ceilf(fullFilteredOcc.GetData()[sourceLevelIt.GetPosition()]);
    }while(targetLevelIt.Next() && sourceLevelIt.Next(downsamplingRate, downsamplingRate));
}

void IDSDDownsampling::DownsampleCasellesDescriptors(int const targetLevel){
    int leap = (int)powf(downsamplingRate, targetLevel); // Could overflow easy.

    InpaintingDataSource * const sourceIds = idsLevels[0];
    InpaintingDataSource * const targetIds = idsLevels[targetLevel];

    PaddedMatrixIt sourceIt = PaddedMatrixIt(sourceIds->GetDim());
    PaddedMatrixIt targetIt = PaddedMatrixIt(targetIds->GetDim());
    do{
        targetIds->GetGradX()->GetData()[targetIt.GetPosition()] = sourceIds->GetGradX()->GetData()[sourceIt.GetPosition()];
        targetIds->GetGradY()->GetData()[targetIt.GetPosition()] = sourceIds->GetGradY()->GetData()[sourceIt.GetPosition()];
    }while(sourceIt.Next(leap, leap) && targetIt.Next());
}

263

sarthou's avatar
sarthou committed
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
/*
    PaddedMatrixDim gaussianKernelDim = PaddedMatrixDim(downsamplingRate + 1, downsamplingRate + 1, 1);

    InpaintingDataSource * const sourceIds = idsLevels[targetLevel - 1];
    InpaintingDataSource * const targetIds = idsLevels[targetLevel];

    // Need approximation for out-of-bound dataSources, so we mirror the edge pixels in the padding.
    sourceIds->GetRed()->MirrorIntoBuffer();
    sourceIds->GetGreen()->MirrorIntoBuffer();
    sourceIds->GetBlue()->MirrorIntoBuffer();
    sourceIds->GetOcc()->MirrorIntoBuffer();

    float* sourceRed = sourceIds->GetRed()->GetData();
    float* sourceBlu = sourceIds->GetBlue()->GetData();
    float* sourceGrn = sourceIds->GetGreen()->GetData();
    float* sourceOcc = sourceIds->GetOcc()->GetData();

    float* dwnSmpRed = targetIds->GetRed()->GetData();
    float* dwnSmpBlu = targetIds->GetBlue()->GetData();
    float* dwnSmpGrn = targetIds->GetGreen()->GetData();
    float* dwnSmpOcc = targetIds->GetOcc()->GetData();

    PaddedMatrixIt sourceLevelIt = PaddedMatrixIt(sourceIds->GetDim());
    PaddedMatrixIt targetLevelIt = PaddedMatrixIt(targetIds->GetDim());
    PaddedMatrixKernelIt gaussianKernelIt = PaddedMatrixKernelIt(sourceIds->GetDim(), &gaussianKernelDim);
    do{
        float redVal = 0.0f, bluVal = 0.0f, grnVal = 0.0f, occVal = 0.0f;
        gaussianKernelIt.Start(sourceLevelIt.GetPosition());
        do{
            redVal += sourceRed[gaussianKernelIt.GetPosition()] * gaussianFilter[gaussianKernelIt.GetIndex()];
            bluVal += sourceBlu[gaussianKernelIt.GetPosition()] * gaussianFilter[gaussianKernelIt.GetIndex()];
            grnVal += sourceGrn[gaussianKernelIt.GetPosition()] * gaussianFilter[gaussianKernelIt.GetIndex()];
            occVal += sourceOcc[gaussianKernelIt.GetPosition()] * gaussianFilter[gaussianKernelIt.GetIndex()];
        }while(gaussianKernelIt.Next());

        // If always occlude more.
        occVal = ceilf(occVal);

        dwnSmpRed[targetLevelIt.GetPosition()] = redVal;
        dwnSmpBlu[targetLevelIt.GetPosition()] = bluVal;
        dwnSmpGrn[targetLevelIt.GetPosition()] = grnVal;
        dwnSmpOcc[targetLevelIt.GetPosition()] = occVal;
    }while(targetLevelIt.Next() && sourceLevelIt.Next(downsamplingRate, downsamplingRate));
307
    */