tp@0: function EDB1plotmodel(eddatafile,plotoptions,plotoptions2,plotoptions3)
tp@0: % EDB1plotmodel - Plots a model which is given in an eddatafile.
tp@0: %
tp@0: % Input parameters:
tp@0: %   eddatafile (optional)   The input file. If not specified, a file
tp@0: %                           opening window will be presented.
tp@0: %   plotoptions (optional)  An integer that can give extra options:
tp@0: %                           if bit0 = 1 (= 1) => plot sources
tp@0: %                           if bit1 = 1 (= 2) => plot receivers
tp@0: %                           if bit2 = 1 (= 4) => plot plane normal vectors
tp@0: %                           if bit3 = 1 (= 8) => print plane numbers
tp@0: %                           if bit4 = 1 (=16) => print edge numbers
tp@0: %                           if bit5 = 1 (=32) => print corner numbers
tp@0: %                           if bit6 = 1 (=64) => print plane numbers using the
tp@0: %                                                CAD file numbering
tp@0: %                           if bit7 = 1 (=128)=> print corner numbers using the
tp@0: %                                                CAD file numbering
tp@0: %                           if bit8 = 1 (=256) => offedges are not included
tp@0: %                           if bit9 = 1 (=512) => only offedges are included
tp@0: %                           if plotoptions == -1 => plot edge by edge
tp@0: %                           if plotoptions == -2 => plot plane by plane
tp@0: %                           Example: the integer 11 = 8+2+1
tp@0: %                           so sources, receivers, and plane numbers will be plotted.
tp@0: %   plotoptions2 (optional) A matrix with two columns of 1 or many integers, the first
tp@0: %                           giving the source numbers to be plotted and the
tp@0: %                           second giving the receiver numbers to be
tp@0: %                           plotted. The first in the sequence will be
tp@0: %                           plotted with a black color.
tp@0: %                           Alternatively, if plotoptions = -1, then
tp@0: %                           plotoptions2 can be used to give a list of
tp@0: %                           edges to plot.
tp@0: %                           Alternatively, if plotoptions = -2, then
tp@0: %                           plotoptions2 can be used to give a list of
tp@0: %                           planes to plot.
tp@0: %   plotoptions3 (optional) A vector with two or three values, giving the
tp@0: %                           view position for the 3D plot view.
tp@0: %
tp@0: % Sources and receivers are taken from an sdatafile and an rdatafile, the file name of
tp@0: % which is assumed to be similar to the eddatafile.
tp@0: %
tp@0: % Uses functions EDB1strpend, EDB1myvalueinput
tp@0: %
tp@0: % ----------------------------------------------------------------------------------------------
tp@0: %   This file is part of the Edge Diffraction Toolbox by Peter Svensson.                       
tp@0: %                                                                                              
tp@0: %   The Edge Diffraction Toolbox is free software: you can redistribute it and/or modify       
tp@0: %   it under the terms of the GNU General Public License as published by the Free Software     
tp@0: %   Foundation, either version 3 of the License, or (at your option) any later version.        
tp@0: %                                                                                              
tp@0: %   The Edge Diffraction Toolbox is distributed in the hope that it will be useful,       
tp@0: %   but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS  
tp@0: %   FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.             
tp@0: %                                                                                              
tp@0: %   You should have received a copy of the GNU General Public License along with the           
tp@0: %   Edge Diffraction Toolbox. If not, see <http://www.gnu.org/licenses/>.                 
tp@0: % ----------------------------------------------------------------------------------------------
tp@0: % Peter Svensson (svensson@iet.ntnu.no) 20110704
tp@0: %
tp@0: % EDB1plotmodel(eddatafile,plotoptions,plotoptions2,plotoptions3);
tp@0: 
tp@0: %--------------------------------------------------------------
tp@0: % Read the eddatafile
tp@0: 
tp@0: if nargin == 0 | (nargin >= 1 & isstr(eddatafile) ~= 1)
tp@0: 	[eddatafile,eddatafilepath] = uigetfile('*eddata.mat','Please select the eddatafile');
tp@0:     [eddatafilepath,temp1,temp2] = fileparts(eddatafilepath);
tp@0: 	if ~isstr(eddatafile) | isempty(eddatafile)
tp@0: 		return
tp@0: 	end
tp@0: else
tp@0: 	[eddatafilepath,eddatafile,fileext] = fileparts(eddatafile);
tp@0:     eddatafile = [eddatafile,fileext];
tp@0: end
tp@0: Filestem = EDB1strpend(eddatafile,'_eddata');
tp@0: 
tp@0: disp(['eddatafile is ',eddatafilepath,filesep,eddatafile])
tp@0: 
tp@0: if nargin == 1
tp@0:     if isstr(eddatafile) ~= 1
tp@0:         plotoptions = eddatafile;
tp@0:     else
tp@0:         plotoptions = 0;    
tp@0:     end
tp@0:     plotoptions2 = [];
tp@0:     plotoptions3 = [];
tp@0: elseif nargin == 0
tp@0:     plotoptions = 0;    
tp@0:     plotoptions2 = [];
tp@0:     plotoptions3 = [];
tp@0: elseif nargin == 2
tp@0:     plotoptions2 = [];
tp@0:     plotoptions3 = [];
tp@0: elseif nargin == 3
tp@0:     nplotoptions2 = size(plotoptions2);
tp@0:     plotoptions3 = [];
tp@0: else
tp@0:     nplotoptions2 = size(plotoptions2);    
tp@0: end
tp@0: 
tp@0: plotedgesequence = 0; 
tp@0: plotplanesequence = 0;
tp@0: if plotoptions < 0    
tp@0:     if plotoptions == -1
tp@0:         plotedgesequence = 1;
tp@0:     elseif plotoptions == -2 
tp@0:         plotplanesequence = 1;
tp@0:     end
tp@0:     plotoptions = 0;
tp@0: end
tp@0: plotsources = bitget(plotoptions,1);
tp@0: plotreceivers = bitget(plotoptions,2);
tp@0: plotnvecs = bitget(plotoptions,3);
tp@0: plotplnumbers = bitget(plotoptions,4);
tp@0: plotednumbers = bitget(plotoptions,5);
tp@0: plotconumbers = bitget(plotoptions,6);
tp@0: plotplCADnumbers = bitget(plotoptions,7);
tp@0: plotcoCADnumbers = bitget(plotoptions,8);
tp@0: plotnooffedges = bitget(plotoptions,9);
tp@0: plotonlyoffedges = bitget(plotoptions,10);
tp@0: 
tp@0: if plotplCADnumbers, plotplnumbers = 1; end
tp@0: if plotcoCADnumbers, plotconumbers = 1; end
tp@0: 
tp@0: %--------------------------------------------------------------------------
tp@0: % Load the needed input files
tp@0: 
tp@0: if ~isempty(eddatafilepath)
tp@0:     eval(['load ',eddatafilepath,filesep,eddatafile])
tp@0: else
tp@0:     eval(['load ',eddatafile])    
tp@0: end
tp@0:     
tp@0: ncornersperplanevec = double(ncornersperplanevec);
tp@0: if plotsources
tp@0:     sdatafile = [eddatafilepath,filesep,Filestem,'_sdata_1.mat'];
tp@0:     if exist(sdatafile) == 2
tp@0:         allsources = [];
tp@0:         sfilecounter = 1;
tp@0:         while exist(sdatafile) == 2
tp@0:             eval(['load ',sdatafile])
tp@0:             allsources = [allsources;sources];
tp@0:             sfilecounter = sfilecounter +1;
tp@0:             sdatafile = [eddatafilepath,filesep,Filestem,'_sdata_',int2str(sfilecounter),'.mat']            
tp@0:         end 
tp@0:         sources = allsources;
tp@0:     else
tp@0:         sdatafile = [eddatafilepath,filesep,Filestem,'_sdata.mat'];
tp@0:         if exist(sdatafile) == 2
tp@0:             eval(['load ',sdatafile])
tp@0:         else
tp@0:             error(['ERROR: The sdata file named ',sdatafile,' could not be opened'])    
tp@0:         end
tp@0:     end
tp@0: end
tp@0: if plotreceivers
tp@0:     rdatafile = [eddatafilepath,filesep,Filestem,'_rdata_1.mat'];
tp@0:     if exist(rdatafile) == 2
tp@0:         allreceivers = [];
tp@0:         rfilecounter = 1;
tp@0:         while exist(rdatafile) == 2
tp@0:             eval(['load ',rdatafile])
tp@0:             allreceivers = [allreceivers;receivers];
tp@0:             rfilecounter = rfilecounter +1;
tp@0:             rdatafile = [eddatafilepath,filesep,Filestem,'_rdata_',int2str(rfilecounter),'.mat']            
tp@0:         end 
tp@0:         receivers = allreceivers;
tp@0:     else
tp@0:         rdatafile = [eddatafilepath,filesep,Filestem,'_rdata.mat'];
tp@0:         if exist(rdatafile) == 2
tp@0:             eval(['load ',rdatafile])
tp@0:         else
tp@0:             error(['ERROR: The rdata file named ',rdatafile,' could not be opened'])    
tp@0:         end
tp@0:     end    
tp@0: end
tp@0: 
tp@0: if plotplCADnumbers | plotcoCADnumbers
tp@0:     cadgeofile = [eddatafilepath,Filestem,'_cadgeo.mat'];
tp@0:     if exist(cadgeofile) == 2
tp@0:         eval(['load ',cadgeofile])
tp@0:     else
tp@0:         error(['ERROR: The cadgeo file named ',cadgeofile,' could not be opened'])            
tp@0:     end
tp@0: end
tp@0: 
tp@0: %--------------------------------------------------------------
tp@0: 
tp@0: [ncorners,slask] = size(corners);
tp@0: [nedges,slask] = size(edgecorners);
tp@0: [nplanes,slask] = size(planenvecs);
tp@0: 
tp@0: planelist = 1:nplanes;
tp@0: 
tp@0: if isempty(plotoptions3)
tp@0:     viewpos = EDB1myvalueinput('From which point do you want to watch the model? (x y z, or az el, with spaces inbetween) ',-1);
tp@0: else
tp@0:     viewpos = plotoptions3;
tp@0: end
tp@0: 
tp@0: hold off
tp@0: for ii = 1:nedges
tp@0: 	co1 = edgecorners(ii,1);
tp@0: 	co2 = edgecorners(ii,2);
tp@0: 	iv = [co1;co2];
tp@0: 	plot3(corners(iv,1),corners(iv,2),corners(iv,3))
tp@0:     if ii ==1
tp@0:         view(viewpos)
tp@0:         hold
tp@0:     end
tp@0: end
tp@0: 
tp@0: if plotedgesequence
tp@0:     for ii = 1:nedges
tp@0:     	co1 = edgecorners(ii,1);
tp@0:         co2 = edgecorners(ii,2);
tp@0:         iv = [co1;co2];
tp@0:         h = plot3(corners(iv,1),corners(iv,2),corners(iv,3));
tp@0:         set(h,'LineWidth',2)
tp@0:         disp(['Edge no. ',int2str(ii),' Closed wedge angle: ',num2str(closwedangvec(ii)*180/pi)])
tp@0:         pause
tp@0:         set(h,'LineWidth',1)    
tp@0:     
tp@0:     end
tp@0: end
tp@0: if plotplanesequence
tp@0:     if ~isempty(plotoptions2)
tp@0:        listofplanes = plotoptions2; 
tp@0:     else
tp@0:         listofplanes = [1:nplanes];
tp@0:     end
tp@0:     for jj = listofplanes
tp@0:         edgelist = edgesatplane(jj,:);
tp@0:         edgelist = edgelist(find(edgelist));
tp@0:         for ii = edgelist
tp@0:             co1 = edgecorners(ii,1);
tp@0:             co2 = edgecorners(ii,2);
tp@0:             iv = [co1;co2];
tp@0:             h = plot3(corners(iv,1),corners(iv,2),corners(iv,3));
tp@0:             set(h,'LineWidth',2)
tp@0:         end
tp@0:         disp(['Plane no. ',int2str(jj)])
tp@0:         pause
tp@0:         
tp@0:         % For some reason it doesnt work to replot the same list of edges
tp@0:         % with a thinner linewidth; even trying out with replotting the
tp@0:         % edges with another color just replots one of the edges??
tp@0:         % Therefore the inelegant method of replotting the whole model.
tp@0:         hold off
tp@0:         for ii = 1:nedges
tp@0:             co1 = edgecorners(ii,1);
tp@0:             co2 = edgecorners(ii,2);
tp@0:             iv = [co1;co2];
tp@0:             plot3(corners(iv,1),corners(iv,2),corners(iv,3))
tp@0:             if ii ==1
tp@0:                 view(viewpos)
tp@0:                 hold
tp@0:             end
tp@0:         end
tp@0:         
tp@0:         
tp@0:         
tp@0:     end
tp@0: end
tp@0: 
tp@0: if plotnvecs
tp@0: 	for ii = 1:nplanes
tp@0:         midpoint = mean(corners(planecorners(ii,1:(ncornersperplanevec(ii))),:));
tp@0:         
tp@0:         endpoint = midpoint + planenvecs(ii,:);
tp@0:         bothpoints = [midpoint;endpoint];
tp@0:         plot3(bothpoints(:,1),bothpoints(:,2),bothpoints(:,3));
tp@0:         plot3(midpoint(1),midpoint(2),midpoint(3),'ro');
tp@0: 	end
tp@0: end
tp@0: 
tp@0: if plotplnumbers
tp@0: 	for ii = 1:nplanes
tp@0:         midpoint = mean(corners(planecorners(ii,1:ncornersperplanevec(ii)),:));
tp@0:         endpoint = midpoint + planenvecs(ii,:)*0.1;
tp@0:         if plotplCADnumbers == 1
tp@0:             text(endpoint(1),endpoint(2),endpoint(3),int2str(planenumbers(ii)));
tp@0:         else
tp@0:             text(endpoint(1),endpoint(2),endpoint(3),int2str(ii));
tp@0:         end
tp@0:     
tp@0: 	end
tp@0: end
tp@0: 
tp@0: if plotednumbers
tp@0:     iv = [1:nedges];
tp@0:     if plotnooffedges, iv(offedges) = []; end
tp@0:     if plotonlyoffedges, iv = offedges; end
tp@0:     iv = iv(:).';
tp@0:     
tp@0:     for ii = iv
tp@0:         midpoint = mean(corners(edgecorners(ii,1:2),:));
tp@0:         text(midpoint(1),midpoint(2),midpoint(3),int2str(ii));
tp@0:     end
tp@0: end
tp@0: 
tp@0: if plotconumbers
tp@0:     for ii = 1:ncorners
tp@0:         if plotcoCADnumbers == 1
tp@0:             text(corners(ii,1),corners(ii,2),corners(ii,3),int2str(cornernumbers(ii)))    
tp@0:         else
tp@0:             text(corners(ii,1),corners(ii,2),corners(ii,3),int2str(ii))    
tp@0:         end
tp@0:     end
tp@0:     
tp@0: end
tp@0: 
tp@0: if plotsources == 1
tp@0:     if isempty(plotoptions2)
tp@0:         plot3(sources(:,1),sources(:,2),sources(:,3),'*');
tp@0:     else
tp@0:         if nplotoptions2 == 1
tp@0:             plot3(sources(plotoptions2(1),1),sources(plotoptions2(1),2),sources(plotoptions2(1),3),'*')
tp@0:         else
tp@0:             plot3(sources(plotoptions2(1,1),1),sources(plotoptions2(1,1),2),sources(plotoptions2(1,1),3),'k*')
tp@0:             plot3(sources(plotoptions2(2:nplotoptions2,1),1),sources(plotoptions2(2:nplotoptions2,1),2),sources(plotoptions2(2:nplotoptions2,1),3),'*')            
tp@0:         end
tp@0:     end
tp@0: end
tp@0: 
tp@0: if plotreceivers == 1
tp@0:     if isempty(plotoptions2)
tp@0:         plot3(receivers(:,1),receivers(:,2),receivers(:,3),'ro');
tp@0:     else
tp@0:         if nplotoptions2 == 1
tp@0:             plot3(receivers(plotoptions2(2),1),receivers(plotoptions2(2),2),receivers(plotoptions2(2),3),'ro');
tp@0:         else
tp@0:             plot3(receivers(plotoptions2(1,2),1),receivers(plotoptions2(1,2),2),receivers(plotoptions2(1,2),3),'ko');
tp@0:             plot3(receivers(plotoptions2(2:nplotoptions2,2),1),receivers(plotoptions2(2:nplotoptions2,2),2),receivers(plotoptions2(2:nplotoptions2,2),3),'ro');
tp@0:             
tp@0:         end
tp@0:     end
tp@0: end