VB测量平差程序的设计讲稿.docx

VB测量平差程序的设计讲稿.docx

《VB测量平差程序的设计讲稿.docx》由会员分享，可在线阅读，更多相关《VB测量平差程序的设计讲稿.docx（18页珍藏版）》请在冰豆网上搜索。

VB测量平差程序的设计讲稿

Case0'读入观测值文件

Text1.Visible=False

CommonDialog1.ShowOpen

fname=CommonDialog1.FileName'将用户在"打开"对话框中选择的文件名对变量fname赋值

Iffname<>""Then'若无此判断当对话框中选择取消时、下面赋值语句将出错

Setts=fso.OpenTextFile（fname）'将fname作为文本文件打开,并设置句柄

j=0:

k=0:

p=0:

h=0

'j是测站数累计变量,k是已知点累计变量,l（j）、ns（j）分别是方向值、边长累积计数

DoWhilets.AtEndOfLine<>True'前测型循环，进入循环的条件是没有读到文件结束尾

B=ts.ReadLine'读一行，置入b

B=Trim（B）:

i=1:

'删除B可能有的前导和尾随空格，i是工作变量，

m（i）=InStr（B,","）'查行中第一个逗号的左数位置,并保存在整形数组变量m（i）

DoWhilem（i）<>0'前测型Do...Loop循环，成立条件是该行字符串中有逗号

tr（i）=Mid（B,m（i-1）+1,m（i）-m（i-1）-1）'提取指定位置开始的指定数目字符。

i=i+1

m（i）=InStr（m（i-1）+1,B,","）'从上一个找到的逗号位置起，查找下一个逗号的位置

Loop

Ifm（i）=0Andi>1Thentr（i）=Right（B,Len（B）-m（i-1））'处理一行中最后一个逗号后的字符串

'以下部分是将存储在数组变量m（i）中的字符分类存放到方向、边长、已知坐标、网型信息等数组中

Ifp=0Then'读到的是文件第一行。

ma=tr

（1）:

ms=tr

（2）:

mk=tr（3）:

p=1'提取观测方向,边先验精度值，并使该句以后不能再执行。

Else

Ifm

（1）=0Thenj=j+1:

ReDimPreservedm（j）:

ReDimPreservenl（j）:

ReDimPreservens（j）:

dm（j）=B:

nl（j）=nl（j-1）:

ns（j）=ns（j-1）

'该行中没有逗号,但又不是结束符，则一定是测站点名。

将读出的字符串赋值到点名数组变量dm（j），

Nextj

Ifp=0Thend=d+1:

ReDimPreservens（d）:

ns（d）=ns（d-1）:

ReDimPreservedm（d）:

dm（d）=lb（i）'如p=0，表明目标点未设过测站,将该点点名加入点名数组

Nexti

zds=d'将总点数存入模块级变量zds

ReDimx（zds）,y（zds）'重新定义坐标数组

x

（1）=10000:

y

（1）=10000'为推算近似坐标,对第一个点赋假设坐标值

k=1

Fori=1Toyds

Iflb

（1）=ym（i）Thenk=k+1

Nexti

ss=sid（1,k）'调出第一点到未知点方向的边长，参数是测站点序号，照准方向号

h=seqn（lb（k））'查k方向照准点的计算序号

x（h）=x

（1）+ss*Cos（0）:

y（h）=y

（1）+ss*Sin（0）'计算第一点上第k方向值的目标点假设坐标

Fori=1Tonl（cds）'遍访所有方向值，将其由角度值转换为弧度值.

Ifl（i）>0.001Thenl（i）=radian（l（i））'零方向值不参加转换

Nexti

n=0

Do

n=n+1'n是循环计数变量，控制循环次数，避免假定坐标计算不出时，进入死循环。

Fori=1Tocds'按测站循环

Ifx（i）>1Then'在该测站假设坐标已计算出的前提下，求照准点假设坐标

p1=0

Forj=nl（i-1）+1Tonl（i）'遍访i测站上所有方向值

h=seqn（lb（j））'查目标点对应的序号

Ifx（h）>1Then'目标点坐标已解出

t2=azimuth（xo

（1）,yo

（1）,xo

（2）,yo

（2））

dt=t2-t1:

x1=x（m

（1））:

y1=y（m

（1））

Fori=1Tozds'将假设坐标转换到实际坐标

dx=x（i）-x1:

dy=y（i）-y1

x（i）=xo

（1）+dx*Cos（dt）-dy*Sin（dt）

y（i）=yo

（1）+dx*Sin（dt）+dy*Cos（dt）

Nexti

Fori=1Toyds'置入已知点坐标

