1 | % * This code was used in the following articles:
|
---|
2 | % * [1] Learning 3-D Scene Structure from a Single Still Image,
|
---|
3 | % * Ashutosh Saxena, Min Sun, Andrew Y. Ng,
|
---|
4 | % * In ICCV workshop on 3D Representation for Recognition (3dRR-07), 2007.
|
---|
5 | % * (best paper)
|
---|
6 | % * [2] 3-D Reconstruction from Sparse Views using Monocular Vision,
|
---|
7 | % * Ashutosh Saxena, Min Sun, Andrew Y. Ng,
|
---|
8 | % * In ICCV workshop on Virtual Representations and Modeling
|
---|
9 | % * of Large-scale environments (VRML), 2007.
|
---|
10 | % * [3] 3-D Depth Reconstruction from a Single Still Image,
|
---|
11 | % * Ashutosh Saxena, Sung H. Chung, Andrew Y. Ng.
|
---|
12 | % * International Journal of Computer Vision (IJCV), Aug 2007.
|
---|
13 | % * [6] Learning Depth from Single Monocular Images,
|
---|
14 | % * Ashutosh Saxena, Sung H. Chung, Andrew Y. Ng.
|
---|
15 | % * In Neural Information Processing Systems (NIPS) 18, 2005.
|
---|
16 | % *
|
---|
17 | % * These articles are available at:
|
---|
18 | % * http://make3d.stanford.edu/publications
|
---|
19 | % *
|
---|
20 | % * We request that you cite the papers [1], [3] and [6] in any of
|
---|
21 | % * your reports that uses this code.
|
---|
22 | % * Further, if you use the code in image3dstiching/ (multiple image version),
|
---|
23 | % * then please cite [2].
|
---|
24 | % *
|
---|
25 | % * If you use the code in third_party/, then PLEASE CITE and follow the
|
---|
26 | % * LICENSE OF THE CORRESPONDING THIRD PARTY CODE.
|
---|
27 | % *
|
---|
28 | % * Finally, this code is for non-commercial use only. For further
|
---|
29 | % * information and to obtain a copy of the license, see
|
---|
30 | % *
|
---|
31 | % * http://make3d.stanford.edu/publications/code
|
---|
32 | % *
|
---|
33 | % * Also, the software distributed under the License is distributed on an
|
---|
34 | % * "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either
|
---|
35 | % * express or implied. See the License for the specific language governing
|
---|
36 | % * permissions and limitations under the License.
|
---|
37 | % *
|
---|
38 | % */
|
---|
39 | function ReportPlaneParaMRF_FreeSup(step, DepthFolder,... |
---|
40 | Sup,MedSup,depthMap,VarMap,RayOri, Ray,MedRay,maskSky,maskG,Algo,k, CornerList, OccluList,... |
---|
41 | MultiScaleSupTable, StraightLineTable, HBrokeBook, VBrokeBook,previoslyStored,... |
---|
42 | baseline); |
---|
43 | % This function runs the RMF over the plane parameter of each superpixels |
---|
44 | |
---|
45 | global GeneralDataFolder ScratchDataFolder LocalFolder ClusterExecutionDirectory... |
---|
46 | ImgFolder VertYNuPatch VertYNuDepth HoriXNuPatch HoriXNuDepth a_default b_default Ox_default Oy_default... |
---|
47 | Horizon_default filename batchSize NuRow_default SegVertYSize SegHoriXSize WeiBatchSize PopUpVertY PopUpHoriX taskName; |
---|
48 | |
---|
49 | if nargin <20 |
---|
50 | baseline = 0; |
---|
51 | end |
---|
52 | |
---|
53 | % initialize parameters |
---|
54 | displayFlag = false; |
---|
55 | RenderVrmlFlag = true; |
---|
56 | StickHori = 5%0.1; % sticking power in horizontal direction |
---|
57 | StickVert = 10; % sticking power in vertical direction |
---|
58 | Center = 5; % Co-Planar weight at the Center of each superpixel |
---|
59 | HoriConf = 1; % set the confidant of the learned depth at the middle in Horizontal direction of the image |
---|
60 | VertConf = 1; % set the confidant of the learned depth at the top of the image |
---|
61 | mapVert = linspace(VertConf,1,VertYNuDepth); % modeling the gravity prior |
---|
62 | mapHori = [linspace(HoriConf,1,round(HoriXNuDepth/2)) fliplr(linspace(HoriConf,1,HoriXNuDepth-round(HoriXNuDepth/2)))]; % modeling the gravity prior |
---|
63 | % ========set the range of depth that our model in |
---|
64 | ClosestDist = 1; |
---|
65 | % set the FarestDist to very 5 times to median depth |
---|
66 | FarestDist = 5*median(depthMap(:)) |
---|
67 | % ================================================ |
---|
68 | ceiling = 0*VertYNuDepth; % set the position of the ceiling, related to No plane coming back constrain |
---|
69 | Name{1} = 'FraWOPri'; |
---|
70 | Name{2} = 'FraCoP'; |
---|
71 | %pause |
---|
72 | if isempty(MultiScaleSupTable) |
---|
73 | Name{3} = 'Var_FraStickCoP'; |
---|
74 | else |
---|
75 | Name{3} = 'Var_FraStickCoPSTasCoP'; |
---|
76 | end |
---|
77 | if ~isempty(MultiScaleSupTable) |
---|
78 | MultiScaleFlag = true; |
---|
79 | WeiV = 2*ones(1,size(MultiScaleSupTable,2)-1); |
---|
80 | else |
---|
81 | MultiScaleFlag = false; |
---|
82 | WeiV = 1; |
---|
83 | end |
---|
84 | WeiV(1,1:2:end) = 6; % emphasize the middle scale three times smaller than large scale |
---|
85 | WeiV =WeiV./sum(WeiV);% normalize if pair of superpixels have same index in all the scale, their weight will be 10 |
---|
86 | ShiftStick = -.1; % between -1 and 0, more means more smoothing. |
---|
87 | ShiftCoP = -.5; % between -1 and 0, more means more smoothing. |
---|
88 | gravity =true; % if true, apply the HoriConf and VertConf linear scale weight |
---|
89 | CoPST = true; % if true, apply the Straight line prior as the Co-Planar constrain |
---|
90 | ConerImprove = true; |
---|
91 | FractionalDepthError = true; |
---|
92 | % get rid of the error point and sky point in the depthMap |
---|
93 | % set every depth bigger than FarestDistmeter to FarestDistmeters |
---|
94 | %CleanedDepthMap = (depthMapif ~previoslyStored >80).*medfilt2(depthMap,[4 4])+(depthMap<=80).*depthMap; |
---|
95 | CleanedDepthMap = depthMap; |
---|
96 | %CleanedDepthMap(depthMap>FarestDist) = FarestDist; % don't clean the point >80 sometimes it occlusion |
---|
97 | disp('Nu of depthMap>FarestDist'); |
---|
98 | sum(sum(depthMap>FarestDist)) |
---|
99 | CleanedDepthMap(depthMap>FarestDist) = NaN; % don't clean the point >80 sometimes it occlusion |
---|
100 | size(CleanedDepthMap) |
---|
101 | size(depthMap) |
---|
102 | Posi3D = im_cr2w_cr(CleanedDepthMap,permute(Ray,[2 3 1])); |
---|
103 | if ~previoslyStored |
---|
104 | |
---|
105 | NewMap = [rand(max(Sup(:)),3); [0 0 0]]; |
---|
106 | % Clean the Sup near sky |
---|
107 | %maskSky = imdilate(maskSky, strel('disk', 3) ); |
---|
108 | %[Sup,MedSup]=CleanedSup(Sup,MedSup,maskSky); |
---|
109 | maskSky = Sup == 0; |
---|
110 | maskSkyEroded = imerode(maskSky, strel('disk', 4) ); |
---|
111 | SupEpand = ExpandSup2Sky(Sup,maskSkyEroded); |
---|
112 | NuPatch = HoriXNuDepth*VertYNuDepth-sum(maskSky(:)); |
---|
113 | |
---|
114 | NuSup = setdiff(unique(Sup)',0); |
---|
115 | NuSup = sort(NuSup); |
---|
116 | NuSupSize = size(NuSup,2); |
---|
117 | |
---|
118 | % Sup index and planeParameter index inverse map |
---|
119 | Sup2Para = sparse(1,max(Sup(:))); |
---|
120 | Sup2Para(NuSup) = 1:NuSupSize; |
---|
121 | |
---|
122 | % constructinf the Straight line prior matrix Will be add in the CoPlane matrix |
---|
123 | NuLine = size(StraightLineTable,2); |
---|
124 | CoPSTList = []; |
---|
125 | for i = 1:NuLine |
---|
126 | L = StraightLineTable{i}; |
---|
127 | X = L(1:(end-1))'; |
---|
128 | Y = L(2:end)'; |
---|
129 | if isempty(X) |
---|
130 | continue; |
---|
131 | end |
---|
132 | for j = 1:size(X,1) |
---|
133 | if X(j)~=Y(j) |
---|
134 | CoPSTList = [CoPSTList; X(j) Y(j)]; |
---|
135 | end |
---|
136 | end |
---|
137 | end |
---|
138 | end |
---|
139 | |
---|
140 | % Generate the Matrix for MRF |
---|
141 | tic |
---|
142 | PosiM = sparse(0,0); |
---|
143 | VarM = sparse(0,0); |
---|
144 | RayMd = sparse(0,0); |
---|
145 | RayAllOriM = sparse(0,0); |
---|
146 | RayAllM = sparse(0,0); |
---|
147 | RayMtilt = sparse(0,0); |
---|
148 | RayMCent = sparse(0,0); |
---|
149 | DepthInverseMCent = []; |
---|
150 | DepthInverseM = []; |
---|
151 | YPointer = []; |
---|
152 | beta = []; |
---|
153 | EmptyIndex = []; |
---|
154 | for i = NuSup |
---|
155 | % mask = Sup ==i; |
---|
156 | mask = SupEpand ==i; % include the Ray that will be use to expand the NonSky |
---|
157 | RayAllOriM = blkdiag( RayAllOriM, RayOri(:,mask)'); |
---|
158 | RayAllM = blkdiag( RayAllM, Ray(:,mask)'); |
---|
159 | mask = Sup ==i; % Not include the Ray that will be use to expand the NonSky |
---|
160 | [yt x] = find(mask); |
---|
161 | CenterX = round(median(x)); |
---|
162 | CenterY = round(median(yt)); |
---|
163 | YPointer = [YPointer; CenterY >= ceiling]; |
---|
164 | mask(isnan(CleanedDepthMap)) = false; |
---|
165 | SupNuPatch(i) = sum(mask(:)); |
---|
166 | % if sum(mask(:)) < 5 |
---|
167 | % EmptyIndex = [EmptyIndex; i]; |
---|
168 | % mask(:) = false; |
---|
169 | % end |
---|
170 | % find center point |
---|
171 | [yt x] = find(mask); |
---|
172 | CenterX = round(median(x)); |
---|
173 | CenterY = round(median(yt)); |
---|
174 | |
---|
175 | if ~all(mask(:)==0) |
---|
176 | if gravity |
---|
177 | Pa2CenterRatio = median(CleanedDepthMap(mask))./CleanedDepthMap(mask); |
---|
178 | if sum(mask(:)) > 0 |
---|
179 | RayMtilt = blkdiag(RayMtilt, ... |
---|
180 | ( Pa2CenterRatio(:,[1 1 1]).*repmat(RayOri(:,CenterY,CenterX)',[ SupNuPatch(i) 1])- RayOri(:,mask)')); |
---|
181 | else |
---|
182 | RayMtilt = blkdiag(RayMtilt, RayOri(:,mask)'); |
---|
183 | end |
---|
184 | RayMCent = blkdiag(RayMCent, RayOri(:,CenterY,CenterX)'*SupNuPatch(i)*mapVert(CenterY)*mapHori(CenterX)); |
---|
185 | PosiM = blkdiag(PosiM,Posi3D(:,mask)'.*repmat( mapVert(yt)',[1 3]).*repmat( mapHori(x)',[1 3])); |
---|
186 | VarM = [VarM; VarMap(mask).*(mapVert(yt)').*( mapHori(x)')]; |
---|
187 | RayMd = blkdiag(RayMd,RayOri(:,mask)'.*repmat( mapVert(yt)',[1 3]).*repmat( mapHori(x)',[1 3])); |
---|
188 | DepthInverseMCent = [DepthInverseMCent; 1./median(CleanedDepthMap(mask)).*SupNuPatch(i)* mapVert(CenterY)'*mapHori(CenterX)]; |
---|
189 | DepthInverseM = [DepthInverseM; 1./CleanedDepthMap(mask).* mapVert(yt)'.*mapHori(x)']; |
---|
190 | else |
---|
191 | Pa2CenterRatio = CleanedDepthMap(CenterY,CenterX)./CleanedDepthMap(mask); |
---|
192 | if sum(mask(:)) > 0 |
---|
193 | RayMtilt = blkdiag(RayMtilt, ... |
---|
194 | ( Pa2CenterRatio(:,[1 1 1]).*repmat(RayOri(:,CenterY,CenterX)',[ SupNuPatch(i) 1])- RayOri(:,mask)')); |
---|
195 | alse |
---|
196 | RayMtilt = blkdiag(RayMtilt, RayOri(:,mask)'); |
---|
197 | end |
---|
198 | RayMCent = blkdiag(RayMCent, RayOri(:,CenterY,CenterX)'*SupNuPatch(i)); |
---|
199 | PosiM = blkdiag(PosiM,Posi3D(:,mask)'); |
---|
200 | VarM = [VarM; VarMap(mask)]; |
---|
201 | RayMd = blkdiag(RayMd,RayOri(:,mask)'); |
---|
202 | DepthInverseMCent = [DepthInverseMCent; 1./median(CleanedDepthMap(mask)).*SupNuPatch(i)]; |
---|
203 | DepthInverseM = [DepthInverseM; 1./CleanedDepthMap(mask)]; |
---|
204 | end |
---|
205 | else |
---|
206 | RayMtilt = blkdiag(RayMtilt, RayOri(:,mask)'); |
---|
207 | RayMCent = blkdiag(RayMCent, RayOri(:,mask)'); |
---|
208 | PosiM = blkdiag(PosiM, Posi3D(:,mask)'); |
---|
209 | VarM = [VarM; VarMap(mask)]; |
---|
210 | RayMd = blkdiag(RayMd, RayOri(:,mask)'); |
---|
211 | end |
---|
212 | end |
---|
213 | YPointer(YPointer==0) = -1; |
---|
214 | |
---|
215 | % buliding CoPlane Matrix========================================================================= |
---|
216 | BounaryPHori = conv2(Sup,[1 -1],'same') ~=0; |
---|
217 | BounaryPHori(:,end) = 0; |
---|
218 | BounaryPVert = conv2(Sup,[1; -1],'same') ~=0; |
---|
219 | BounaryPVert(end,:) = 0; |
---|
220 | ClosestNList = [ Sup(find(BounaryPHori==1)) Sup(find(BounaryPHori==1)+VertYNuDepth);... |
---|
221 | Sup(find(BounaryPVert==1)) Sup(find(BounaryPVert==1)+1)]; |
---|
222 | ClosestNList = sort(ClosestNList,2); |
---|
223 | ClosestNList = unique(ClosestNList,'rows'); |
---|
224 | ClosestNList(ClosestNList(:,1) == 0,:) = []; |
---|
225 | BoundaryAll = BounaryPHori + BounaryPHori(:,[1 1:(end-1)])... |
---|
226 | +BounaryPVert + BounaryPVert([1 1:(end-1)],:); |
---|
227 | BoundaryAll([1 end],:) = 1; |
---|
228 | BoundaryAll(:,[1 end]) = 1; |
---|
229 | NuSTList = 0; |
---|
230 | if CoPST |
---|
231 | ClosestNList = [ClosestNList; CoPSTList]; |
---|
232 | NuSTList = size(CoPSTList,1) |
---|
233 | end |
---|
234 | NuNei = size(ClosestNList,1); |
---|
235 | CoPM1 = sparse(0,3*NuSupSize); |
---|
236 | CoPM2 = sparse(0,3*NuSupSize); |
---|
237 | CoPEstDepth = sparse(0,0); |
---|
238 | CoPNorM = []; |
---|
239 | WeiCoP = []; |
---|
240 | for i = 1: NuNei |
---|
241 | % if ~CornerList(i) |
---|
242 | mask = Sup == ClosestNList(i,1); |
---|
243 | SizeMaskAll = sum(mask(:)); |
---|
244 | [y x] = find(mask); |
---|
245 | CenterX = round(median(x)); |
---|
246 | CenterY = round(median(y)); |
---|
247 | y = find(mask(:,CenterX)); |
---|
248 | if ~isempty(y) |
---|
249 | CenterY = round(median(y)); |
---|
250 | end |
---|
251 | % CenterX = round(median(x)); |
---|
252 | % CenterY = round(median(yt)); |
---|
253 | |
---|
254 | temp1 = sparse(1, 3*NuSupSize); |
---|
255 | temp2 = sparse(1, 3*NuSupSize); |
---|
256 | temp1(:,(Sup2Para(ClosestNList(i,1))*3-2): Sup2Para(ClosestNList(i,1))*3) = Ray(:,CenterY,CenterX)'; |
---|
257 | temp2(:,(Sup2Para(ClosestNList(i,2))*3-2): Sup2Para(ClosestNList(i,2))*3) = Ray(:,CenterY,CenterX)'; |
---|
258 | NuNei-NuSTList |
---|
259 | i |
---|
260 | if i < NuNei-NuSTList % only immediate connecting superpixels are weighted using MultiScaleSup |
---|
261 | % wei = WeiV*10;%*(MultiScaleSupTable(Sup2Para(ClosestNList(i,1)),2:end) == MultiScaleSupTable(Sup2Para(ClosestNList(i,2)),2:end))'; |
---|
262 | if MultiScaleFlag |
---|
263 | vector = (MultiScaleSupTable(Sup2Para(ClosestNList(i,1)),2:end) == MultiScaleSupTable(Sup2Para(ClosestNList(i,2)),2:end)); |
---|
264 | expV = exp(-10*(WeiV*vector' + ShiftCoP) ); |
---|
265 | wei = 1/(1+expV); |
---|
266 | else |
---|
267 | wei = 1; |
---|
268 | end |
---|
269 | else |
---|
270 | wei = 1; |
---|
271 | end |
---|
272 | CoPM1 = [CoPM1; temp1*wei]; |
---|
273 | CoPM2 = [CoPM2; temp2*wei]; |
---|
274 | tempWeiCoP = [SizeMaskAll]; |
---|
275 | CoPEstDepth = [CoPEstDepth; max(median(CleanedDepthMap(mask)),ClosestDist)]; |
---|
276 | |
---|
277 | mask = Sup == ClosestNList(i,2); |
---|
278 | SizeMaskAll = sum(mask(:)); |
---|
279 | [y x] = find(mask); |
---|
280 | CenterX = round(median(x)); |
---|
281 | CenterY = round(median(y)); |
---|
282 | y = find(mask(:,CenterX)); |
---|
283 | if ~isempty(y) |
---|
284 | CenterY = round(median(y)); |
---|
285 | end |
---|
286 | % [yt x] = find(mask); |
---|
287 | % CenterX = round(median(x)); |
---|
288 | % CenterY = round(median(yt)); |
---|
289 | |
---|
290 | temp1 = sparse(1, 3*NuSupSize); |
---|
291 | temp2 = sparse(1, 3*NuSupSize); |
---|
292 | temp1(:,(Sup2Para(ClosestNList(i,1))*3-2): Sup2Para(ClosestNList(i,1))*3) = Ray(:,CenterY,CenterX)'; |
---|
293 | temp2(:,(Sup2Para(ClosestNList(i,2))*3-2): Sup2Para(ClosestNList(i,2))*3) = Ray(:,CenterY,CenterX)'; |
---|
294 | CoPM1 = [CoPM1; temp1*wei]; |
---|
295 | CoPM2 = [CoPM2; temp2*wei]; |
---|
296 | tempWeiCoP = [tempWeiCoP; SizeMaskAll]; |
---|
297 | WeiCoP = [WeiCoP; tempWeiCoP]; |
---|
298 | CoPEstDepth = [CoPEstDepth; max(median(CleanedDepthMap(mask)),ClosestDist)]; |
---|
299 | % end |
---|
300 | end%========================================================================================================= |
---|
301 | |
---|
302 | % find the boundary point that might need to be stick ot each other========================================== |
---|
303 | HoriStickM_i = sparse(0,3*NuSupSize); |
---|
304 | HoriStickM_j = sparse(0,3*NuSupSize); |
---|
305 | HoriStickPointInd = []; |
---|
306 | EstDepHoriStick = []; |
---|
307 | for i = find(BounaryPHori==1)' |
---|
308 | j = i+VertYNuDepth; |
---|
309 | if Sup(i) == 0 || Sup(j) == 0 |
---|
310 | continue; |
---|
311 | end |
---|
312 | % if ~OccluList(sum( ClosestNList == repmat(sort([Sup(i) Sup(j)]), [NuNei 1]),2) == 2) |
---|
313 | % size(OccluList) |
---|
314 | % if ~any(sum( OccluList == repmat(sort([Sup(i) Sup(j)]), [size(OccluList,1) 1]),2) == 2) |
---|
315 | Target(1) = Sup2Para(Sup(i)); |
---|
316 | Target(2) = Sup2Para(Sup(j)); |
---|
317 | rayBoundary(:,1) = Ray(:,i); |
---|
318 | rayBoundary(:,2) = Ray(:,i); |
---|
319 | % betaTemp = StickHori;%*(DistStickLengthNormWei.^2)*beta(Target(I)); |
---|
320 | % betaTemp = StickHori*WeiV;%*(MultiScaleSupTable(Sup2Para(Sup(i)),2:end) == MultiScaleSupTable(Sup2Para(Sup(j)),2:end))';%*(DistStickLengthNormWei.^2)*beta(Target(I)); |
---|
321 | if MultiScaleFlag |
---|
322 | vector = (MultiScaleSupTable(Sup2Para(Sup(i)),2:end) == MultiScaleSupTable(Sup2Para(Sup(j)),2:end)); |
---|
323 | expV = exp(-10*(WeiV*vector' + ShiftStick) ); |
---|
324 | betaTemp = StickHori*(0.5+1/(1+expV)); %*(DistStickLengthNormWei.^2)*beta(Target(I)); |
---|
325 | % therr should always be sticking (know don't care about occlusion) |
---|
326 | else |
---|
327 | betaTemp = StickHori; |
---|
328 | end |
---|
329 | temp = sparse(3,NuSupSize); |
---|
330 | temp(:,Target(1)) = rayBoundary(:,1); |
---|
331 | HoriStickM_i = [HoriStickM_i; betaTemp*temp(:)']; |
---|
332 | temp = sparse(3,NuSupSize); |
---|
333 | temp(:,Target(2)) = rayBoundary(:,2); |
---|
334 | HoriStickM_j = [HoriStickM_j; betaTemp*temp(:)']; |
---|
335 | EstDepHoriStick = [EstDepHoriStick; sqrt(max(CleanedDepthMap(i),ClosestDist)*max(CleanedDepthMap(j),ClosestDist))]; |
---|
336 | HoriStickPointInd = [HoriStickPointInd i ]; |
---|
337 | % else |
---|
338 | % disp('Occlu'); |
---|
339 | % end |
---|
340 | end |
---|
341 | VertStickM_i = sparse(0,3*NuSupSize); |
---|
342 | VertStickM_j = sparse(0,3*NuSupSize); |
---|
343 | VertStickPointInd = []; |
---|
344 | EstDepVertStick = []; |
---|
345 | for i = find(BounaryPVert==1)' |
---|
346 | j = i+1; |
---|
347 | if Sup(i) == 0 || Sup(j) == 0 |
---|
348 | continue; |
---|
349 | end |
---|
350 | % if ~OccluList(sum( ClosestNList == repmat(sort([Sup(i) Sup(j)]), [NuNei 1]),2) == 2) |
---|
351 | % if ~any(sum( OccluList == repmat(sort([Sup(i) Sup(j)]), [size(OccluList,1) 1]),2) == 2) |
---|
352 | Target(1) = Sup2Para(Sup(i)); |
---|
353 | Target(2) = Sup2Para(Sup(j)); |
---|
354 | rayBoundary(:,1) = Ray(:,i); |
---|
355 | rayBoundary(:,2) = Ray(:,i); |
---|
356 | % betaTemp = StickVert;%*(DistStickLengthNormWei.^2)*beta(Target(I)); |
---|
357 | if MultiScaleFlag |
---|
358 | vector = (MultiScaleSupTable(Sup2Para(Sup(i)),2:end) == MultiScaleSupTable(Sup2Para(Sup(j)),2:end)); |
---|
359 | expV = exp(-10*(WeiV*vector' + ShiftStick) ); |
---|
360 | betaTemp = StickVert*(0.5+1/(1+expV)); |
---|
361 | % therr should always be sticking (know don't care about occlusion) |
---|
362 | else |
---|
363 | betaTemp = StickVert; |
---|
364 | end |
---|
365 | temp = sparse(3,NuSupSize); |
---|
366 | temp(:,Target(1)) = rayBoundary(:,1); |
---|
367 | VertStickM_i = [VertStickM_i; betaTemp*temp(:)']; |
---|
368 | temp = sparse(3,NuSupSize); |
---|
369 | temp(:,Target(2)) = rayBoundary(:,2); |
---|
370 | VertStickM_j = [VertStickM_j; betaTemp*temp(:)']; |
---|
371 | EstDepVertStick = [EstDepVertStick; sqrt(max(CleanedDepthMap(i),ClosestDist)*max(CleanedDepthMap(j),ClosestDist))]; |
---|
372 | VertStickPointInd = [VertStickPointInd i ]; |
---|
373 | % else |
---|
374 | % disp('Occlu'); |
---|
375 | % end |
---|
376 | end |
---|
377 | % ======================================Finish building up matrix=====================hard work====================== |
---|
378 | % ================================================================================================================ |
---|
379 | depthfile = strrep(filename{k},'img','depth_learned'); % |
---|
380 | |
---|
381 | for i = step |
---|
382 | opt = sdpsettings('solver','sedumi'); |
---|
383 | ParaPPCP = sdpvar(3*NuSupSize,1); |
---|
384 | F = set(ParaPPCP(3*(1:NuSupSize)-1).*YPointer<=0) ... |
---|
385 | +set(RayAllM*ParaPPCP >=1/FarestDist)... |
---|
386 | +set(RayAllM*ParaPPCP <=1/ClosestDist)... |
---|
387 | +set(RayAllOriM*ParaPPCP >=1/FarestDist)... |
---|
388 | +set(RayAllOriM*ParaPPCP <=1/ClosestDist); |
---|
389 | % First fit the plane to find the estimated plane parameters |
---|
390 | % If PosiM contain NaN data the corresponding Plane Parameter will also be NaN |
---|
391 | if i == 1 % W/O priors |
---|
392 | solvesdp(F,norm( (PosiM*ParaPPCP-ones(size(PosiM,1),1))./(abs(VarM)),1)... |
---|
393 | , opt); |
---|
394 | elseif i ==2 % Coner |
---|
395 | solvesdp(F,norm( (PosiM*ParaPPCP-ones(size(PosiM,1),1))./(abs(VarM)),1)... |
---|
396 | +Center*norm(((CoPM1 - CoPM2)*ParaPPCP).*CoPEstDepth, 1)... |
---|
397 | , opt); |
---|
398 | else % all the priors |
---|
399 | sum(sum(isnan(PosiM))) |
---|
400 | disp('Whole MRF') |
---|
401 | solvesdp(F,norm( (PosiM*ParaPPCP-ones(size(PosiM,1),1))./(abs(VarM)),1)... |
---|
402 | +Center*norm(((CoPM1 - CoPM2)*ParaPPCP).*CoPEstDepth, 1)... |
---|
403 | +norm(((HoriStickM_i-HoriStickM_j)*ParaPPCP).*EstDepHoriStick,1)... |
---|
404 | +norm(((VertStickM_i-VertStickM_j)*ParaPPCP).*EstDepVertStick,1)... |
---|
405 | , opt); |
---|
406 | end |
---|
407 | |
---|
408 | ParaPPCP = double(ParaPPCP); |
---|
409 | isnan(ParaPPCP) |
---|
410 | SepPointMeasureHori = 1./(HoriStickM_i*ParaPPCP)-1./(HoriStickM_j*ParaPPCP); |
---|
411 | SepPointMeasureVert = 1./(VertStickM_i*ParaPPCP)-1./(VertStickM_j*ParaPPCP); |
---|
412 | ParaPPCP = reshape(ParaPPCP,3,[]); |
---|
413 | %any(any(isnan(ParaPPCP))) |
---|
414 | % porject the ray on planes to generate the ProjDepth |
---|
415 | FitDepthPPCP = FarestDist*ones(1,VertYNuDepth*HoriXNuDepth); |
---|
416 | FitDepthPPCP(~maskSkyEroded) = (1./sum(ParaPPCP(:,Sup2Para(SupEpand(~maskSkyEroded))).*Ray(:,~maskSkyEroded),1))'; |
---|
417 | FitDepthPPCP = reshape(FitDepthPPCP,VertYNuDepth,[]); |
---|
418 | % storage for comparison purpose ====================== |
---|
419 | depthMap = FarestDist*ones(1,VertYNuDepth*HoriXNuDepth); |
---|
420 | depthMap(~maskSkyEroded) = (1./sum(ParaPPCP(:,Sup2Para(SupEpand(~maskSkyEroded))).*RayOri(:,~maskSkyEroded),1))'; |
---|
421 | depthMap = reshape(depthMap,VertYNuDepth,[]); |
---|
422 | % if baseline == 0 |
---|
423 | system(['mkdir ' ScratchDataFolder '/' Name{i} '_' DepthFolder '/']); |
---|
424 | save([ScratchDataFolder '/' Name{i} '_' DepthFolder '/' depthfile '.mat'],'depthMap'); |
---|
425 | % elseif baseline == 1 |
---|
426 | % system(['mkdir ' ScratchDataFolder '/' Name{i} '_' DepthFolder '_baseline/']); |
---|
427 | % save([ScratchDataFolder '/' Name{i} '_' DepthFolder '_baseline/' depthfile '.mat'],'depthMap'); |
---|
428 | % else |
---|
429 | % system(['mkdir ' ScratchDataFolder '/' Name{i} '_' DepthFolder '_baseline2/']); |
---|
430 | % save([ScratchDataFolder '/' Name{i} '_' DepthFolder '_baseline2/' depthfile '.mat'],'depthMap'); |
---|
431 | % end |
---|
432 | |
---|
433 | % ===================================================== |
---|
434 | [Position3DFitedPPCP] = im_cr2w_cr(FitDepthPPCP,permute(Ray,[2 3 1])); |
---|
435 | if false %RenderVrmlFlag && i==3 |
---|
436 | [VrmlName] = vrml_test_faceset_goodSkyBoundary( filename{k}, Position3DFitedPPCP, FitDepthPPCP, permute(Ray,[2 3 1]), [Name{i} DepthFolder], ... |
---|
437 | [], [], 0, maskSkyEroded, maskG, 1, 0, a_default, b_default, Ox_default, Oy_default); |
---|
438 | system(['gzip -9 -c ' ScratchDataFolder '/vrml/' VrmlName ' > ' ScratchDataFolder '/vrml/' VrmlName '.gz']); |
---|
439 | delete([ScratchDataFolder '/vrml/' VrmlName]); |
---|
440 | end |
---|
441 | % save([ScratchDataFolder '/data/PreDepth/FitDepthPPCP' num2str(k) '.mat'],'FitDepthPPCP','SepPointMeasureHori',... |
---|
442 | % 'SepPointMeasureVert','VertStickPointInd','HoriStickPointInd'); |
---|
443 | end |
---|
444 | |
---|
445 | %return; |
---|
446 | %==================Finished for one step MRF========================================================================================================== |
---|
447 | |
---|
448 | |
---|
449 | % ====================================following are 2nd and 3rd step MRF to give more visually pleasing result======================================= |
---|
450 | ;% generating new PosiMPPCP using the new position |
---|
451 | |
---|
452 | save([ScratchDataFolder '/data/VertShiftVrml' mat2str(k) '.mat' ] ); |
---|
453 | groundThreshold = cos(15*pi/180); |
---|
454 | verticalThreshold = cos(50*pi/180); % before Nov 29 11:30Pm I used 30 which is too big |
---|
455 | normPara = norms(ParaPPCP); |
---|
456 | normalizedPara = ParaPPCP ./ repmat( normPara, [3 1]); |
---|
457 | groundPara = abs(normalizedPara(2,:)) >= groundThreshold; |
---|
458 | groundParaInd = find(groundPara); |
---|
459 | verticalPara = abs(normalizedPara(2,:)) <= verticalThreshold; |
---|
460 | verticalParaInd = find(verticalPara); |
---|
461 | |
---|
462 | if false |
---|
463 | [CoPEstDepth CoPM1 CoPM2 EstDepHoriStick HoriStickM_i HoriStickM_j EstDepVertStick VertStickM_i VertStickM_j WeiCoP NewSup NewNuSup NewNuSupSize NewSup2Para FixPara VertStickPointInd HoriStickPointInd] ... |
---|
464 | = FreeSupSharpCorner( Sup, Ray, OccluList, WeiV, ShiftStick, StickHori, StickVert, ClosestDist, depthMap,... |
---|
465 | MultiScaleSupTable, verticalParaInd, graoundParaInd, MultiScaleFlag); |
---|
466 | |
---|
467 | end |
---|
468 | |
---|
469 | indexVertical = find( verticalPara)*3-1; |
---|
470 | indexGroundX = find( groundPara)*3-2; |
---|
471 | indexGroundZ = find( groundPara)*3; |
---|
472 | |
---|
473 | PosiMPPCP = sparse(0,0); |
---|
474 | VarM2 = sparse(0,0); |
---|
475 | % VertVar = 10^(-2); |
---|
476 | for i = NuSup |
---|
477 | mask = Sup == i; |
---|
478 | if any(verticalParaInd == Sup2Para(i)) |
---|
479 | mask = logical(zeros(size(Sup))); |
---|
480 | end |
---|
481 | mask(ipermute( isnan( sum(Posi3D,1) ),[2 3 1] ) ) = false; % new added to avoid NaN |
---|
482 | % PosiMPPCP = blkdiag(PosiMPPCP, Posi3D(:,mask)'); |
---|
483 | % VarM2 = [VarM2; VertVar*ones(sum(mask(:)), 1)]; |
---|
484 | % mask = logical(zeros(size(Sup))); |
---|
485 | % elseif any(groundParaInd == Sup2Para(i)) |
---|
486 | % PosiMPPCP = blkdiag(PosiMPPCP, Position3DFitedPPCP(:,mask)'); |
---|
487 | % VarM2 = [VarM2; VarMap(mask)]; |
---|
488 | % else |
---|
489 | PosiMPPCP = blkdiag(PosiMPPCP, Posi3D(:,mask)'); |
---|
490 | VarM2 = [VarM2; VarMap(mask)]; |
---|
491 | % end |
---|
492 | end |
---|
493 | |
---|
494 | opt = sdpsettings('solver','sedumi'); |
---|
495 | Para = sdpvar(3*NuSupSize,1); |
---|
496 | F = set(Para(3*(1:NuSupSize)-1).*YPointer<=0) + set(RayAllM*Para >=1/FarestDist)... |
---|
497 | +set(RayAllM*Para <=1/ClosestDist)... |
---|
498 | +set(Para(indexVertical) == 0); |
---|
499 | % +set(Para(indexGroundX) == 0)... |
---|
500 | % +set(Para(indexGroundZ) == 0); |
---|
501 | if FractionalDepthError |
---|
502 | size(PosiMPPCP) |
---|
503 | solvesdp(F,norm( ( PosiMPPCP*Para-ones(size(PosiMPPCP,1),1))./VarM2,1)... |
---|
504 | +Center*norm(((CoPM1 - CoPM2)*Para).*CoPEstDepth, 1)... |
---|
505 | +norm(((HoriStickM_i-HoriStickM_j)*Para).*EstDepHoriStick,1)... |
---|
506 | +norm(((VertStickM_i-VertStickM_j)*Para).*EstDepVertStick,1)... |
---|
507 | , opt); |
---|
508 | % +1000*norm( Para(indexGroundX)./ normPara( find(groundPara) )', 1) ... |
---|
509 | % +1000*norm( Para(indexGroundZ)./ normPara( find(groundPara) )', 1) ... |
---|
510 | % sqrt( max(CoPM1*ParaPPCP(:), 1/FarestDist).*max(CoPM2*ParaPPCP(:), 1/FarestDist)), 1)... |
---|
511 | % +norm(((VertStickM_i-VertStickM_j)*Para)./... |
---|
512 | % sqrt( max(VertStickM_i*ParaPPCP(:), 1/FarestDist).*max(VertStickM_j*ParaPPCP(:), 1/FarestDist)),1)... |
---|
513 | % pause; |
---|
514 | % +1000*norm( Para(indexVertical)./ normPara( find(verticalPara) )',1) ... |
---|
515 | % +norm(((HoriStickM_i-HoriStickM_j)*Para)./sqrt((VertStickM_i*ParaPPCP(:)).*(VertStickM_j*ParaPPCP(:))),1)... |
---|
516 | else |
---|
517 | solvesdp(F,norm( RayMd*Para - DepthInverseM,1)... |
---|
518 | +Center*norm((CoPM1 - CoPM2)*Para,1)... |
---|
519 | +norm((VertStickM_i-VertStickM_j)*Para,1)... |
---|
520 | +norm((HoriStickM_i-HoriStickM_j)*Para,1)... |
---|
521 | +1000*norm( Para(indexVertical)./ normPara( find(verticalPara) )',1) ... |
---|
522 | +1000*norm( Para(indexGroundX)./ normPara( find(groundPara) )', 1) ... |
---|
523 | +1000*norm( Para(indexGroundZ)./ normPara( find(groundPara) )', 1) ... |
---|
524 | , opt); |
---|
525 | % +0.01*norm( RayMtilt*ParaPPCP,1)... |
---|
526 | % +norm(spdiags(WeiCoP,0,size(CoPM1,1),size(CoPM1,1))*(CoPM1 - CoPM2)*ParaPPCP,1)... |
---|
527 | % +norm( (CoPM1 - CoPM2)*ParaPPCP,1 )... |
---|
528 | end |
---|
529 | %save([ScratchDataFolder '/data/All.mat']); |
---|
530 | Para = double(Para); |
---|
531 | SepPointMeasureHori = 1./(HoriStickM_i*Para)-1./(HoriStickM_j*Para); |
---|
532 | SepPointMeasureVert = 1./(VertStickM_i*Para)-1./(VertStickM_j*Para); |
---|
533 | Para = reshape(Para,3,[]); |
---|
534 | |
---|
535 | FitDepth = 80*ones(1,VertYNuDepth*HoriXNuDepth); |
---|
536 | FitDepth(~maskSkyEroded) = (1./sum(Para(:,Sup2Para(SupEpand(~maskSkyEroded))).*Ray(:,~maskSkyEroded),1))'; |
---|
537 | FitDepth = reshape(FitDepth,VertYNuDepth,[]); |
---|
538 | sum(isnan(FitDepth(:))) |
---|
539 | [Position3DFited] = im_cr2w_cr(FitDepth,permute(Ray,[2 3 1])); |
---|
540 | % [VrmlName] = vrml_test_faceset_goodSkyBoundary( filename{k}, Position3DFited, FitDepth, permute(Ray,[2 3 1]), 'VNonSupport', ... |
---|
541 | % [], [], 0, maskSkyEroded, maskG, 1, 0, a_default, b_default, Ox_default, Oy_default); |
---|
542 | % system(['gzip -9 -c ' ScratchDataFolder '/vrml/' VrmlName ' > ' ScratchDataFolder '/vrml/' VrmlName '.gz']); |
---|
543 | % delete([ScratchDataFolder '/vrml/' VrmlName]); |
---|
544 | |
---|
545 | %end |
---|
546 | % ============================================================ |
---|
547 | if ConerImprove |
---|
548 | save([ScratchDataFolder '/data/All.mat']); |
---|
549 | |
---|
550 | groundThreshold = cos(5*pi/180); % a little bigger than 15 |
---|
551 | verticalThreshold = cos(50*pi/180); |
---|
552 | normPara = norms(Para); |
---|
553 | normalizedPara = Para ./ repmat( normPara, [3 1]); |
---|
554 | groundPara = abs(normalizedPara(2,:)) >= groundThreshold; |
---|
555 | groundParaInd = find(groundPara); |
---|
556 | verticalPara = abs(normalizedPara(2,:)) <= verticalThreshold; |
---|
557 | verticalParaInd = find(verticalPara); |
---|
558 | [CoPEstDepth CoPM1 CoPM2 EstDepHoriStick HoriStickM_i HoriStickM_j EstDepVertStick VertStickM_i VertStickM_j ... |
---|
559 | WeiCoP NewSup NewNuSup NewNuSupSize NewSup2Para FixPara VertStickPointInd HoriStickPointInd] ... |
---|
560 | = FreeSupSharpCorner( Sup, Ray, OccluList, WeiV, ShiftStick, StickHori, StickVert, ClosestDist, depthMap,... |
---|
561 | MultiScaleSupTable, verticalParaInd, groundParaInd, MultiScaleFlag, Sup2Para); |
---|
562 | |
---|
563 | % [CoPM1 CoPM2 HoriStickM_i HoriStickM_j VertStickM_i VertStickM_j WeiCoP NewSup NewNuSup NewNuSupSize NewSup2Para FixPara VertStickPointInd HoriStickPointInd] ... |
---|
564 | % = FreeSupSharpCorner( Sup, Ray, OccluList, WeiV, ShiftStick, StickHori, StickVert,... |
---|
565 | % MultiScaleSupTable, verticalParaInd, groundParaInd, MultiScaleFlag, Sup2Para); |
---|
566 | %[CoPM1 CoPM2 WeiCoP Sup Sup2ParaNew FixPara NewNuSup VertStickPointInd HoriStickPointInd] ... |
---|
567 | % = FixParaFreeCorner( Sup, Ray, HoriStickPointInd, VertStickPointInd,... |
---|
568 | % SepPointMeasureHori, SepPointMeasureVert, OccluList, 1); % may change 1 to logistic later |
---|
569 | |
---|
570 | %[CoPM1 CoPM2 WeiCoP Sup Sup2ParaNew FixPara NewNuSup VertStickPointInd HoriStickPointInd] ... |
---|
571 | % = FixParaFreeCornerNew( Sup, Ray, verticalParaInd, groundParaInd, NuSup, NuSupSize,OccluList, 1); % may change 1 to logistic later |
---|
572 | |
---|
573 | indexVertical = NewSup2Para(NuSup(verticalParaInd)); |
---|
574 | indexGround = NewSup2Para(NuSup(groundParaInd)); |
---|
575 | % indexGroundZ = Sup2ParaNew(NuSup(groundParaInd)); |
---|
576 | ZeroMaskVertical = indexVertical ==0; |
---|
577 | ZeroMaskGround = indexGround ==0; |
---|
578 | % ZeroMaskGroundZ = indexGroundZ ==0; |
---|
579 | indexVertical = indexVertical(~ZeroMaskVertical); |
---|
580 | indexGround = indexGround(~ZeroMaskGround); |
---|
581 | % indexGroundZ = indexGroundZ(~ZeroMaskGroundZ); |
---|
582 | indexVertical = indexVertical*3-1; |
---|
583 | indexGroundX = indexGround*3-2; |
---|
584 | indexGroundZ = indexGround*3; |
---|
585 | |
---|
586 | |
---|
587 | RayAllM = sparse(0,0); |
---|
588 | %PosiM = sparse(0,0); |
---|
589 | RayM = sparse(0,0); |
---|
590 | DepthIveEst = []; |
---|
591 | YPointer = []; |
---|
592 | for i = NewNuSup |
---|
593 | mask = NewSup ==i; |
---|
594 | RayAllM = blkdiag( RayAllM, Ray(:,mask)'); |
---|
595 | [yt x] = find(mask); |
---|
596 | CenterX = round(median(x)); |
---|
597 | CenterY = round(median(yt)); |
---|
598 | YPointer = [YPointer; CenterY >= ceiling]; |
---|
599 | if ~any( FixPara ==i ) |
---|
600 | mask = logical(zeros(size(Sup))); |
---|
601 | end |
---|
602 | RayM = blkdiag( RayM, Ray(:,mask)'); |
---|
603 | DepthIveEst = [DepthIveEst; 1./FitDepth(mask)]; |
---|
604 | end |
---|
605 | PosiM = spdiags(1./DepthIveEst,0,size(DepthIveEst,1),size(DepthIveEst,1))*RayM; |
---|
606 | |
---|
607 | NewNuSupSize = size(NewNuSup,2); |
---|
608 | opt = sdpsettings('solver','sedumi'); |
---|
609 | ParaCorner = sdpvar(3*NewNuSupSize,1); |
---|
610 | F = set(ParaCorner(3*(1:NewNuSupSize)-1).*YPointer<=0) + set(RayAllM*ParaCorner >=1/FarestDist)... |
---|
611 | +set(RayAllM*ParaCorner <=1/ClosestDist)... |
---|
612 | +set(ParaCorner([(NewSup2Para(FixPara)*3-2); (NewSup2Para(FixPara)*3-1); (NewSup2Para(FixPara)*3)]) == ... |
---|
613 | Para([(NewSup2Para(FixPara)*3-2); (NewSup2Para(FixPara)*3-1); (NewSup2Para(FixPara)*3)]));; |
---|
614 | % +set(ParaCorner(indexVertical) == 0); |
---|
615 | % +set(ParaCorner(indexGroundX) == 0); |
---|
616 | % +set(ParaCorner(indexGroundZ) == 0)... |
---|
617 | % +set(ParaCorner([(Sup2ParaNew(FixPara)*3-2); (Sup2ParaNew(FixPara)*3-1); (Sup2ParaNew(FixPara)*3)]) == ... |
---|
618 | % Para([(Sup2Para(FixPara)*3-2); (Sup2Para(FixPara)*3-1); (Sup2Para(FixPara)*3)])); |
---|
619 | solvesdp(F,Center*norm( ((CoPM1 - CoPM2)*ParaCorner), 1 )... |
---|
620 | +1000*norm( ParaCorner(indexGroundX)./ normPara( groundParaInd(~ZeroMaskGround) )', 1) ... |
---|
621 | +1000*norm( ParaCorner(indexGroundZ)./ normPara( groundParaInd(~ZeroMaskGround) )', 1) ... |
---|
622 | , opt); |
---|
623 | %+norm(((HoriStickM_i-HoriStickM_j)*ParaCorner).*EstDepHoriStick ,1)... |
---|
624 | %+norm(((VertStickM_i-VertStickM_j)*ParaCorner).*EstDepVertStick ,1)... |
---|
625 | % CoPEstDepth |
---|
626 | % PosiM*ParaCorner - 1,1) |
---|
627 | % norm(RayM*ParaCorner - DepthIveEst,1) |
---|
628 | |
---|
629 | %pause; |
---|
630 | ParaCorner = double(ParaCorner); |
---|
631 | ParaCorner = reshape(ParaCorner,3,[]); |
---|
632 | FitDepthCorner = 100*ones(1,VertYNuDepth*HoriXNuDepth); |
---|
633 | FitDepthCorner(NewSup~=0) = (1./sum(ParaCorner(:,NewSup2Para(NewSup(NewSup~=0))).*Ray(:,NewSup~=0),1))'; |
---|
634 | FitDepthCorner = reshape(FitDepthCorner,VertYNuDepth,[]); |
---|
635 | [Position3DFitedCorner] = im_cr2w_cr(FitDepthCorner,permute(Ray,[2 3 1])); |
---|
636 | [VrmlName] = vrml_test_faceset_goodSkyBoundary( filename{k}, Position3DFitedCorner, FitDepthCorner, permute(Ray,[2 3 1]), ['Corner'], ... |
---|
637 | [], [], 0, maskSkyEroded, maskG, 1, 0, a_default, b_default, Ox_default, Oy_default); |
---|
638 | system(['gzip -9 -c ' ScratchDataFolder '/vrml/' VrmlName ' > ' ScratchDataFolder '/vrml/' VrmlName '.gz']); |
---|
639 | delete([ScratchDataFolder '/vrml/' VrmlName]); |
---|
640 | end |
---|
641 | %return; |
---|
642 | return; |
---|