matlab_合成孔径雷达干涉测量技术 INSAR

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % BranchCuts.m generates branch cuts based on the phase residues. This is % done using the Goldstein method, as described in "Two-dimensional phase % unwrapping: theory, algorithms and software" by Dennis Ghiglia and % Mark Pritt. % "residue_charge" is a matrix wherein positive residues are 1 and % negative residues are 0. % "max_box_radius" defines the maximum search radius for the balancing of % residues. If this is too large, areas will be isolated by the branch % cuts. % "IM_mask" is a binary matrix. This serves as an artificial border for the % branch cuts to connect to. % Created by B.S. Spottiswoode on 15/10/2004 % Last modified on 18/10/2004 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% function branch_cuts=BranchCuts(residue_charge, max_box_radius, IM_mask); [rowdim, coldim]=size(residue_charge); branch_cuts=~IM_mask; %Define initial branch cuts borders as the mask. residue_charge_masked=residue_charge; residue_charge(logical(~IM_mask))=0; %Remove all residues except those in the mask cluster_counter=1; %Keep track of the number of residues in each cluster satellite_residues=0; %Keep track of the number of satellite residues accounted for residue_binary=(residue_charge~=0); %Logical matrix indicating the position of the residues residue_balanced=zeros(rowdim, coldim); %Initially designate all residues as unbalanced [rowres,colres] = find(residue_binary); %Find the coordinates of the residues adjacent_residues=zeros(rowdim, coldim); %Defines the positions of additional residues found in the search box missed_residues=0; %Keep track of the effective number of residues left unbalanced because of disp('Calculating branch cuts ...'); tic; temp=size(rowres); for i=1:temp(1); %Loop through the residues radius=1; %Set the initial box size r_active=rowres(i); %Coordinates of the active residue c_active=colres(i); count_nearby_residues_flag=1; %Flag to indicate whether or not to keep track of the nearby residues cluster_counter=1; %Reset the cluster counter adjacent_residues=zeros(rowdim, coldim); %Reset the adjacent residues indicator charge_counter=residue_charge_masked(r_active, c_active); %Store the initial residue charge if residue_balanced(r_active, c_active)~=1 %Has this residue already been balanced? while (charge_counter~=0) %Loop until balanced %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %This portion of code serves to search the box perimeter, %place branch cuts, and keep track of the summed residue charge %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% for m=r_active-radius:r_active+radius %Coordinates of the box border pixels ***I COULD MAKE THIS MORE EFFICIENT*** for n=c_active-radius:c_active+radius if (abs(m - r_active)==radius | abs(n - c_active)==radius) & charge_counter~=0 %Ensure that only the border pixels are being scrutinised if m<=1 | m>=rowdim | n<=1 | n>=coldim %Is the current pixel on the image border? if m>=rowdim m=rowdim; end %Make sure that the indices aren't too large for the matrix if n>coldim n=coldim; end if n<1 n=1; end if m<1 m=1; end branch_cuts=PlaceBranchCutInternal(branch_cuts, r_active, c_active, m, n); %Place a branch cut between the active point and the border cluster_counter=cluster_counter+1; %Keep track of how many residues have been clustered charge_counter=0; %Label the charge as balanced residue_balanced(r_active, c_active)=1; %Mark the centre residue as balanced end if IM_mask(m,n)==0 branch_cuts=PlaceBranchCutInternal(branch_cuts, r_active, c_active, m, n); %Place a branch cut between the active point and the mask border cluster_counter=cluster_counter+1; %Keep track of how many residues have been clustered charge_counter=0; %Label the charge as balanced residue_balanced(r_active, c_active)=1; %Mark the centre residue as balanced end if residue_binary(m,n)==1 %Is the current pixel a residue? if count_nearby_residues_flag==1 %Are we keeping track of the residues encountered? adjacent_residues(m,n)=1; %Mark the current residue as adjacent end branch_cuts=PlaceBranchCutInternal(branch_cuts, r_active, c_active, m, n); %Place a branch cut regardless of the charge_counter value cluster_counter=cluster_counter+1; %Keep track of how many residues have been clustered if residue_balanced(m,n)==0; residue_balanced(m,n)=1; %Mark the current residue as balanced charge_counter=charge_counter + residue_charge_masked(m,n); %Update the charge counter end if charge_counter==0 %Is the active residue balanced? residue_balanced(r_active, c_active)=1; %Mark the centre (active) residue as balanced end end end end end %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %This next portion of code centres the box on the adjacent %residues. If the charge is still not balanced after moving %through all adjacent residues, increase the box radius and %centre the box around the initial residue. %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% if sum(sum(adjacent_residues))==0 %If there are no adjacent residues: radius=radius+1; %Enlarge the box r_active=rowres(i); %Centre the larger box about the original active residue c_active=colres(i); else %If there are adjacent residues: if count_nearby_residues_flag==1; %Run this bit once per box being searched [r_adjacent,c_adjacent] = find(adjacent_residues); %F