Mercurial > hg > human-echolocation
view private/EDB1findpathsISESx.m @ 18:2d5f50205527 jabuilder_int tip
Escape the trailing backslash as well
author | Chris Cannam |
---|---|
date | Tue, 30 Sep 2014 16:23:00 +0100 |
parents | 90220f7249fc |
children |
line wrap: on
line source
function edpathsfile = EDB1findpathsISESx(eddatafile,lengthNspecmatrix,lengthNdiffmatrix,... singlediffcol,startindicessinglediff,endindicessinglediff,... ndecimaldivider,PointertoIRcombs,IRoriginsfrom,S,R,isou,irec,... directsound,specorder,difforder,nedgesubs,visplanesfromS,visplanesfromR,... vispartedgesfromS,vispartedgesfromR,desiredfilename) % EDB1findpathsISESx - Finds all possible paths that include direct sound, specular, diffraction. % Finds all possible paths of types direct sound, specular, diffraction % combinations of specular and diffraction. % % Input parameters: % eddatafile Taken directly from the setup-file or created automatically % lengthNspecmatrix, lengthNdiffmatrix % singlediffcol, startindicessinglediff, endindicessinglediff % ndecimaldivider, PointertoIRcombs, IRoriginsfrom % Taken directly from the ISEStree-file % S, R, isou, irec Coordinates and counter number of source and receiver, taken % from the setup file. % directsound, specorder, difforder % Taken directly from the setup-file. % nedgesubs Taken directly from the setup-file. % visplanesfromS, visplanesfromR, vispartedgesfromS, vispartedgesfromR % Taken directly from the srdatafile. % desiredfilename (optional) The desired name of the output file. % ISCOORDS (global) The matrix ISCOORDS from the ISEStreefile % POTENTIALISES (global) The matrix POTENTIALISES from the ISEStreefile % ORIGINSFROM (global) % ISESVISIBILITY (global) % IVNDIFFMATRIX (global) The list IVNDIFFMATRIX from the ISEStreefile % IVNSPECMATRIX (global) % REFLORDER (global) % % See descriptions in EDB1findISEStree, EDB1srgeo and EDB1mainISES % % Output parameters: % edpathsfile The output file name, where the output data is % stored. If the optional input parameter desiredfilename, % is given, then edpathsfile will simple repeat % this name. If desiredfilename was not specified % a file name is created by extracting the file % stem from the input parameter eddatafile, and % with '_edpaths' added. % % Output data in the edpathsfile: % pathtypevec A matrix, [ncombs,specorder], describing the % type of reflection, using 'f' for the direct % sound, 's' for specular reflection and 'd' for % diffraction. % reflpaths A matrix, [ncombs,specorder], giving the plane % and edge numbers involved in each refl. comb. % specextradata A sparse matrix, [ncombs,(specorder+1)*3] % containing the coordinates of the image source % (col 1-3), the image receiver (col 4-6) % and the coordinates of all specular hit points % (only for purely specular combinations). % edgeextradata A sparse matrix, [ncombs,(specorder*2)], with % the visibility of the involved edge, denoted % by two values between 0 and 1. % S The coordinates of the source. % R The coordinates of the receiver. % mainlistguide A matrix, [nposscombs,3], which for each row gives % 1. the number of rows in "reflpaths" and % "pathtypevec" that have one type of combination. % 2. the first row and 3. the last row % mainlistguidepattern A matrix, [nposscombs,specorder] which for each % row gives the description of each possible combination % using 'f', 's' and 'd'. % directsoundrow The value 0 or 1 indicating whether the first row in % mainlistguide & mainlistguidepattern & reflpaths & pathtypevec % contains the direct sound or not. % allspecrows A list, [1,2], containing the first and last row numbers % of mainlistguide & mainlistguidepattern that contain specular % reflections. If there are no specular reflections, allspecrows % contains [0 0]. % firstdiffrow A number indicating which is the first row of mainlistguide & % mainlistguidepattern that contain a diffraction component. If % there are no diffraction components, firstdiffrow = 0. % Sinsideplanenumber The plane number that the source is directly at. If the % source is not placed directly at a plane, then % souinsideplanenumber = []; % Rinsideplanenumber The plane number that the receiver is directly at. If the % receiver is not placed directly at a plane, then % recinsideplanenumber = []; % % Uses the functions EDB1strpend EDB1directsound % EDB1speculISES EDB1diffISESx EDB1diff2ISES EDB1diffNISES % % ---------------------------------------------------------------------------------------------- % This file is part of the Edge Diffraction Toolbox by Peter Svensson. % % The Edge Diffraction Toolbox is free software: you can redistribute it and/or modify % it under the terms of the GNU General Public License as published by the Free Software % Foundation, either version 3 of the License, or (at your option) any later version. % % The Edge Diffraction Toolbox is distributed in the hope that it will be useful, % but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS % FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details. % % You should have received a copy of the GNU General Public License along with the % Edge Diffraction Toolbox. If not, see <http://www.gnu.org/licenses/>. % ---------------------------------------------------------------------------------------------- % Peter Svensson (svensson@iet.ntnu.no) 20050404 % % edpathsfile = EDB1findpathsISES(eddatafile,... % lengthNspecmatrix,lengthNdiffmatrix,singlediffcol,startindicessinglediff,... % endindicessinglediff,ndecimaldivider,PointertoIRcombs,IRoriginsfrom,S,R,isou,irec,... % directsound,specorder,difforder,nedgesubs,visplanesfromS,visplanesfromR,... % vispartedgesfromS,vispartedgesfromR,desiredfilename); global SHOWTEXT global POTENTIALISES ISCOORDS IVNDIFFMATRIX global IVNSPECMATRIX ORIGINSFROM ISESVISIBILITY REFLORDER eval(['load ',eddatafile]) clear cornerinfrontofplane edgeseesplane planeseesplane if nargin < 22 [filepath,filestem,fileext] = fileparts(eddatafile); edpathsfile = [[filepath,filesep],EDB1strpend(filestem,'_eddata'),'_',int2str(isou),'_',int2str(irec),'_edpaths.mat']; else edpathsfile = desiredfilename; end nplanes = length(planeisthin); pathtypevec = []; reflpaths = []; specextradata = []; edgeextradata = []; mainlistguide = []; mainlistguidepattern = []; %------------------------------------------------------- % ############################# % ## DIRECT SOUND ## % ############################# if directsound == 1 if SHOWTEXT >= 2 disp(' Direct sound') end dirsoundok = EDB1directsound(eddatafile,S,R,visplanesfromS,visplanesfromR); if dirsoundok == 1 pathtypevec = [pathtypevec;['f',zeros(1,specorder-1)]]; reflpaths = [reflpaths;zeros(1,specorder)]; if specorder > 0 specextradata = [specextradata;[S R zeros(1,(specorder-1)*3)]]; else specextradata = [specextradata;[S R]]; end edgeextradata = [edgeextradata;zeros(1,specorder*2)]; mainlistguide = [1 1 1]; mainlistguidepattern = ['f' ' '*ones(1,specorder-1)]; directsoundrow = 1; else directsoundrow = 0; end else directsoundrow = 0; end %------------------------------------------------------- % ############################# % ## SPECULAR ## % ############################# if specorder >= 1 if SHOWTEXT >= 2 disp(' SPECULAR reflections') end [validISlist,validIScoords,allreflpoints,listguide,listofreflorder] = EDB1speculISES(eddatafile,... S,R,lengthNspecmatrix,specorder,visplanesfromR); [nvalidreflorders,slask] = size(listguide); for ii = 1:nvalidreflorders ncombs = listguide(ii,1); norder = listofreflorder(ii); if ncombs > 0 pathtypevec = [pathtypevec;['s'*ones(ncombs,norder) zeros(ncombs,specorder-norder)]]; reflpaths = [reflpaths;validISlist(listguide(ii,2):listguide(ii,3),:)]; if size(pathtypevec,2) > size(reflpaths,2) reflpaths = [reflpaths zeros(size(reflpaths,1),size(pathtypevec,2)-size(reflpaths,2))]; end specextradata = [specextradata;[validIScoords(listguide(ii,2):listguide(ii,3),1:3) allreflpoints(listguide(ii,2):listguide(ii,3),:)]]; if size(specextradata,2) < (size(pathtypevec,2)+1)*3 specextradata = [specextradata zeros(size(specextradata,1),(size(pathtypevec,2)+1)*3-size(specextradata,2))]; end edgeextradata = [edgeextradata;zeros(ncombs,specorder*2)]; end mainlistguidepattern = [mainlistguidepattern;['s'*ones(1,norder) ' '*ones(1,specorder-norder)]]; end [nprevious,slask] = size(mainlistguide); if nprevious > 0 listguide(:,2:3) = listguide(:,2:3)+mainlistguide(nprevious,3); end mainlistguide = [mainlistguide;listguide]; if nvalidreflorders == 0 allspecrows = [0 0]; else allspecrows = [min([1 nvalidreflorders]) nvalidreflorders]+directsoundrow; end firstdiffrow = nvalidreflorders + directsoundrow + 1; else allspecrows = [0 0]; firstdiffrow = directsoundrow + 1; end %------------------------------------------------------- % #################################### % ## SINGLE DIFFRACTION ## % #################################### if difforder >=1 & ~isempty(lengthNdiffmatrix) if SHOWTEXT >= 2 disp(' SINGLE DIFFRACTION combined with any order specular') end [edgedifflist,startandendpoints,prespeclist,postspeclist,validISEDcoords,validEDIRcoords,listguide,listoforders,... bigedgeweightlist] = EDB1diffISESx(eddatafile,S,R,... IVNDIFFMATRIX(1:lengthNdiffmatrix(1),1),singlediffcol,startindicessinglediff,endindicessinglediff,... specorder,visplanesfromR,vispartedgesfromS,vispartedgesfromR,nedgesubs,ndecimaldivider,PointertoIRcombs,IRoriginsfrom); [nvalidcombs,slask] = size(listguide); nprespecs = listoforders(:,2) - 1; npostspecs = listoforders(:,1) - nprespecs - 1; [nrows,nsubsegmentcols] = size(startandendpoints); nsubsegments = uint8( (startandendpoints(:,2:2:nsubsegmentcols))~=0 ); if nrows > 1 & nsubsegmentcols > 2 nsubsegments = sum(nsubsegments.').'; end nmaxsubsegments = max(nsubsegments); expandedlistguide = listguide; for ii = 1:nvalidcombs ncombs = listguide(ii,1); iv = [listguide(ii,2):listguide(ii,3)]; pathtypevec = [pathtypevec; ['s'*ones(ncombs,nprespecs(ii)) 'd'*ones(ncombs,1) 's'*ones(ncombs,npostspecs(ii)) zeros(ncombs,specorder-listoforders(ii)) ]]; reflpaths = [reflpaths; [prespeclist(iv,1:nprespecs(ii)) edgedifflist(iv) postspeclist(iv,1:npostspecs(ii)) zeros(ncombs,specorder-listoforders(ii)) ]]; specextradata = [specextradata;[ validISEDcoords(iv,1:3) validEDIRcoords(iv,1:3) zeros(ncombs,(specorder-1)*3) ]]; edgeextradata = [edgeextradata;[startandendpoints(iv,1:2) zeros(ncombs,(specorder-1)*2) ]]; mainlistguidepattern = [mainlistguidepattern;['s'*ones(1,nprespecs(ii)) 'd' 's'*ones(1,npostspecs(ii)) zeros(1,specorder-listoforders(ii))]]; if nmaxsubsegments > 2 for jj = 2:nmaxsubsegments ivsub = find(nsubsegments(iv)>=jj); ncombssubset = length(ivsub); if ncombssubset > 0 pathtypevec = [pathtypevec; ['s'*ones(ncombssubset,nprespecs(ii)) 'd'*ones(ncombssubset,1) 's'*ones(ncombssubset,npostspecs(ii)) zeros(ncombssubset,specorder-listoforders(ii)) ]]; reflpaths = [reflpaths; [prespeclist(iv(ivsub),1:nprespecs(ii)) edgedifflist(iv(ivsub)) postspeclist(iv(ivsub),1:npostspecs(ii)) zeros(ncombssubset,specorder-listoforders(ii)) ]]; specextradata = [specextradata;[ validISEDcoords(iv(ivsub),1:3) validEDIRcoords(iv(ivsub),1:3) zeros(ncombssubset,(specorder-1)*3) ]]; edgeextradata = [edgeextradata;[startandendpoints(iv(ivsub),(jj-1)*2+1:(jj-1)*2+2) zeros(ncombssubset,(specorder-1)*2) ]]; expandedlistguide(ii,1) = expandedlistguide(ii,1) + ncombssubset; expandedlistguide(ii,3) = expandedlistguide(ii,3) + ncombssubset; for kk = ii+1:nvalidcombs expandedlistguide(kk,2:3) = expandedlistguide(kk,2:3) + ncombssubset; end end end end end [n1,n2] = size(mainlistguide); if n1 > 0 listoffset = mainlistguide(n1,3); else listoffset = 0; end if size(expandedlistguide,1) > 0 mainlistguide = [mainlistguide;[expandedlistguide(:,1) expandedlistguide(:,2:3)+listoffset]]; end end %------------------------------------------------------- % ################################################ % ## DOUBLE AND HIGHER DIFFRACTION ## % ################################################ for kk = 2:min([specorder difforder]) if ~isempty(lengthNdiffmatrix) if SHOWTEXT >= 2 disp([' DIFFRACTION, ORDER ',int2str(kk),' combined with any order specular']) end nmaxsubsegments = 0; % We will use the spec-edge-spec combs that were found valid for % first-order diffraction, because if higher-order diffraction combs % end with the same edge-spec-spec-... comb we know the visibility % already. % First, pick out only the combs that have postspecs. Second, keep % only a unique set because there will be many repeated combinations. if kk == 2 iv = find(sum(postspeclist.')>0 ); edgedifflist = edgedifflist(iv,:); postspeclist = postspeclist(iv,:); bigedgeweightlist = bigedgeweightlist(iv,:); validEDIRcoords = validEDIRcoords(iv,:); patternthatisOK = [edgedifflist postspeclist]; [patternthatisOK,iv,slask] = unique(patternthatisOK,'rows'); edgedifflistin = edgedifflist(iv,:); postspeclistin = postspeclist(iv,:); bigedgeweightlistin = bigedgeweightlist(iv,:); validEDIRcoordsin = validEDIRcoords(iv,:); end if kk == 2 maxrownumber = max(lengthNdiffmatrix(1:2)); [edgedifflist,startandendpoints,prespeclist,midspeclist,postspeclist,validISEDcoords,validEDIRcoords,listguide,listofallspecs] = EDB1diff2ISES(eddatafile,S,R,... IVNDIFFMATRIX(1:maxrownumber,1:2),lengthNdiffmatrix(1:2),... specorder,visplanesfromR,vispartedgesfromS,vispartedgesfromR,nedgesubs,ndecimaldivider,edgedifflistin,postspeclistin,... bigedgeweightlistin,validEDIRcoordsin,edgeplaneperptoplane1,edgeplaneperptoplane2); [nrows,nsubsegmentcols] = size(startandendpoints); nsubsegments = uint8((startandendpoints(:,2:4:nsubsegmentcols))~=0); if nrows > 1 & nsubsegmentcols > 4 nsubsegments = sum(nsubsegments.').'; end nmaxsubsegments = max(nsubsegments); elseif kk >= 3 maxrownumber = max(lengthNdiffmatrix(1:kk)); [edgedifflist,startandendpoints,prespeclist,postspeclist,validISEDcoords,validEDIRcoords,listguide,listoforders] = EDB1diffNISES(eddatafile,S,R,... IVNDIFFMATRIX(1:maxrownumber,1:kk),lengthNdiffmatrix(1:kk),kk,... specorder,visplanesfromR,vispartedgesfromS,vispartedgesfromR,nedgesubs,ndecimaldivider,edgedifflistin,postspeclistin,... bigedgeweightlistin,validEDIRcoordsin); end [nvalidcombs,slask] = size(listguide); if kk >= 3 nprespecs = listoforders(:,2) - 1; npostspecs = listoforders(:,1) - nprespecs - kk; else nprespecs = listofallspecs(:,1); nmidspecs = listofallspecs(:,2); npostspecs = listofallspecs(:,3); listoforders = nprespecs+nmidspecs+npostspecs+2; end for ii = 1:nvalidcombs ncombs = listguide(ii,1); i1 = listguide(ii,2); i2 = listguide(ii,3); if kk == 2 pathtypevec = [pathtypevec; ['s'*ones(ncombs,nprespecs(ii)) 'd'*ones(ncombs,1) 's'*ones(ncombs,nmidspecs(ii)) 'd'*ones(ncombs,1) 's'*ones(ncombs,npostspecs(ii)) zeros(ncombs,specorder-listoforders(ii)) ]]; else pathtypevec = [pathtypevec; ['s'*ones(ncombs,nprespecs(ii)) 'd'*ones(ncombs,kk) 's'*ones(ncombs,npostspecs(ii)) zeros(ncombs,specorder-listoforders(ii)) ]]; end if kk == 2 reflpaths = [reflpaths; [prespeclist(i1:i2,1:nprespecs(ii)) edgedifflist(i1:i2,1) midspeclist(i1:i2,1:nmidspecs(ii)) edgedifflist(i1:i2,2) postspeclist(i1:i2,1:npostspecs(ii)) zeros(ncombs,specorder-listoforders(ii)) ]]; else reflpaths = [reflpaths; [prespeclist(i1:i2,1:nprespecs(ii)) edgedifflist(i1:i2,1:kk) postspeclist(i1:i2,1:npostspecs(ii)) zeros(ncombs,specorder-listoforders(ii)) ]]; end specextradata = [specextradata;[ validISEDcoords(i1:i2,1:3) validEDIRcoords(i1:i2,1:3) zeros(ncombs,(specorder-1)*3) ]]; edgeextradata = [edgeextradata;[startandendpoints(i1:i2,1:kk*2) zeros(ncombs,(specorder-kk)*2) ]]; if kk == 2 mainlistguidepattern = [mainlistguidepattern;['s'*ones(1,nprespecs(ii)) 'd' 's'*ones(1,nmidspecs(ii)) 'd' 's'*ones(1,npostspecs(ii)) zeros(1,specorder-listoforders(ii))]]; else mainlistguidepattern = [mainlistguidepattern;['s'*ones(1,nprespecs(ii)) 'd'*ones(1,kk) 's'*ones(1,npostspecs(ii)) zeros(1,specorder-listoforders(ii))]]; end for jj = 2:nmaxsubsegments ivsub = find(nsubsegments(iv)>=jj); ncombssubset = length(ivsub); if ncombssubset > 0 pathtypevec = [pathtypevec; ['s'*ones(ncombssubset,nprespecs(ii)) 'd'*ones(ncombssubset,1) 's'*ones(ncombssubset,nmidspecs(ii)) 'd'*ones(ncombssubset,1) 's'*ones(ncombssubset,npostspecs(ii)) zeros(ncombssubset,specorder-listoforders(ii)) ]]; reflpaths = [reflpaths; [prespeclist(iv(ivsub),1:nprespecs(ii)) edgedifflist(iv(ivsub),1) midspeclist(iv(ivsub),1:nmidspecs(ii)) edgedifflist(iv(ivsub),2) postspeclist(iv(ivsub),1:npostspecs(ii)) zeros(ncombssubset,specorder-listoforders(ii)) ]]; specextradata = [specextradata;[ validISEDcoords(iv(ivsub),1:3) validEDIRcoords(iv(ivsub),1:3) zeros(ncombssubset,(specorder-1)*3) ]]; edgeextradata = [edgeextradata;[startandendpoints(iv(ivsub),(jj-1)*4+1:(jj-1)*4+4) zeros(ncombssubset,(specorder-1)*4) ]]; expandedlistguide(ii,1) = expandedlistguide(ii,1) + ncombssubset; expandedlistguide(ii,3) = expandedlistguide(ii,3) + ncombssubset; for kk = ii+1:nvalidcombs expandedlistguide(kk,2:3) = expandedlistguide(kk,2:3) + ncombssubset; end end end end [n1,n2] = size(mainlistguide); if n1 > 0 listoffset = mainlistguide(n1,3); else listoffset = 0; end mainlistguide = [mainlistguide;[listguide(:,1) listguide(:,2:3)+listoffset]]; end end [n1,n2] = size(mainlistguide); if firstdiffrow > n1 firstdiffrow = 0; end %------------------------------------------------------- % Is the source or receiver directly at a plane? Sinsideplanenumber = find(visplanesfromS==4); Rinsideplanenumber = find(visplanesfromR==4); %------------------------------------------------------- % Get rid of empty columns if ~isempty(pathtypevec) checksum = sum(pathtypevec); ncols = length(checksum); while checksum(ncols) == 0 pathtypevec(:,ncols) = []; checksum = sum(pathtypevec); ncols = length(checksum); end end %------------------------------------------------------- % Save the output data specextradata = sparse(specextradata); edgeextradata = sparse(edgeextradata); pathtypevec = uint8(pathtypevec); [ncombs,slask] = size(reflpaths); if ncombs+1 < 256 mainlistguide = uint8(mainlistguide); elseif ncombs+1 < 65536 mainlistguide = uint16(mainlistguide); else mainlistguide = uint32(mainlistguide); end Varlist = [' pathtypevec reflpaths specextradata edgeextradata S R mainlistguide mainlistguidepattern directsoundrow allspecrows firstdiffrow Sinsideplanenumber Rinsideplanenumber']; eval(['save ',edpathsfile,Varlist])