x（m（i））=xo（i）:

y（m（i））=yo（i）

Nexti

Case2'组法方程

Text1.Visible=False

Diml1AsDouble,ppAsDouble,n2AsLong'定义过程级变量

q=206265:

ll=0

n1=2*（zds-yds）'未知数数目

n2=n1*（n1+1）/2'一维存储法方程系数数组上限

ReDimNX（n2）,UX（n1）'重新定义法方程系数、常数数组

Callorder（m（）,yds）'对保存已知点序号的m（）数组排序

Fori=1Tocds'按测站循环

z=0'将按测站累积的变量清零

'下面开始处理一个测站的方向观测值

k1=nl（i-1）+1:

k2=nl（i）'一测站上最小和最大方向号

Forj=k1Tok2'在i测站上按方向循环，求定向角未知数

h=seqn（lb（j））

t=azimuth（x（i）,y（i）,x（h）,y（h））

f=t-l（j）:

Iff<0Thenf=f+2*pi'f是解算的零方角方位角

pp=1'方向观测值的权为1

Callequation（nb（）,pp,l1）'组法方程，参数分别是误差方程系数数组、权、误差方程常数项

Nextj

pp=-1/（k2-k1+1）:

l1=ln'pp是和方程的权

Callequation（nc（）,pp,l1）'和方程组法方程,nc（）是一测站方向误差方程和方程数组

'一个测站方向观测值误差方程处理完毕,下面开始处理该测站的边观测值

Ifns（i）-ns（i-1）-1>=0Then'i测站有观测边

Forj=ns（i-1）+1Tons（i）'依次遍访i测站上各观测边

ReDimnb（n1）'重新定义误差方程系数数组,并且每循环到一新边长前清零

h=seqn（sb（j））

t=azimuth（x（i）,y（i）,x（h）,y（h））

A=Cos（t）:

B=Sin（t）

cha=charact（i,k）'自定义函数，查测站点i是否已知点，如不是,用k返回i前面有几个已知点

Ifcha="n"Then'测站点i不是已知点

d=2*（i-k-1）+1'计算测站i点x坐标未知数在未知数点集中的序号

nb（d）=-A:

nb（d+1）=-B

EndIf

cha=charact（h,k）

Ifcha="n"Then'照准方向点h不是已知点

d=2*（h-k-1）+1'计算照准方向h点x坐标未知数在未知数点集中的序号

nb（d）=A:

nb（d+1）=B

EndIf

dx=x（h）-x（i）:

dy=y（h）-y（i）

ss=Sqr（dx^2+dy^2）'反算边长，置于判断式外是因为两已知点之间不会测边

pp=（ma/（0.1*ms+mk*ss*10^-4））^2'边长观测值定权

l1=（ss-s（j））*100'求边误差方程常数（单位是厘米）

'并对该测站最大方向号，最大边长号数组变量nl（j）、ns（j）赋值累计起始值

Ifm

（1）<>0Then'不是第一行，并且该行中有逗号分割的多个字串

tr

（2）=LCase（tr

（2））

Iftr

（2）<>"l"Andtr

（2）<>"s"Then'这一行不是方向或边长观测值，而是已知点坐标值

k=k+1:

ReDimPreserveym（k）:

ReDimPreservexo（k）:

ReDimPreserveyo（k）:

ym（k）=tr

（1）:

xo（k）=Val（tr

（2））:

yo（k）=Val（tr（3））'输入已知点坐标

Else

Iftr

（2）="l"Thennl（j）=nl（j）+1:

ReDimPreservelb（nl（j））:

ReDimPreservel（nl（j））:

lb（nl（j））=tr

（1）:

l（nl（j））=Val（tr（3））'累计测站方向数、提取照准点、方向值

Iftr

（2）="s"Thenns（j）=ns（j）+1:

ReDimPreservesb（ns（j））:

ReDimPreserves（ns（j））:

sb（ns（j））=tr

（1）:

s（ns（j））=Val（tr（3））'提取边观测数、提取照准点、观测边

EndIf

EndIf

EndIf

Loop

ts.Close

cds=j:

yds=k'用模块级变量cds、yds保存测站点总数、已知点总数

MsgBox"数据已成功读入",0+64+0,"数据输入"

EndIf

Case1'解算近似坐标

d=cds'd是测站数

Fori=1Tonl（cds）'依次访问所有的方向值

p=0'设识别变量

Forj=1Tod'依次访问所有测站

Ifdm（j）=lb（i）Thenp=1'查看目标点是否设过测站，是则对识别变量p赋值1。

t=azimuth（x（i）,y（i）,x（h）,y（h））:

