theodisplay.m

基于matlab的反演程序,用于地球物理勘探中射线追踪及偏移成像程序.

	it=find(dat<0);
	dt=get(abs(dat(it)),'userdata');

	%determine the labeling interval
	dat=get(htimelabel,'userdata');
	it=find(dat<0);
	ievery=get(abs(dat(it)),'userdata');

	%first tick at the first even interval smaller than t(1)
	tnot=floor(t(1)/dt)*dt;
	if( abs((tnot-t(1)))>.2*dt )
		tnot=ceil(t(1)/dt)*dt;
	end

	ticks=tnot:dt:t(length(t));
	%generate tick labels
	lblmat=[];
	for k=1:length(ticks)
		if( rem(k-1,ievery)== 0)
			lblmat=strmat(lblmat,num2str(ticks(k)));
		else
			lblmat=strmat(lblmat,' ');
		end
	end
	ind=find(lblmat==1);
   if(~isempty(ind))
      lblmat(ind)=32;
   end	
	set(gca,'ytick',ticks,'yticklabel',lblmat);

	%turn on grid if needed
	if(strcmp(get(hgridon,'checked'),'on'))
		grid on;
	end
	
	%put the title at the bottom
	if(~isempty(title))
		xlim=get(gca,'xlim');
		ylim=get(gca,'ylim');
		xtit=mean(xlim);
		ytit=max(ylim)+.05*abs(diff(ylim));
		htit=text(xtit,ytit,title);
		ext=get(htit,'extent');
		xtit=xtit-.5*ext(3);
		set(htit,'position',[xtit,ytit]);
	end

	%yaxis label
	ylabel('Time in Seconds');
	
	return;
	
end

%
%	*******REPOSITION THE TOPS LABELS*********
%
if(strcmp(action,'repositops'))
	h=get(gcf,'userdata');
  	hstore1=h(12);
		hlocate=h(10);
	info=get(hstore1,'userdata');
	if(~isempty(info))
	 xlim=get(gca,'xlim');
	 ntops=info(1);
	 if(ntops)
			hline=info(2:ntops+1);
			htext=info(ntops+2:2*ntops+1);
			
			for k=1:ntops
				pos=get(htext(k),'position');
				set(htext(k),'position',[xlim(2) pos(2)]);
				set(hline(k),'xdata',xlim);
			end
		end
	end

	%if locate was on, turn it back on
	if( strcmp(get(hlocate,'checked'),'on'))
		set(hlocate,'checked','off');
		theodisplay('locate');
	end

	return;
end


%
% manipulate the grid
%
if(strcmp(action,'grid'))
	hmenu=gcbo;
	flag=get(hmenu,'checked');

	if(strcmp(flag,'on'))%turn it off
		set(hmenu,'checked','off');
		grid off
	else
		set(hmenu,'checked','on');
		grid on
	end

	return;
end

% change timing line interval
if(strcmp(action,'timeint'))
 h=get(gcf,'userdata');
 hmsg=h(1);
	hmenu=gcbo;
	htimeint=get(hmenu,'parent');
	htimelabel=h(9);

 dat=get(htimeint,'userdata');

 %turn off previous
 it=find(dat<0);
 dat(it)=-dat(it);
 set(dat(it),'checked','off');

 %turn on
 set(hmenu,'checked','on');
 it=find(dat==hmenu);
 dat(it)=-dat(it);

 set(htimeint,'userdata',dat);

	%determine the yticks
	% get the trace etc
	trcs=get(hmsg,'userdata');
	t=trcs{2};

	%determine the tick interval
	dat=get(htimeint,'userdata');
	it=find(dat<0);
	dt=get(abs(dat(it)),'userdata');

	%determine the labeling interval
	dat=get(htimelabel,'userdata');
	it=find(dat<0);
	ievery=get(abs(dat(it)),'userdata');

	%first tick at the first even interval smaller than t(1)
	tnot=floor(t(1)/dt)*dt;
	if( abs((tnot-t(1)))>.2*dt )
		tnot=ceil(t(1)/dt)*dt;
	end

	ticks=tnot:dt:t(length(t));
	%generate tick labels
	lblmat=[];
	for k=1:length(ticks)
		if( rem(k-1,ievery)== 0)
			lblmat=strmat(lblmat,num2str(ticks(k)));
		else
			lblmat=strmat(lblmat,' ');
		end
	end
	ind=find(lblmat==1);
   if(~isempty(ind))
      lblmat(ind)=32;
   end		
	set(gca,'ytick',ticks,'yticklabel',lblmat);


