zip2lonlat.m

物流分析工具包。Facility location: Continuous minisum facility location, alternate location-allocation (ALA)

function [XY,idx,flg] = zip2lonlat(zip)
%ZIP2LONLAT Determine longitude and latitude of ZIP codes.
% [XY,idx,flg] = zip2lonlat(zip)
%     zip = m-element vector of 5-digit zip codes
%      XY = m x 2 matrix of lon-lat (in decimal degrees)
%           ([NaN NaN] returned if 5- and 3-digit code not found)
%     idx = m-element index vector of data used, where idx(flg==5) are the
%           indices for USZIP5 and idx(flg==3) are the indices for USZIP3
%           (0 returned if 5- and 3-digit code not found)
%     flg = 5, if 5-digit code from USZIP5 used
%              (not all USZIP5 codes are unique, so code with max pop used)
%         = 3, if 3-digit code from USZIP3 used because 5-digit code not
%              found (not all possible codes are in USZIP5, e.g., PO Boxes)
%         = 0, if 3-digit code is not found
%
% See also LONLAT2CITY

% Copyright (c) 1994-2006 by Michael G. Kay
% Matlog Version 9 13-Jan-2006 (http://www.ie.ncsu.edu/kay/matlog)

% Input Error Checking ****************************************************
zip = zip(:);
if any(zip < 601 | zip > 99950)
   error('Valid zip codes only in the range of 601 to 99950');
end
% End (Input Error Checking) **********************************************

m = length(zip);
XY = NaN*ones(m,2);
idx = zeros(m,1);
flg = zeros(m,1);
[Code5,XY5,Pop] = uszip5('Code5','XY','Pop');
[Code3,XY3] = uszip3('Code3','XY');
for i = 1:length(zip)
   idxi = find(Code5 == zip(i));
   if ~isempty(idxi)
      if length(idxi) > 1
         idxi = idxi(argmax(Pop(idxi)));  % Not all Code5's are unique
      end                                 % so use code with max pop
      XY(i,:) = XY5(idxi,:);
      idx(i) = idxi;
      flg(i) = 5;
   else
      idxi = find(Code3 == floor(zip(i)/100));
      if ~isempty(idxi)
         XY(i,:) = XY3(idxi,:);
         idx(i) = idxi;
         flg(i) = 3;
      end
   end
end