p1=j'反算方位角并记下方向号

ExitFor'退出For－Next循环

EndIf

Nextj

Forj=nl（i-1）+1Tonl（i）'再次遍访i测站上的所有方向值

h=seqn（lb（j））'查目标点对应的序号

Ifx（h）<1Then'该照准方向点坐标未求出

tt=t+l（j）-l（p1）'计算j方向的方位角

ss=sid（i,j）'根据测站号、方向号去查找边长

Ifss=0Thenss=side（i,h）

x（h）=x（i）+ss*Cos（tt）:

y（h）=y（i）+ss*Sin（tt）'计算近似坐标

EndIf

Nextj

EndIf

Nexti

p=0

Fork=1Tozds

Ifx（k）<1Thenp=1'查看是否还有没解算出坐标的点,有则对p赋值，使进入循环再次搜索计算

Nextk

LoopUntilp=0Orn>zds'坐标已全部算出或虽还有未算出的，但根据循环次数已不能算出时结束循环。

Fori=1Toyds'查点名数组中那些点是已知点,用m（）数组存其序号

Forj=1Tozds

Ifym（i）=dm（j）Thenm（i）=j

Nextj

Nexti

t1=azimuth（x（m

（1））,y（m

（1））,x（m

（2））,y（m

（2）））'利用一对公共点求旋转角

z=z+f

Nextj

zo=z/（k2-k1+1）'zo是定向角未知数（零方向的坐标方位角）

ReDimnc（n1）:

ln=0'重新定义和方程系数数组，并且每循环到一新测站前清零，ln是和方程常数项

Forj=k1Tok2'在i测站上按方向循环，求误差方程系数、常数

ReDimnb（n1）'重新定义误差方程系数数组,并且每循环到一新方向前清零

h=seqn（lb（j））

dx=x（h）-x（i）:

dy=y（h）-y（i）

ss=Sqr（dx^2+dy^2）'反算边长

t=azimuth（x（i）,y（i）,x（h）,y（h））'反算坐标方位角

A=q*Sin（t）*10^-2/ss:

B=-q*Cos（t）*10^-2/ss'求方向误差方程系数

cha=charact（i,k）'自定义函数，查测站点i是否已知点，如不是,用k返回i前面有几个已知点

Ifcha="n"Then'序号为i的测站点不是已知点

d=2*（i-k-1）+1'计算测站i点x坐标未知数在未知数点集中的序号

nb（d）=A:

nb（d+1）=B'对误差方程系数数组赋值

nc（d）=nc（d）+nb（d）:

nc（d+1）=nc（d+1）+nb（d+1）'累积到和方程数组中去

EndIf

cha=charact（h,k）'查照准方向点h是否已知点，如不是,用k返回h前面有几个已知点

Ifcha="n"Then'照准方向点不是已知点

d=2*（h-k-1）+1'计算照准方向点x坐标未知数在未知数点集中的序号

nb（d）=-A:

nb（d+1）=-B'对误差方程系数数组赋值

nc（d）=nc（d）+nb（d）:

nc（d+1）=nc（d+1）+nb（d+1）'累积到和方程数组中去

EndIf

Ift-zo<-0.01Thent=t+2*pi'当t是零方向方位角时，可能比平均值zo略小，所以不用<0的比较方法

l1=q*（t-zo-l（j））'误差方程常数项,以秒为单位

ln=ln+l1:

ll=ll+l1^2'累积和方程常数项及求[pvv]的[Pll]值

Callequation（nb（）,pp,l1）

ll=ll+pp*l1^2'累积用于求[pvv]的[Pll]值

Nextj

EndIf

'一个测站的误差方程处理完毕，进入下一测站

Nexti

Text1.Visible=False

MsgBox（"法方程组成完毕"）

Case3'高斯约化法解法方程

Dimm2AsLong

n=2*（zds-yds）

ReDimnb（n）,nc（n）'数组清零

Fori=1Ton

nb（i）=-UX（i）'将法方程常数项赋予一工作数组,保留原值用于求pvv

Nexti

Fork=1Ton-1'从法方程的第一行到倒数第二行,k实际上控制约化次数

m1=k*（k+1）/2'm1是k行的自乘元素在法方程系数阵中序号（一维按列上三角阵）

Forj=k+1Ton'在j循环内完成k次约化

m2=（（j-1）*j/2）+k'm2在j循环内依次是k行自乘系数右边各元素序号

Fori=jTon'i循环内完成法方程系数阵j行元素的k次约化

d=（i-1）/2*i

n1=d+k'与约化元素同列的k行元素

n2=d+j'约化元素在法方程系数阵中编号