end

%change the labeling interval
if(strcmp(action,'label'))
	h=get(gcf,'userdata');
	hmsg=h(1);
	htimeint=h(5);
	hmenu=gcbo;
	htimelabel=get(hmenu,'parent');

 dat=get(htimelabel,'userdata');

 %turn off previous
 it=find(dat<0);
 dat(it)=-dat(it);
 set(dat(it),'checked','off');

 %turn on
 set(hmenu,'checked','on');
 it=find(dat==hmenu);
 dat(it)=-dat(it);

 set(htimelabel,'userdata',dat);

	%determine the yticks
	% get the trace etc
	trcs=get(hmsg,'userdata');
	t=trcs{2};

	%determine the tick interval
	dat=get(htimeint,'userdata');
	it=find(dat<0);
	dt=get(abs(dat(it)),'userdata');

	%determine the labeling interval
	dat=get(htimelabel,'userdata');
	it=find(dat<0);
	ievery=get(abs(dat(it)),'userdata');

	%first tick at the first even interval smaller than t(1)
	tnot=floor(t(1)/dt)*dt;
	if( abs((tnot-t(1)))>.2*dt )
		tnot=ceil(t(1)/dt)*dt;
	end

	ticks=tnot:dt:t(length(t));
	%generate tick labels
	lblmat=[];
	for k=1:length(ticks)
		if( rem(k-1,ievery)== 0)
			lblmat=strmat(lblmat,num2str(ticks(k)));
		else
			lblmat=strmat(lblmat,' ');
		end
	end
	ind=find(lblmat==1);
   if(~isempty(ind))
      lblmat(ind)=32;
   end		
	set(gca,'ytick',ticks,'yticklabel',lblmat);

end


if(strcmp(action,'quit'))
	close(gcf);
	return;
end

%compute a new phase rotation
if(strcmp(action,'newphs'))
	% ask for the phase angle
	askthingsinit('theodisplay(''newphs2'')','Phase angle (degrees):');
	return;
end

if(strcmp(action,'newphs2'))
	h=get(gcf,'userdata');
	hmsg=h(1);
	hoptions=h(2);
	hgrid=h(3);
	hphase=h(6);
	hhard=h(7);
	hquit=h(8);
	
	%get the angle
	a=askthingsfini;
	if(a==-1) %test for a cancel
		return;
	end
	
	ang=sscanf(a,'%f');
	if( isempty(ang) | ang<-180 | ang>180 )
		set(hmsg,'string','Angle must be between -180 and 180');
		askthingsinit('theodisplay(''newphs2'')','Phase angle (degrees):',a);
		return;
	end
	
	% get the trace etc
	trcs=get(hmsg,'userdata');
	theo=trcs{1};
	t=trcs{2};
	
	trcs=get(hquit,'userdata');
	
	phs=get(hphase,'userdata');
	
	%rotate
	newtrc=phsrot(theo,ang);
	
	%insert in the proper place
	ind=surround(phs,ang);
	if(isempty(ind))
		if(ang<phs(1))
			trcs=[newtrc trcs];
			phs=[ang phs];
		elseif( ang>phs(length(phs)))
			trcs=[trcs newtrc];
			phs=[phs ang];
		end
	else
		trcs=[trcs(:,1:ind) newtrc trcs(:,ind+1:length(phs))];
		phs=[phs(1:ind) ang phs(ind+1:length(phs))];
	end
	
	%set the userdata
	set(hquit,'userdata',trcs);
	set(hphase,'userdata',phs);
	
	%plot
	
	theodisplay('plot');
	
	return;
	
end

% hardcopy
if(strcmp(action,'hardcopy') )
	h=get(gcf,'userdata');
	hhard=h(7);
	
	scales=get(hhard,'userdata');

	shardcopy(gcf,'theodisplay(''hardcopy2'')',scales(1),scales(2),...
		'temp.ps');
	return;
end

