[37] | 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(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 = 2; % 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 = 0.01; % set the confidant of the learned depth at the top of the image |
---|
| 61 | BandWith = 1; % Nov29 1 Nov30 0.1 check result change to 0.1 12/1 0.1 lost detail |
---|
| 62 | mapVert = linspace(VertConf,1,VertYNuDepth); % modeling the gravity prior |
---|
| 63 | mapHori = [linspace(HoriConf,1,round(HoriXNuDepth/2)) fliplr(linspace(HoriConf,1,HoriXNuDepth-round(HoriXNuDepth/2)))]; % modeling the gravity prior |
---|
| 64 | % ========set the range of depth that our model in |
---|
| 65 | ClosestDist = 1; |
---|
| 66 | % set the FarestDist to very 5 times to median depth |
---|
| 67 | FarestDist = 1.5*median(depthMap(:)); % tried on university % nogood effect but keep it since usually it the rangelike this % change to 1.5 for church |
---|
| 68 | % ================================================ |
---|
| 69 | ceiling = 0*VertYNuDepth; % set the position of the ceiling, related to No plane coming back constrain % changed for newchurch |
---|
| 70 | Name{1} = 'FraWOPri'; |
---|
| 71 | Name{2} = 'FraCoP'; |
---|
| 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 = false; |
---|
| 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 | % =========================================================================================================================================== |
---|
| 143 | groundThreshold = cos([ zeros(1, VertYNuDepth - ceil(VertYNuDepth/2)+10) linspace(0,15,ceil(VertYNuDepth/2)-10)]*pi/180); |
---|
| 144 | % v1 15 v2 20 too big v3 20 to ensure non misclassified as ground. |
---|
| 145 | % verticalThreshold = cos(linspace(5,55,VertYNuDepth)*pi/180); % give a vector of size 55 in top to down : |
---|
| 146 | verticalThreshold = cos([ 5*ones(1,VertYNuDepth - ceil(VertYNuDepth/2)) linspace(5,55,ceil(VertYNuDepth/2))]*pi/180); % give a vector of size 55 in top to down : |
---|
| 147 | % 50 means suface norm away from y axis more than 50 degree |
---|
| 148 | % =========================================================================================================================================== |
---|
| 149 | PosiM = sparse(0,0); |
---|
| 150 | VarM = sparse(0,0); |
---|
| 151 | RayMd = sparse(0,0); |
---|
| 152 | RayAllOriM = sparse(0,0); |
---|
| 153 | RayAllM = sparse(0,0); |
---|
| 154 | RayMtilt = sparse(0,0); |
---|
| 155 | RayMCent = sparse(0,0); |
---|
| 156 | DepthInverseMCent = []; |
---|
| 157 | DepthInverseM = []; |
---|
| 158 | YPointer = []; |
---|
| 159 | YPosition = []; |
---|
| 160 | beta = []; |
---|
| 161 | EmptyIndex = []; |
---|
| 162 | for i = NuSup |
---|
| 163 | % mask = Sup ==i; |
---|
| 164 | mask = SupEpand ==i; % include the Ray that will be use to expand the NonSky |
---|
| 165 | RayAllOriM = blkdiag( RayAllOriM, RayOri(:,mask)'); |
---|
| 166 | RayAllM = blkdiag( RayAllM, Ray(:,mask)'); |
---|
| 167 | mask = Sup ==i; % Not include the Ray that will be use to expand the NonSky |
---|
| 168 | [yt x] = find(mask); |
---|
| 169 | CenterX = round(median(x)); |
---|
| 170 | CenterY = round(median(yt)); |
---|
| 171 | YPointer = [YPointer; CenterY >= ceiling]; |
---|
| 172 | YPosition = [YPosition; CenterY]; |
---|
| 173 | mask(isnan(CleanedDepthMap)) = false; |
---|
| 174 | SupNuPatch(i) = sum(mask(:)); |
---|
| 175 | % if sum(mask(:)) < 5 |
---|
| 176 | % EmptyIndex = [EmptyIndex; i]; |
---|
| 177 | % mask(:) = false; |
---|
| 178 | % end |
---|
| 179 | % find center point |
---|
| 180 | [yt x] = find(mask); |
---|
| 181 | CenterX = round(median(x)); |
---|
| 182 | CenterY = round(median(yt)); |
---|
| 183 | |
---|
| 184 | if ~all(mask(:)==0) |
---|
| 185 | if gravity |
---|
| 186 | Pa2CenterRatio = median(CleanedDepthMap(mask))./CleanedDepthMap(mask); |
---|
| 187 | if sum(mask(:)) > 0 |
---|
| 188 | RayMtilt = blkdiag(RayMtilt, ... |
---|
| 189 | ( Pa2CenterRatio(:,[1 1 1]).*repmat(RayOri(:,CenterY,CenterX)',[ SupNuPatch(i) 1])- RayOri(:,mask)')); |
---|
| 190 | else |
---|
| 191 | RayMtilt = blkdiag(RayMtilt, RayOri(:,mask)'); |
---|
| 192 | end |
---|
| 193 | RayMCent = blkdiag(RayMCent, RayOri(:,CenterY,CenterX)'*SupNuPatch(i)*mapVert(CenterY)*mapHori(CenterX)); |
---|
| 194 | PosiM = blkdiag(PosiM,Posi3D(:,mask)');%.*repmat( mapVert(yt)',[1 3]).*repmat( mapHori(x)',[1 3])); |
---|
| 195 | VarM = [VarM; VarMap(mask).*(mapVert(yt)').*( mapHori(x)')]; |
---|
| 196 | RayMd = blkdiag(RayMd,RayOri(:,mask)'.*repmat( mapVert(yt)',[1 3]).*repmat( mapHori(x)',[1 3])); |
---|
| 197 | DepthInverseMCent = [DepthInverseMCent; 1./median(CleanedDepthMap(mask)).*SupNuPatch(i)* mapVert(CenterY)'*mapHori(CenterX)]; |
---|
| 198 | DepthInverseM = [DepthInverseM; 1./CleanedDepthMap(mask).* mapVert(yt)'.*mapHori(x)']; |
---|
| 199 | else |
---|
| 200 | Pa2CenterRatio = CleanedDepthMap(CenterY,CenterX)./CleanedDepthMap(mask); |
---|
| 201 | if sum(mask(:)) > 0 |
---|
| 202 | RayMtilt = blkdiag(RayMtilt, ... |
---|
| 203 | ( Pa2CenterRatio(:,[1 1 1]).*repmat(RayOri(:,CenterY,CenterX)',[ SupNuPatch(i) 1])- RayOri(:,mask)')); |
---|
| 204 | else |
---|
| 205 | RayMtilt = blkdiag(RayMtilt, RayOri(:,mask)'); |
---|
| 206 | end |
---|
| 207 | RayMCent = blkdiag(RayMCent, RayOri(:,CenterY,CenterX)'*SupNuPatch(i)); |
---|
| 208 | PosiM = blkdiag(PosiM,Posi3D(:,mask)'); |
---|
| 209 | VarM = [VarM; VarMap(mask)]; |
---|
| 210 | RayMd = blkdiag(RayMd,RayOri(:,mask)'); |
---|
| 211 | DepthInverseMCent = [DepthInverseMCent; 1./median(CleanedDepthMap(mask)).*SupNuPatch(i)]; |
---|
| 212 | DepthInverseM = [DepthInverseM; 1./CleanedDepthMap(mask)]; |
---|
| 213 | end |
---|
| 214 | else |
---|
| 215 | RayMtilt = blkdiag(RayMtilt, RayOri(:,mask)'); |
---|
| 216 | RayMCent = blkdiag(RayMCent, RayOri(:,mask)'); |
---|
| 217 | PosiM = blkdiag(PosiM, Posi3D(:,mask)'); |
---|
| 218 | VarM = [VarM; VarMap(mask)]; |
---|
| 219 | RayMd = blkdiag(RayMd, RayOri(:,mask)'); |
---|
| 220 | end |
---|
| 221 | end |
---|
| 222 | YPointer(YPointer==0) = -1; |
---|
| 223 | |
---|
| 224 | % buliding CoPlane Matrix========================================================================= |
---|
| 225 | BounaryPHori = conv2(Sup,[1 -1],'same') ~=0; |
---|
| 226 | BounaryPHori(:,end) = 0; |
---|
| 227 | BounaryPVert = conv2(Sup,[1; -1],'same') ~=0; |
---|
| 228 | BounaryPVert(end,:) = 0; |
---|
| 229 | ClosestNList = [ Sup(find(BounaryPHori==1)) Sup(find(BounaryPHori==1)+VertYNuDepth);... |
---|
| 230 | Sup(find(BounaryPVert==1)) Sup(find(BounaryPVert==1)+1)]; |
---|
| 231 | ClosestNList = sort(ClosestNList,2); |
---|
| 232 | ClosestNList = unique(ClosestNList,'rows'); |
---|
| 233 | ClosestNList(ClosestNList(:,1) == 0,:) = []; |
---|
| 234 | BoundaryAll = BounaryPHori + BounaryPHori(:,[1 1:(end-1)])... |
---|
| 235 | +BounaryPVert + BounaryPVert([1 1:(end-1)],:); |
---|
| 236 | BoundaryAll([1 end],:) = 1; |
---|
| 237 | BoundaryAll(:,[1 end]) = 1; |
---|
| 238 | NuSTList = 0; |
---|
| 239 | if CoPST |
---|
| 240 | ClosestNList = [ClosestNList; CoPSTList]; |
---|
| 241 | NuSTList = size(CoPSTList,1); |
---|
| 242 | end |
---|
| 243 | NuNei = size(ClosestNList,1); |
---|
| 244 | CoPM1 = sparse(0,3*NuSupSize); |
---|
| 245 | CoPM2 = sparse(0,3*NuSupSize); |
---|
| 246 | CoPEstDepth = sparse(0,0); |
---|
| 247 | CoPNorM = []; |
---|
| 248 | WeiCoP = []; |
---|
| 249 | for i = 1: NuNei |
---|
| 250 | % if ~CornerList(i) |
---|
| 251 | mask = Sup == ClosestNList(i,1); |
---|
| 252 | SizeMaskAll = sum(mask(:)); |
---|
| 253 | [y x] = find(mask); |
---|
| 254 | CenterX = round(median(x)); |
---|
| 255 | CenterY = round(median(y)); |
---|
| 256 | y = find(mask(:,CenterX)); |
---|
| 257 | if ~isempty(y) |
---|
| 258 | CenterY = round(median(y)); |
---|
| 259 | end |
---|
| 260 | % CenterX = round(median(x)); |
---|
| 261 | % CenterY = round(median(yt)); |
---|
| 262 | |
---|
| 263 | temp1 = sparse(1, 3*NuSupSize); |
---|
| 264 | temp2 = sparse(1, 3*NuSupSize); |
---|
| 265 | temp1(:,(Sup2Para(ClosestNList(i,1))*3-2): Sup2Para(ClosestNList(i,1))*3) = Ray(:,CenterY,CenterX)'; |
---|
| 266 | temp2(:,(Sup2Para(ClosestNList(i,2))*3-2): Sup2Para(ClosestNList(i,2))*3) = Ray(:,CenterY,CenterX)'; |
---|
| 267 | %NuNei-NuSTList |
---|
| 268 | % i |
---|
| 269 | if i < NuNei-NuSTList % only immediate connecting superpixels are weighted using MultiScaleSup |
---|
| 270 | % wei = WeiV*10;%*(MultiScaleSupTable(Sup2Para(ClosestNList(i,1)),2:end) == MultiScaleSupTable(Sup2Para(ClosestNList(i,2)),2:end))'; |
---|
| 271 | if MultiScaleFlag |
---|
| 272 | vector = (MultiScaleSupTable(Sup2Para(ClosestNList(i,1)),2:end) == MultiScaleSupTable(Sup2Para(ClosestNList(i,2)),2:end)); |
---|
| 273 | expV = exp(-10*(WeiV*vector' + ShiftCoP) ); |
---|
| 274 | wei = 1/(1+expV); |
---|
| 275 | else |
---|
| 276 | wei = 1; |
---|
| 277 | end |
---|
| 278 | else |
---|
| 279 | wei = 1; |
---|
| 280 | end |
---|
| 281 | CoPM1 = [CoPM1; temp1*wei]; |
---|
| 282 | CoPM2 = [CoPM2; temp2*wei]; |
---|
| 283 | tempWeiCoP = [SizeMaskAll]; |
---|
| 284 | CoPEstDepth = [CoPEstDepth; max(median(CleanedDepthMap(mask)),ClosestDist)]; |
---|
| 285 | |
---|
| 286 | mask = Sup == ClosestNList(i,2); |
---|
| 287 | SizeMaskAll = sum(mask(:)); |
---|
| 288 | [y x] = find(mask); |
---|
| 289 | CenterX = round(median(x)); |
---|
| 290 | CenterY = round(median(y)); |
---|
| 291 | y = find(mask(:,CenterX)); |
---|
| 292 | if ~isempty(y) |
---|
| 293 | CenterY = round(median(y)); |
---|
| 294 | end |
---|
| 295 | % [yt x] = find(mask); |
---|
| 296 | % CenterX = round(median(x)); |
---|
| 297 | % CenterY = round(median(yt)); |
---|
| 298 | |
---|
| 299 | temp1 = sparse(1, 3*NuSupSize); |
---|
| 300 | temp2 = sparse(1, 3*NuSupSize); |
---|
| 301 | temp1(:,(Sup2Para(ClosestNList(i,1))*3-2): Sup2Para(ClosestNList(i,1))*3) = Ray(:,CenterY,CenterX)'; |
---|
| 302 | temp2(:,(Sup2Para(ClosestNList(i,2))*3-2): Sup2Para(ClosestNList(i,2))*3) = Ray(:,CenterY,CenterX)'; |
---|
| 303 | CoPM1 = [CoPM1; temp1*wei]; |
---|
| 304 | CoPM2 = [CoPM2; temp2*wei]; |
---|
| 305 | tempWeiCoP = [tempWeiCoP; SizeMaskAll]; |
---|
| 306 | WeiCoP = [WeiCoP; tempWeiCoP]; |
---|
| 307 | CoPEstDepth = [CoPEstDepth; max(median(CleanedDepthMap(mask)),ClosestDist)]; |
---|
| 308 | % end |
---|
| 309 | end%========================================================================================================= |
---|
| 310 | |
---|
| 311 | % find the boundary point that might need to be stick ot each other========================================== |
---|
| 312 | HoriStickM_i = sparse(0,3*NuSupSize); |
---|
| 313 | HoriStickM_j = sparse(0,3*NuSupSize); |
---|
| 314 | HoriStickPointInd = []; |
---|
| 315 | EstDepHoriStick = []; |
---|
| 316 | for i = find(BounaryPHori==1)' |
---|
| 317 | j = i+VertYNuDepth; |
---|
| 318 | if Sup(i) == 0 || Sup(j) == 0 |
---|
| 319 | continue; |
---|
| 320 | end |
---|
| 321 | % if ~OccluList(sum( ClosestNList == repmat(sort([Sup(i) Sup(j)]), [NuNei 1]),2) == 2) |
---|
| 322 | % size(OccluList) |
---|
| 323 | % if ~any(sum( OccluList == repmat(sort([Sup(i) Sup(j)]), [size(OccluList,1) 1]),2) == 2) |
---|
| 324 | Target(1) = Sup2Para(Sup(i)); |
---|
| 325 | Target(2) = Sup2Para(Sup(j)); |
---|
| 326 | rayBoundary(:,1) = Ray(:,i); |
---|
| 327 | rayBoundary(:,2) = Ray(:,i); |
---|
| 328 | % betaTemp = StickHori;%*(DistStickLengthNormWei.^2)*beta(Target(I)); |
---|
| 329 | % betaTemp = StickHori*WeiV;%*(MultiScaleSupTable(Sup2Para(Sup(i)),2:end) == MultiScaleSupTable(Sup2Para(Sup(j)),2:end))';%*(DistStickLengthNormWei.^2)*beta(Target(I)); |
---|
| 330 | if MultiScaleFlag |
---|
| 331 | vector = (MultiScaleSupTable(Sup2Para(Sup(i)),2:end) == MultiScaleSupTable(Sup2Para(Sup(j)),2:end)); |
---|
| 332 | expV = exp(-10*(WeiV*vector' + ShiftStick) ); |
---|
| 333 | betaTemp = StickHori*(0.5+1/(1+expV)); %*(DistStickLengthNormWei.^2)*beta(Target(I)); |
---|
| 334 | % therr should always be sticking (know don't care about occlusion) |
---|
| 335 | else |
---|
| 336 | betaTemp = StickHori; |
---|
| 337 | end |
---|
| 338 | temp = sparse(3,NuSupSize); |
---|
| 339 | temp(:,Target(1)) = rayBoundary(:,1); |
---|
| 340 | HoriStickM_i = [HoriStickM_i; betaTemp*temp(:)']; |
---|
| 341 | temp = sparse(3,NuSupSize); |
---|
| 342 | temp(:,Target(2)) = rayBoundary(:,2); |
---|
| 343 | HoriStickM_j = [HoriStickM_j; betaTemp*temp(:)']; |
---|
| 344 | EstDepHoriStick = [EstDepHoriStick; sqrt(max(CleanedDepthMap(i),ClosestDist)*max(CleanedDepthMap(j),ClosestDist))]; |
---|
| 345 | HoriStickPointInd = [HoriStickPointInd i ]; |
---|
| 346 | % else |
---|
| 347 | % disp('Occlu'); |
---|
| 348 | % end |
---|
| 349 | end |
---|
| 350 | VertStickM_i = sparse(0,3*NuSupSize); |
---|
| 351 | VertStickM_j = sparse(0,3*NuSupSize); |
---|
| 352 | VertStickPointInd = []; |
---|
| 353 | EstDepVertStick = []; |
---|
| 354 | for i = find(BounaryPVert==1)' |
---|
| 355 | j = i+1; |
---|
| 356 | if Sup(i) == 0 || Sup(j) == 0 |
---|
| 357 | continue; |
---|
| 358 | end |
---|
| 359 | % if ~OccluList(sum( ClosestNList == repmat(sort([Sup(i) Sup(j)]), [NuNei 1]),2) == 2) |
---|
| 360 | % if ~any(sum( OccluList == repmat(sort([Sup(i) Sup(j)]), [size(OccluList,1) 1]),2) == 2) |
---|
| 361 | Target(1) = Sup2Para(Sup(i)); |
---|
| 362 | Target(2) = Sup2Para(Sup(j)); |
---|
| 363 | rayBoundary(:,1) = Ray(:,i); |
---|
| 364 | rayBoundary(:,2) = Ray(:,i); |
---|
| 365 | % betaTemp = StickVert;%*(DistStickLengthNormWei.^2)*beta(Target(I)); |
---|
| 366 | if MultiScaleFlag |
---|
| 367 | vector = (MultiScaleSupTable(Sup2Para(Sup(i)),2:end) == MultiScaleSupTable(Sup2Para(Sup(j)),2:end)); |
---|
| 368 | expV = exp(-10*(WeiV*vector' + ShiftStick) ); |
---|
| 369 | betaTemp = StickVert*(0.5+1/(1+expV)); |
---|
| 370 | % therr should always be sticking (know don't care about occlusion) |
---|
| 371 | else |
---|
| 372 | betaTemp = StickVert; |
---|
| 373 | end |
---|
| 374 | temp = sparse(3,NuSupSize); |
---|
| 375 | temp(:,Target(1)) = rayBoundary(:,1); |
---|
| 376 | VertStickM_i = [VertStickM_i; betaTemp*temp(:)']; |
---|
| 377 | temp = sparse(3,NuSupSize); |
---|
| 378 | temp(:,Target(2)) = rayBoundary(:,2); |
---|
| 379 | VertStickM_j = [VertStickM_j; betaTemp*temp(:)']; |
---|
| 380 | EstDepVertStick = [EstDepVertStick; sqrt(max(CleanedDepthMap(i),ClosestDist)*max(CleanedDepthMap(j),ClosestDist))]; |
---|
| 381 | VertStickPointInd = [VertStickPointInd i ]; |
---|
| 382 | % else |
---|
| 383 | % disp('Occlu'); |
---|
| 384 | % end |
---|
| 385 | end |
---|
| 386 | % ======================================Finish building up matrix=====================hard work====================== |
---|
| 387 | % ================================================================================================================ |
---|
| 388 | depthfile = strrep(filename{k},'img','depth_learned'); % |
---|
| 389 | |
---|
| 390 | for i = step |
---|
| 391 | opt = sdpsettings('solver','sedumi'); |
---|
| 392 | ParaPPCP = sdpvar(3*NuSupSize,1); |
---|
| 393 | F = set(ParaPPCP(3*(1:NuSupSize)-1).*YPointer<=0) ... |
---|
| 394 | +set(RayAllM*ParaPPCP >=1/FarestDist)... |
---|
| 395 | +set(RayAllM*ParaPPCP <=1/ClosestDist)... |
---|
| 396 | +set(RayAllOriM*ParaPPCP >=1/FarestDist)... |
---|
| 397 | +set(RayAllOriM*ParaPPCP <=1/ClosestDist); |
---|
| 398 | % First fit the plane to find the estimated plane parameters |
---|
| 399 | % If PosiM contain NaN data the corresponding Plane Parameter will also be NaN |
---|
| 400 | if i == 1 % W/O priors |
---|
| 401 | % solvesdp(F,norm( (PosiM*ParaPPCP-ones(size(PosiM,1),1))./(abs(VarM)),1)... |
---|
| 402 | solvesdp(F,norm( (PosiM*ParaPPCP-ones(size(PosiM,1),1))./exp(abs(VarM)/BandWith),1)... |
---|
| 403 | , opt); |
---|
| 404 | elseif i ==2 % Coner |
---|
| 405 | solvesdp(F,norm( (PosiM*ParaPPCP-ones(size(PosiM,1),1))./exp(abs(VarM)/BandWith),1)... |
---|
| 406 | +Center*norm(((CoPM1 - CoPM2)*ParaPPCP).*CoPEstDepth, 1)... |
---|
| 407 | , opt); |
---|
| 408 | else % all the priors |
---|
| 409 | % sum(sum(isnan(PosiM))) |
---|
| 410 | disp('Whole MRF') |
---|
| 411 | solvesdp(F,norm( (PosiM*ParaPPCP-ones(size(PosiM,1),1))./exp(abs(VarM)/BandWith),1)... |
---|
| 412 | +Center*norm(((CoPM1 - CoPM2)*ParaPPCP).*CoPEstDepth, 1)... |
---|
| 413 | +norm(((HoriStickM_i-HoriStickM_j)*ParaPPCP).*EstDepHoriStick,1)... |
---|
| 414 | +norm(((VertStickM_i-VertStickM_j)*ParaPPCP).*EstDepVertStick,1)... |
---|
| 415 | , opt); |
---|
| 416 | end |
---|
| 417 | |
---|
| 418 | %pause |
---|
| 419 | ParaPPCP = double(ParaPPCP); |
---|
| 420 | sum(isnan(ParaPPCP)) |
---|
| 421 | %pause |
---|
| 422 | SepPointMeasureHori = 1./(HoriStickM_i*ParaPPCP)-1./(HoriStickM_j*ParaPPCP); |
---|
| 423 | SepPointMeasureVert = 1./(VertStickM_i*ParaPPCP)-1./(VertStickM_j*ParaPPCP); |
---|
| 424 | ParaPPCP = reshape(ParaPPCP,3,[]); |
---|
| 425 | %any(any(isnan(ParaPPCP))) |
---|
| 426 | % porject the ray on planes to generate the ProjDepth |
---|
| 427 | FitDepthPPCP = FarestDist*ones(1,VertYNuDepth*HoriXNuDepth); |
---|
| 428 | FitDepthPPCP(~maskSkyEroded) = (1./sum(ParaPPCP(:,Sup2Para(SupEpand(~maskSkyEroded))).*Ray(:,~maskSkyEroded),1))'; |
---|
| 429 | FitDepthPPCP = reshape(FitDepthPPCP,VertYNuDepth,[]); |
---|
| 430 | % storage for comparison purpose ====================== |
---|
| 431 | depthMap = FarestDist*ones(1,VertYNuDepth*HoriXNuDepth); |
---|
| 432 | depthMap(~maskSkyEroded) = (1./sum(ParaPPCP(:,Sup2Para(SupEpand(~maskSkyEroded))).*RayOri(:,~maskSkyEroded),1))'; |
---|
| 433 | depthMap = reshape(depthMap,VertYNuDepth,[]); |
---|
| 434 | % if baseline == 0 |
---|
| 435 | % system(['mkdir ' ScratchDataFolder '/' Name{i} '_' DepthFolder '/']); |
---|
| 436 | % save([ScratchDataFolder '/' Name{i} '_' DepthFolder '/' depthfile '.mat'],'depthMap'); |
---|
| 437 | % elseif baseline == 1 |
---|
| 438 | % system(['mkdir ' ScratchDataFolder '/' Name{i} '_' DepthFolder '_baseline/']); |
---|
| 439 | % save([ScratchDataFolder '/' Name{i} '_' DepthFolder '_baseline/' depthfile '.mat'],'depthMap'); |
---|
| 440 | % else |
---|
| 441 | % system(['mkdir ' ScratchDataFolder '/' Name{i} '_' DepthFolder '_baseline2/']); |
---|
| 442 | % save([ScratchDataFolder '/' Name{i} '_' DepthFolder '_baseline2/' depthfile '.mat'],'depthMap'); |
---|
| 443 | % end |
---|
| 444 | |
---|
| 445 | % ===================================================== |
---|
| 446 | [Position3DFitedPPCP] = im_cr2w_cr(FitDepthPPCP,permute(Ray,[2 3 1])); |
---|
| 447 | if false %RenderVrmlFlag && i==3 |
---|
| 448 | [VrmlName] = vrml_test_faceset_goodSkyBoundary( filename{k}, Position3DFitedPPCP, FitDepthPPCP, permute(Ray,[2 3 1]), [Name{i} DepthFolder 'ExpVar'], ... |
---|
| 449 | [], [], 0, maskSkyEroded, maskG, 1, 0, a_default, b_default, Ox_default, Oy_default); |
---|
| 450 | system(['gzip -9 -c ' ScratchDataFolder '/vrml/' VrmlName ' > ' ScratchDataFolder '/vrml/' VrmlName '.gz']); |
---|
| 451 | %delete([ScratchDataFolder '/vrml/' VrmlName]); |
---|
| 452 | end |
---|
| 453 | % save([ScratchDataFolder '/data/PreDepth/FitDepthPPCP' num2str(k) '.mat'],'FitDepthPPCP','SepPointMeasureHori',... |
---|
| 454 | % 'SepPointMeasureVert','VertStickPointInd','HoriStickPointInd'); |
---|
| 455 | %save([ScratchDataFolder '/data/All.mat']); |
---|
| 456 | end |
---|
| 457 | |
---|
| 458 | %return; |
---|
| 459 | %==================Finished for one step MRF========================================================================================================== |
---|
| 460 | |
---|
| 461 | |
---|
| 462 | % ====================================following are 2nd and 3rd step MRF to give more visually pleasing result======================================= |
---|
| 463 | if ConerImprove |
---|
| 464 | [CoPM1 CoPM2 WeiCoP Sup Sup2ParaNew FixPara NuSup VertStickPointInd HoriStickPointInd] ... |
---|
| 465 | = FixParaFreeCorner( Sup, Ray, HoriStickPointInd, VertStickPointInd,... |
---|
| 466 | SepPointMeasureHori, SepPointMeasureVert, OccluList, WeiV); |
---|
| 467 | |
---|
| 468 | RayAllM = sparse(0,0); |
---|
| 469 | YPointer = []; |
---|
| 470 | for i = NuSup |
---|
| 471 | mask = Sup ==i; |
---|
| 472 | RayAllM = blkdiag( RayAllM, Ray(:,mask)'); |
---|
| 473 | [yt x] = find(mask); |
---|
| 474 | CenterX = round(median(x)); |
---|
| 475 | CenterY = round(median(yt)); |
---|
| 476 | YPointer = [YPointer; CenterY >= ceiling]; |
---|
| 477 | end |
---|
| 478 | |
---|
| 479 | NuSupSize = size(NuSup,2); |
---|
| 480 | opt = sdpsettings('solver','sedumi'); |
---|
| 481 | ParaCorner = sdpvar(3*NuSupSize,1); |
---|
| 482 | F = set(ParaCorner(3*(1:NuSupSize)-1).*YPointer<=0) + set(RayAllM*ParaCorner >=1/FarestDist)... |
---|
| 483 | +set(RayAllM*ParaCorner <=1/ClosestDist)... |
---|
| 484 | +set(ParaCorner([(Sup2ParaNew(FixPara)*3-2); (Sup2ParaNew(FixPara)*3-1); (Sup2ParaNew(FixPara)*3)]) == ... |
---|
| 485 | ParaPPCP([(Sup2Para(FixPara)*3-2); (Sup2Para(FixPara)*3-1); (Sup2Para(FixPara)*3)])); |
---|
| 486 | solvesdp(F,Center*norm((CoPM1 - CoPM2)*ParaCorner, 1)... |
---|
| 487 | , opt); |
---|
| 488 | ParaCorner = double(ParaCorner); |
---|
| 489 | ParaCorner = reshape(ParaCorner,3,[]); |
---|
| 490 | %FitDepthCorner = 80*ones(1,VertYNuDepth*HoriXNuDepth); |
---|
| 491 | FitDepthCorner = FarestDist*ones(1,VertYNuDepth*HoriXNuDepth); |
---|
| 492 | FitDepthCorner(Sup~=0) = (1./sum(ParaCorner(:,Sup2ParaNew(Sup(Sup~=0))).*Ray(:,Sup~=0),1))'; |
---|
| 493 | FitDepthCorner = reshape(FitDepthCorner,VertYNuDepth,[]); |
---|
| 494 | [Position3DFitedCorner] = im_cr2w_cr(FitDepthCorner,permute(Ray,[2 3 1])); |
---|
| 495 | [VrmlName] = vrml_test_faceset_goodSkyBoundary( filename{k}, Position3DFitedCorner, FitDepthCorner, permute(Ray,[2 3 1]), [Name 'Corner'], ... |
---|
| 496 | [], [], 0, maskSkyEroded, maskG, 1, 0, a_default, b_default, Ox_default, Oy_default); |
---|
| 497 | system(['gzip -9 -c ' ScratchDataFolder '/vrml/' VrmlName ' > ' ScratchDataFolder '/vrml/' VrmlName '.gz']); |
---|
| 498 | %delete([ScratchDataFolder '/vrml/' VrmlName]); |
---|
| 499 | end |
---|
| 500 | %return; |
---|
| 501 | % ============================================================ |
---|
| 502 | ;% generating new PosiMPPCP using the new position |
---|
| 503 | |
---|
| 504 | save([ScratchDataFolder '/data/VertShiftVrml.mat' ] ); |
---|
| 505 | % groundThreshold = cos(5*pi/180); % make small range |
---|
| 506 | % verticalThreshold = cos(50*pi/180); |
---|
| 507 | normPara = norms(ParaPPCP); |
---|
| 508 | normalizedPara = ParaPPCP ./ repmat( normPara, [3 1]); |
---|
| 509 | groundPara = abs(normalizedPara(2,:)) >= groundThreshold(YPosition); |
---|
| 510 | groundParaInd = find(groundPara); |
---|
| 511 | verticalPara = abs(normalizedPara(2,:)) <= verticalThreshold(YPosition); % change to have different range of vertical thre ================ |
---|
| 512 | verticalParaInd = find(verticalPara); |
---|
| 513 | |
---|
| 514 | if false |
---|
| 515 | [CoPM1 CoPM2 HoriStickM_i HoriStickM_j VertStickM_i VertStickM_j WeiCoP NewSup NewNuSup NewNuSupSize NewSup2Para FixPara VertStickPointInd HoriStickPointInd] ... |
---|
| 516 | = FreeSupSharpCorner( Sup, Ray, OccluList, WeiV, ShiftStick, StickHori, StickVert,... |
---|
| 517 | MultiScaleSupTable, verticalParaInd, graoundParaInd, MultiScaleFlag); |
---|
| 518 | |
---|
| 519 | end |
---|
| 520 | |
---|
| 521 | indexVertical = find( verticalPara)*3-1; |
---|
| 522 | indexGroundX = find( groundPara)*3-2; |
---|
| 523 | indexGroundZ = find( groundPara)*3; |
---|
| 524 | |
---|
| 525 | PosiMPPCP = sparse(0,0); |
---|
| 526 | VarM2 = sparse(0,0); |
---|
| 527 | % VertVar = 10^(-2); |
---|
| 528 | for i = NuSup |
---|
| 529 | mask = Sup == i; |
---|
| 530 | if any(verticalParaInd == Sup2Para(i)) |
---|
| 531 | mask = logical(zeros(size(Sup))); |
---|
| 532 | end |
---|
| 533 | % PosiMPPCP = blkdiag(PosiMPPCP, Posi3D(:,mask)'); |
---|
| 534 | % VarM2 = [VarM2; VertVar*ones(sum(mask(:)), 1)]; |
---|
| 535 | % mask = logical(zeros(size(Sup))); |
---|
| 536 | % elseif any(groundParaInd == Sup2Para(i)) |
---|
| 537 | % PosiMPPCP = blkdiag(PosiMPPCP, Position3DFitedPPCP(:,mask)'); |
---|
| 538 | % VarM2 = [VarM2; VarMap(mask)]; |
---|
| 539 | % else |
---|
| 540 | PosiMPPCP = blkdiag(PosiMPPCP, Posi3D(:,mask)'); |
---|
| 541 | VarM2 = [VarM2; VarMap(mask)]; |
---|
| 542 | % end |
---|
| 543 | end |
---|
| 544 | |
---|
| 545 | opt = sdpsettings('solver','sedumi'); |
---|
| 546 | Para = sdpvar(3*NuSupSize,1); |
---|
| 547 | F = set(Para(3*(1:NuSupSize)-1).*YPointer<=0) + set(RayAllM*Para >=1/FarestDist)... |
---|
| 548 | +set(RayAllM*Para <=1/ClosestDist)... |
---|
| 549 | +set(Para(indexVertical) == 0); |
---|
| 550 | % +set(Para(indexGroundX) == 0)... |
---|
| 551 | % +set(Para(indexGroundZ) == 0); |
---|
| 552 | if FractionalDepthError |
---|
| 553 | size(PosiMPPCP) |
---|
| 554 | solvesdp(F,norm( ( PosiMPPCP*Para-ones(size(PosiMPPCP,1),1))./exp(abs(VarM2)/BandWith),1)... |
---|
| 555 | +Center*norm(((CoPM1 - CoPM2)*Para).*CoPEstDepth, 1)... |
---|
| 556 | +norm(((HoriStickM_i-HoriStickM_j)*Para).*EstDepHoriStick,1)... |
---|
| 557 | +norm(((VertStickM_i-VertStickM_j)*Para).*EstDepVertStick,1)... |
---|
| 558 | +10*norm( Para(indexGroundX)./ normPara( find(groundPara) )', 1) ... |
---|
| 559 | +10*norm( Para(indexGroundZ)./ normPara( find(groundPara) )', 1) ... |
---|
| 560 | , opt); |
---|
| 561 | % sqrt( max(CoPM1*ParaPPCP(:), 1/FarestDist).*max(CoPM2*ParaPPCP(:), 1/FarestDist)), 1)... |
---|
| 562 | % +norm(((VertStickM_i-VertStickM_j)*Para)./... |
---|
| 563 | % sqrt( max(VertStickM_i*ParaPPCP(:), 1/FarestDist).*max(VertStickM_j*ParaPPCP(:), 1/FarestDist)),1)... |
---|
| 564 | % pause; |
---|
| 565 | % +1000*norm( Para(indexVertical)./ normPara( find(verticalPara) )',1) ... |
---|
| 566 | % +norm(((HoriStickM_i-HoriStickM_j)*Para)./sqrt((VertStickM_i*ParaPPCP(:)).*(VertStickM_j*ParaPPCP(:))),1)... |
---|
| 567 | else |
---|
| 568 | solvesdp(F,norm( RayMd*Para - DepthInverseM,1)... |
---|
| 569 | +Center*norm((CoPM1 - CoPM2)*Para,1)... |
---|
| 570 | +norm((VertStickM_i-VertStickM_j)*Para,1)... |
---|
| 571 | +norm((HoriStickM_i-HoriStickM_j)*Para,1)... |
---|
| 572 | +1000*norm( Para(indexVertical)./ normPara( find(verticalPara) )',1) ... |
---|
| 573 | +1000*norm( Para(indexGroundX)./ normPara( find(groundPara) )', 1) ... |
---|
| 574 | +1000*norm( Para(indexGroundZ)./ normPara( find(groundPara) )', 1) ... |
---|
| 575 | , opt); |
---|
| 576 | % +0.01*norm( RayMtilt*ParaPPCP,1)... |
---|
| 577 | % +norm(spdiags(WeiCoP,0,size(CoPM1,1),size(CoPM1,1))*(CoPM1 - CoPM2)*ParaPPCP,1)... |
---|
| 578 | % +norm( (CoPM1 - CoPM2)*ParaPPCP,1 )... |
---|
| 579 | end |
---|
| 580 | Para = reshape(Para,3,[]); |
---|
| 581 | Para = double(Para); |
---|
| 582 | sum(sum(isnan(Para))) |
---|
| 583 | %pause |
---|
| 584 | % FitDepth = 80*ones(1,VertYNuDepth*HoriXNuDepth); |
---|
| 585 | FitDepth = FarestDist*ones(1,VertYNuDepth*HoriXNuDepth); |
---|
| 586 | FitDepth(~maskSkyEroded) = (1./sum(Para(:,Sup2Para(SupEpand(~maskSkyEroded))).*Ray(:,~maskSkyEroded),1))'; |
---|
| 587 | FitDepth = reshape(FitDepth,VertYNuDepth,[]); |
---|
| 588 | %sum(isnan(FitDepth(:))) |
---|
| 589 | [Position3DFited] = im_cr2w_cr(FitDepth,permute(Ray,[2 3 1])); |
---|
| 590 | Position3DFited(3,:) = -Position3DFited(3,:); |
---|
| 591 | Position3DFited = permute(Position3DFited,[2 3 1]); |
---|
| 592 | temp = ray(:,:,1:2)./repmat(ray(:,:,3),[1 1 2]); |
---|
| 593 | PositionTex = permute(temp./repmat(cat(3,a,b),[VertYNuDepth HoriXNuDepth 1])+repmat(cat(3,Ox,Oy),[VertYNuDepth HoriXNuDepth 1]),[3 1 2]); |
---|
| 594 | PositionTex = permute(PositionTex,[2 3 1]); |
---|
| 595 | WrlFacest(Position3DFited,ones(3,3,2),Sup,filename{1},'./') |
---|
| 596 | % [VrmlName] = vrml_test_faceset_goodSkyBoundary( filename{k}, Position3DFited, FitDepth, permute(Ray,[2 3 1]), 'VNonSupport_ExpVar_NewVert5_50', ... |
---|
| 597 | % [], [], 0, maskSkyEroded, maskG, 1, 0, a_default, b_default, Ox_default, Oy_default); |
---|
| 598 | % system(['gzip -9 -c ' ScratchDataFolder '/vrml/' VrmlName ' > ' ScratchDataFolder '/vrml/' VrmlName '.gz']); |
---|
| 599 | % system(['mv ' ScratchDataFolder '/vrml/' VrmlName '.gz ' ScratchDataFolder '/vrml/' VrmlName]); |
---|
| 600 | %delete([ScratchDataFolder '/vrml/' VrmlName]); |
---|
| 601 | |
---|
| 602 | % save depth Map ++++++++++++++++ |
---|
| 603 | % depthMap = FitDepth; |
---|
| 604 | % system(['mkdir ' ScratchDataFolder '/VNonSupport_' DepthFolder '/']); |
---|
| 605 | % save([ScratchDataFolder '/VNonSupport_' DepthFolder '/' depthfile '.mat'],'depthMap'); |
---|
| 606 | % ============================= |
---|
| 607 | %end |
---|
| 608 | toc |
---|
| 609 | return; |
---|