乌尔夫网

%画纬线弧
function gridx(a)
i=0:-.001:-pi;
for j=0:180/a:180
plot(cos(i)*sin(j/180*pi)./(1-sin(i)*sin(j/180*pi)),cos(j/180*pi)./(1-sin(i)*sin(j/180*pi)));
end;
//释放内存和关闭文件
fclose(pDummyFile);
fclose(pWritingFile);
free(pPixelData);
}
Test6源文件
#include <GL/glut.h>
#include<math.h>
#define PI 3.1415926
void gridy(float,long double);
二、乌尔夫网的计算机实现
在有了乌尔夫网的数学描述之后，便可以用计算机程序语言或建模工具来实现它。在这里我用到了本学期地学信息建模与软件工程中学到的OpenGL图形库和当今最为优秀的建模工具matlab来分别实现它。
1、OpenGL
该工程由7个源文件及一个1×1的位图文件组成，最终会生成一个位图文件来保存投影结果。
{
long double i;
glMatrixMode(GL_MODELVIEW);
glRotatef(a,0,0,-1);
glColor3f(.2,0.5,1);
glBegin(GL_LINE_LOOP);
for(i=0;i<=PI;i+=.005)
{
glVertex2d(cos(i)/(1+sin(i)*sin(b/180*PI)),sin(i)*cos(b/180*PI)/(1+sin(i)*sin(b/180*PI)));
//计算像素数据的实际长度
i = WindowWidth * 3; //得到每一行的像素数据长度
while( i%4 != 0 ) //补充数据，直到i是的倍数
++i; //本来还有更快的wenku.baidu.com法，
PixelDataLength = i * WindowHeight;
//分配内存和打开文件
pPixelData = (GLubyte*)malloc(PixelDataLength);
glutDisplayFunc(&myDisplay);
glutMainLoop();
return 0;
}
数据源：可以用记事本对数据进行编辑和保存。
显示投影效果：
2、Matlab
M文件如下：
%画基圆
x=0:.01:2*pi;
plot(cos(x),sin(x));
hold on
axis('equal');
}
glEnd();
glLoadIdentity();
}
Test4源文件
//定义了一个从C文本文件获取系列面状地质要素产状数据的函数，便于//数据的保存和使用
#include<stdio.h>
#include<stdlib.h>
void mian(float,long double);
void wu(void)
function mian(k,j)
i=0:.001:pi;
m=cos(i)./(1+sin(i)*sin(j/180*pi));
n=sin(i)*cos(j/180*pi)./(1+sin(i)*sin(j/180*pi));
glRotatef(a,0,0,-1);
glBegin(GL_LINE_STRIP);
for(i=0;i<=PI;i+=.005)
{glVertex2d(cos(i)/(1+sin(i)*sin(b/180*PI)),sin(i)*cos(b/180*PI)/(1+sin(i)*sin(b/180*PI)));
fread(BMP_Header, sizeof(BMP_Header), 1, pDummyFile);
fwrite(BMP_Header, sizeof(BMP_Header), 1, pWritingFile);
fseek(pWritingFile, 0x0012, SEEK_SET);
i = WindowWidth;
{
FILE *fp;
float a,b;
if((fp=fopen("b1","r"))==NULL)
{
printf("open file error\n");
exit(0);
}
while(!feof(fp))
{
fscanf(fp,"%f,%f",&a,&b);
mian(a,b);
}
fclose(fp);
if( pPixelData == 0 )
exit(0);
pDummyFile = fopen("dummy.bmp", "rb");
if( pDummyFile == 0 )
exit(0);
pWritingFile = fopen("grab.bmp", "wb");
if( pWritingFile == 0 )
glBegin(GL_LINE_STRIP);
for(i=-PI;i<=0;i+=.005)
{glVertex2d(cos(i)*sin(a/180*PI)/(1-sin(i)*sin(a/180*PI)),cos(a/180*PI)/(1-sin(i)*sin(a/180*PI)));
}
glEnd();
void gridx(long double);
void wu(void);
void grab(void);
void myDisplay(void)
{
int i;
long double j;
glClear(GL_COLOR_BUFFER_BIT);
glColor3f(.1,.1,.1);
glBegin(GL_LINES);
#define BMP_Header_Length 54
void grab(void)
{
FILE* pDummyFile;
FILE* pWritingFile;
GLubyte* pPixelData;
GLubyte BMP_Header[BMP_Header_Length];
GLint i, j;
GLint PixelDataLength;
%画经线弧
function gridy(a)
i=0:.001:pi;
for j=0:180/a:180
plot(sin(i)*cos(j/180*pi)./(1+sin(i)*sin(j/180*pi)),cos(i)./(1+sin(i)*sin(j/180*pi)));
end;
%投影产状为k°/j°的面状要素
}
glEnd();
glLoadIdentity();
}
Test2源文件
//定义了绘制乌尔夫网纬线的函数
#include<GL/glut.h>
#include<stdio.h>
#include<math.h>
#define PI 3.1415926
void gridx(long double a)
{
long double i;
gridx(i);
}
wu();
glFlush();
grab();
}
Test7源文件
#include <GL/glut.h>
#define WindowWidth 600
#define WindowHeight 600
void myDisplay(void);
int main(int argc, char *argv[])
{
glutInit(&argc, argv);
glutInitDisplayMode(GLUT_RGBA | GLUT_SINGLE);
glutInitWindowPosition(100, 100);
glutInitWindowSize(WindowWidth,WindowHeight);
glutCreateWindow("Lee's OpenGL program");
}
Test5源文件
//定义一个抓屏函数，将投影结果以一个600×600的24位位图保存，以//便于投影结果的保存
#define WindowWidth 600
#define WindowHeight 600
#include<gl/glut.h>
#include <stdio.h>
#include <stdlib.h>
乌尔夫网
一、乌尔夫网的数学描述
1、经线弧
经线弧相当于一系列走向为90ª，等倾角间隔的一系列面状要素的极射赤平影。
这一系列面状要素的球面投影弧线分别相当于基圆绕x周作其对应倾角的旋转。易知，基圆的右半圆可以用点集A={（cos(i),sin(i)）|0=<i<=pi}来描述，则由坐标系的旋转原理可知，产状为90°/jª(90为倾向，j为倾角)的负半球球面投影弧线点集可以由A乘以坐标系旋转矩阵得到，即为B={（cos(i),sin(i)*cos(j),-sin(i)*sin(j)）|0=<i<=pi}。
j = WindowHeight;
fwrite(&i, sizeof(i), 1, pWritingFile);
fwrite(&j, sizeof(j), 1, pWritingFile);
//写入像素数据
fseek(pWritingFile, 0, SEEK_END);
fwrite(pPixelData, PixelDataLength, 1, pWritingFile);
Test1源文件
//定义了绘制乌尔夫网经线的函数
#include<GL/glut.h>
#include<stdio.h>
#include<math.h>
#define PI 3.1415926
void gridy(float a,long double b)
{
long double i;
glMatrixMode(GL_MODELVIEW);
由赤平投影投影定义可知，点（0,0,1）与A上点集的连线与基圆的交点集即为纬线弧，即为
B=(cos(i)*sin(j/180*pi)/(1+sin(i)*sin(j/180*pi)),cos(j/180*pi)/(1+sin(i)*sin(j/180*pi))|0=<i<=pi);
至于纬线弧的几何性质，这里并未深究，但可以想见它不会是圆弧。
glLoadIdentity();
}
Test3源文件
//定义了对面状地质要素进行乌尔夫网投影的函数
#include<GL/glut.h>
#include<stdio.h>
#include<math.h>
#define PI 3.1415926
void mian(float a,long double b)
exit(0);
//读取像素
glPixelStorei(GL_UNPACK_ALIGNMENT, 4);
glReadPixels(0, 0, WindowWidth, WindowHeight,
GL_BGR_EXT, GL_UNSIGNED_BYTE, pPixelData);
//把dummy.bmp的文件头复制为新文件的文件头
glVertex2f(-1,0);
glVertex2f(1,0);
glVertex2f(0,-1);
glVertex2f(0,1);
glEnd();
//glLineWidth(3);
glBegin(GL_LINE_LOOP);
for(j=0;j<2*PI;j+=.005)
{
glVertex2d(cos(j),sin(j));
由赤平投影投影定义可知，点（0,0,1）与B上点集的连线与基圆的交点集即为经线弧，即为
C={(sin(i)*cos(j/180*pi)/(1+sin(i)*sin(j/180*pi)),cos(i)/(1+sin(i)*sin(j/180*pi)))|0=<i<=pi}。
经解析几何证明，该弧线为一段椭圆弧线，而非圆弧。
}
glEnd();
for(i=0;i<=180;i+=10)
{
//glColor3f((float)(i&0x04),(float)(i&0x02),(float)(i&0x01));
gridy(90,i);
}
for(i=0;i<=180;i+=10)
{
//glColor3f((float)(i&0x04),(float)(i&0x02),(float)(i&0x01));
2、纬线弧
纬线圈相当于一系列产状为0°/90ª的互相平行行但不经过基圆心的面状要素的极射赤平投影。这一系列面状要素切球面所成的球冠所对的锥角是等间隔的。所对锥角为j°的面的负半球球面投影弧线点集可以用点集A={(cos(j),sin(j)*cos(i),-sin(j)*sin(i)) |0=<i<=pi}来描述。