if(strcmp(action,'hardcopy2') )
	h=get(gcf,'userdata');
	hhard=h(7);
	hmessage = h(1);
	
	% get the dialog answers and test for reasonableness
	[ps_scale,xlength,ylength,xscale,yscale]=shardcopyfini;

	if( ps_scale== -999.) %test for a cancel
		set(hmessage,'string','Plot cancelled');
		return;
	end
	
%put out a message
set(hmessage,'string',...
      ['plotsize is ' num2str(xlength) ' by ' num2str(ylength) ...
		' inches. use ps_scale = ' num2str(ps_scale) ' in CHVSUB']);
      
	%remember the scales
	set(hhard,'userdata',[xscale yscale]);
	
	return
end
		
%turn on locate info
if(strcmp(action,'locate'))
	h=get(gcf,'userdata');
	hmsg=h(1);

	hlocate=h(10);

	flag=get(hlocate,'checked');

	if(strcmp(flag,'off'))
		set(hlocate,'checked','on');

		set(gcf,'windowbuttonmotionfcn','theodisplay(''showloc'')');

		set(hmsg,'string','Move mouse to display location');

	else
		set(hlocate,'checked','off');
		set(gcf,'windowbuttonmotionfcn','');
		set(hmsg,'string','Locate turned off');
	end

	return;

end

if(strcmp(action,'showloc'))
	h=get(gcf,'userdata');
	%make sure we have a theogram display window
	doit=0;
	if(length(h)==12)
		% if( get(h(3),'type','uimenu') )
	
			if(strcmp(get(h(3),'label'),'Timing Lines'))
				doit=1;
			end
		% end
	end
	if(doit)
		hmsg=h(1);

		pt=get(gca,'currentpoint');

		set(hmsg,'string',['Time = ' num2str(pt(1,2))]);
	end

	return;

end

if(strcmp(action,'writelas'))
	h=get(gcf,'userdata');
	hmsg=h(1);
	hoptions=h(2);
	hphase=h(6);
	hquit=h(8);
	hwritelas=h(11);
	% hphase = h(6) ... vector of theogram phases
	% hquit = h(8) ... matrix of phase rotated theograms

	%put up a file dialog
	[filename,path]=uiputfile('*','Select LAS Output File');

	if(isempty(filename)|filename==0)
		set(hmsg,'string','Output aborted,no file name given');
		return;
	end

	fullfilename=[path filename];

	%set the watch cursor
	set(gcf,'pointer','watch');

	%get the lasheader
	lasheader=get(hwritelas,'userdata');

	%get the theogram
	trcmat=get(hmsg,'userdata');
	theo=trcmat{1};
	t=trcmat{2};
	
	%get the phase rotated traces
	phstrcs=get(hquit,'userdata');
	
	%get the phases
	phs=get(hphase,'userdata');
	
	%generate names for the lasheader
	nphs=length(phs);
	phsnames=32*ones(nphs,4);
	phsdesc='Phase rotated theo';
	phsdesc=ones(nphs,1)*phsdesc;
	for k=1:nphs
		sphs=int2str(abs(phs(k)));
		bb=blanks(3-length(sphs));
		if(phs(k)>0)
			phsnames(k,:)=['+' sphs bb];
		else
			phsnames(k,:)=['-' sphs bb];
		end
	end
	%modify the header
	u=setstr(32*ones(nphs,1));
	lasheader=lashcurves(lasheader,setstr(phsnames),u,u,phsdesc);

	%get the wavelet
% 	w=get(hoptions,'userdata');
% 	w=w(:,1);
% 	%pad wavelet with -999.25 which las will interprete as nan's
% 	[ns,m]=size(trcmat);
% 	nw=length(w);
% 	if(ns>nw)
% 		w=[w;-999.25*ones(ns-nw,1)];
% 	elseif(nw>ns)
% 		t=trcmat(:,1);
% 		t=xcoord(t(1),t(2)-t(1),nw);
% 		trcmat=[t [trcmat(2:4);-999.25*ones(nw-ns,3)]];
% 		phstrcs=[phstrcs; -999.25*ones(nw-ns,nphs)];
% 	end
		

	%reorder the trcmat
	%trcmat=trcmat(:,[2 1 3 4]);
	
	%attach the the phases
	lasmat=[t theo phstrcs];

	%write out
	writelas(fullfilename, lasheader, lasmat);

	set(hmsg,'string',['Output successful to ' fullfilename]);

	set(gcf,'pointer','arrow');

	return;
end