NX（n2）=NX（n2）-NX（n1）*NX（m2）/NX（m1）

Nexti

Nextj

Fori=k+1Ton'在i循环内完成法方程常数项的k次约化

nub=（i-1）*i/2+k

nb（i）=nb（i）-nb（k）*NX（nub）/NX（m1）

Nexti

Nextk'法方程约化完毕

pvv1=ll'll是[pll]

Fork=1Ton'在循环内按[pvv]=[pll]+（l）*（l）求[pvv]

nub=k*（k+1）/2

pvv1=pvv1-nb（k）*nb（k）/NX（nub）

Nextk

Fork=nTo1Step-1'在k循环内回代求解未知数

u=nb（k）

Ifk

Fori=k+1Ton

nub=（i-1）*i/2+k

u=u-NX（nub）*nb（i）

Nexti

nb（k）=u

EndIf

nub=k*（k+1）/2

nb（k）=nb（k）/NX（nub）

Nextk

pvv2=ll

Fork=1Ton'按[pvv]=[pll]+W*X求[[pvv]

pvv2=pvv2+UX（k）*nb（k）

MsgBox（"法方程解算完毕"）

EndSelect

EndSub

PrivateFunctionazimuth（x1AsDouble,y1AsDouble,x2AsDouble,y2AsDouble）'反算坐标方位角

dx=x2-x1:

dy=y2-y1

azimuth=Atn（dy/dx）

Ifdx<0Then

azimuth=azimuth+pi

Else

Ifdy<0Thenazimuth=azimuth+2*pi

EndIf

EndFunction

PrivateFunctionseqn（strAsString）AsInteger'由点名查计算序号函数

Fori=1Tozds

Ifstr=dm（i）Thenseqn=i'将查到的序号赋给函数名,返回调用处

Nexti

EndFunction

PrivateFunctionsid（nAsLong,mAsLong）AsDouble'提取边长函数，n是测站号，m是方向号

ss=0

Ifns（n）-ns（n-1）-1>=0Then'该测站有观测边，确定m方向的边长

Fork=ns（n-1）+1Tons（n）

Ifsb（k）=lb（m）Thenss=s（k）'依次将边的照准点名与m方向照准点名对比，相等则赋予ss。

Nextk

EndIf

Ifns（n）-ns（n-1）-1<0Orss<1Then'如测站n未找到m方向的边，则转到lb（m）站去找

h=seqn（lb（m））'由点名查计算编号

NX（h）=NX（h）+B（j）*B（i）*p'组法方程系数阵

EndIf

Nextj

EndIf

Nexti

EndSub

PrivateSuborder（m（）AsInteger,nAsInteger）'排序通用过程

Fori=1Ton-1

Forj=i+1Ton

Ifm（j）

A=m（j）:

m（j）=m（i）:

m（i）=A

EndIf

Nextj

Nexti

EndSub

PrivateFunctionside（nAsLong,mAsLong）'反算已知点间边长函数过程

Fori=1Toyds

Ifym（i）=dm（n）Thenk1=i

Ifym（i）=dm（m）Thenk2=i

Nexti

dx=xo（k1）-xo（k2）:

dy=yo（k1）-yo（k2）

side=Sqr（dx^2+dy^2）

EndFunction

PrivateSubinversion（A（）AsDouble,W（）AsDouble）'法方程系数阵求逆子程序

n=2*（zds-yds）'形参数组A（）是约化后法方程系数数组,W（）是工作数组。

Fori=nTo1Step-1'i控制所求逆阵元素的行数，由下而上。

EndIf

Nexti

EndSub'求逆完毕,此时约化法方程系数阵中元素已是逆阵元素

8888888

77777

0000

9999

00000

00000

9999

Nextk

Fori=1Tozds'求坐标平差值

cha=charact（i,k）'自定义函数，查点i是否已知点，如不是,用k返回i前面有几个已知点

Ifcha="n"Then'序号i的点不是已知点

h=i-k

d=2*（h-1）+1'计算i点x坐标未知数在未知数点集中的序号

x（i）=x（i）+nb（d）/100:

y（i）=y（i）+nb（d+1）/100

Else

EndIf

Nexti

ma=Sqr（pvv1/（nl（cds）+ns（cds）-n-cds））'求单位权中误差

Callinversion（NX（）,nb（））'将约化后的法方程系数阵送入通用过程inversion（）求逆

ReDimnb（n）'工作数组清零

Fori=1Tozds'按点号循环,求点位中误差

cha=charact（i,k）'自定义函数，查点i是否已知点，如不是,用k返回i前面有几个已知点

Ifcha="n"Then'不是已知点

n1=2*（i-k-1）