SAS：如何将 sashelp.zipcode 坐标投影到 maps.states 上？

Question

我想在美国地图上标出数千个邮政编码。 sashelp.zipcode 的坐标是 latitude/longitude 的度数，maps.states 的坐标是未投影的 latitude/longitude 弧度。如果我使用以下方法将邮政编码坐标转换为弧度：

x=atan(1)/45*x;
y=atan(1)/45*y;

然后我用

画了一张地图

data states;
   set maps.states;
   where state not in (2, 15, 72);
   x = -x; 
run;

邮政编码绘制准确，但当然地图有（轻微）扭曲。我想要得到你在投影 maps.states 时得到的美国地图的丰满曲线版本。

我承认 gproject 过程对我来说有点不透明。我试过：

data zips; 
    set <dataset with 86,000 zip codes>; 
    where state not in (2, 15, 72, .); 
    function='symbol'; 
    size = .02; 
    text = "dot"; 
    xsys='2'; 
    ysys='2'; 
    hsys='3'; 
    when='a';
    color= "%hsv(100,100,100)";
    x=atan(1)/45*x;
    y=atan(1)/45*y;
    x=-x; 
run;

data states; 
    set maps.states;
    where state not in (2, 15, 72);  
run;

data combo; 
    set zips states; 
run; 

proc gproject data=combo
              out=comboproj;
   id state;
run;

data zipsproj mapproj; 
    set comboproj; 
    if when = "a" then output zipsproj; 
    if when = "" then output mapproj; 
run; 

proc gmap map=mapproj 
          data=mapproj 
          anno = zipsproj all;
   id state;
   choro state / nolegend levels=1;
run;
quit;

这给了我一个横向的美国地图，所有的 zip 坐标都超出了范围（编辑我也通过将 zip 坐标乘以 -1 来解决范围问题，但是地图仍然是横向的）。

我做错了什么？请帮忙！

Answer 1

你的问题是你没有告诉注释系统你有什么样的坐标。请记住，注释不必只位于图形本身的坐标系中——它们还可以放在 off-graph 位置（例如自定义轴标签、自定义标题等），并且可以给定以百分比而不是实际坐标表示（例如，将右边 20% 的颜色着色为红色，将左边的 20% 着色为绿色）。您使用 xsys、ysys，如果您指定字符的高度，则 hsys。 2 表示 'absolute, data values'； 3 表示 'absolute, percent of graph area'。有关坐标系的详细信息，请参阅 this paper。

您也不应该将状态 x 乘以 -1，至少在我的设置中不是这样。只有邮政编码 x 应乘以 -1。但是，如果您的设置需要这样做，当然可以继续这样做。

最后，当 GPROJECT 在事物的 negative-x 方面看到一些值时，地图会被 GPROJECT 转向一边。这导致它在选择 PARALLEL1/2 时做出奇怪的选择（您可以在日志中看到这种情况）。消除所有非 US48 州（关岛、密克罗尼西亚等）以达到较低的 48 个。

这是适合我的最终代码：

data states;
   set maps.states;
   where state not in (2, 15, 72);
run;

data zips;
  set sashelp.zipcode;
     x=-1*x;

x=atan(1)/45*x;
y=atan(1)/45*y;
if state ne 2 and state ne 15 and state lt 60;
  when='a';
  function='symbol';
  text='dot';
  xsys='2';
  ysys='2';
  hsys='3';

run;

data combo; 
    set zips  states; 
run; 

proc gproject data=combo
              out=comboproj;
   id state;
run;

data zipsproj mapproj; 
    set comboproj; 
    if when = "a" then output zipsproj; 
    if when = "" then output mapproj; 
run; 

proc gmap map=mapproj 
          data=mapproj 
          anno = zipsproj all density=4 ;
   id state;
   choro state / nolegend levels=1;
run;
quit;

SAS：如何将 sashelp.zipcode 坐标投影到 maps.states 上？

SAS: How to project sashelp.zipcode coordinates onto maps.states?

sas

maps