基于Matlab的TEQC绘图程序代码doc.doc

基于Matlab的TEQC绘图程序代码【实用文档】doc
文档可直接使用可编辑，欢迎下载
％在 ｍａtｌab下新建一个m文件,将以下代码直接拷贝进去，即可执行。
％需要一个ＴＥＱC生成的ｐlｏt文件作为参数
fｕｎctioｎ out=teｑcploｔ３（fileｓ);
%读取TEQＣ生成的Plot文件，绘制数据图表,支持Copｍact、Comｐａｃt2、Cｏmpａct３格式
％选取一个TＥQＣ的Ｐlｏt文件
% 　格式说明
％ ＊。ｓn1 　 载波 L1的信噪比 Signal to nｏise ratio (Ｓ/N)
%　　 ＊。ｓn2 载波　Ｌ２的信噪比 Ｓｉｇnal　to noｉｓe ratｉo （Ｓ/N) Caｒrier L2
％　 ＊．iod *.d12　＊．d21 电离层延迟观测值变化率(米／秒)　Ｄeｒivａｔivｅ of ioｎosphｅric dｅlay　oｂseｒvable （m／ｓ）
％　 ＊。ｉoｎ ＊。i１２ ＊.i2１　电离层延迟观测值（米)Ionosphｅrｉｃ delay　obｓeｒvabｌe (ｍ)
％ 　*。ｍp1　＊。m12 载波L１的多路径误差　 Mulｔiｐａtｈ Cａrrier L１
% *。mp2 *.m2１　　载波L2的多路径误差 　Mｕltipａth Cａrｒier L2
% *．ａzi 　　 卫星方位角(°） Ｓatelｌite azimｕtｈal data (degrees)
% ＊.ele　　　 卫星高度角(°） Ｓaｔｅｌliｔｅ ｅlevaｔｉoｎ　data （degrｅeｓ)
ｉｆ nargiｎ==0
［filen,pａｔｈ］=uｉgetｆiｌe（’*．ｓn１;*.sｎ２;＊.iod;＊．ion;*。ｍp1；＊。ｍp２；＊．azi；*。ele;＊．i１2;*i21；＊.m12；m21;*.d１2；*d21’，。..
　 　 ＇请选择TEＱC报告文件：＇)；
elｓｅ
　［ｐaｔｈｆｉlen,ext］=filｅｐaｒｔs（ｆiｌes）；
　 paｔh=[path ’\＇］；
　ｆilen＝{[fiｌeｎ exｔ］};
end
％读取这个文件
file＝char（filen);
％按行读取文件至数组A
［A]=impｏrｔdata(［ｐath　fｉle］,’\t'）;
％定义ＳAT,存放卫星数据
％GPS有３2颗卫星,存放序号1-3２，
％GLONASS有32颗卫星，存放序号33－64，
％BEIDOU有３5颗卫星,存放序号65—99
SAＴ（１:ｌength(A）,1:９9）=NａN;
％sats（1:leｎgth(A)，1：99)＝NaN；
%存放采样时间，单位秒
ｔｓec（1：lｅngｔh（A）)=NaN；
%读取文件的第一行
fｉlｅlx＝A{1};
％判断是哪种格式
sｗitch　fｉlelx
　 　 ｃaｓe ＇COＭPＡCＴ'
　　 %读取数据采样间隔
　 　 t_sａｍｐ=chaｒ(Ａ（3));
％读取开始时间
　 　 mjl=chａr（A(4）);
　 %读取数据采样间隔
　 　　　 Ｔ＿ＳAMＰ=str2nｕm（t_ｓamp（max（ｆiｎd（t＿sａmp＝=’ ’）)：eｎd））;
　 　 ％读取数据采样开始时间
　 　MJL_START=sｔr2ｎuｍ（ｍjl(mａｘ（fiｎd（mjl=＝＇ ')):end））;
　 　%转成时间序列数字,dａte serｉal　number,从0000年１月1日０时０分０秒开始计算的十进制天数
　　 MJD_ＳＴART=MJL_STARＴ+678941．999９９9７41；
　 　 %i为行号
　n＝1;ｉ=5；
　 casｅ ’ＣOMPACT2'
　 %读取数据采用间隔
　 　　 ｔ_sａmp=char（A(2）);
％读取开始时间
　 ｍjl＝cｈar(Ａ（3）)；
　 %读取数据采样间隔
　 　　　 T＿SＡMP＝str2ｎum（t＿sａmp(max（finｄ（t_ｓaｍp=＝＇ ’）)：end)）；
％读取数据采样开始时间
　 　 　ＭＪL＿START＝sｔｒ2nｕm(mｊl（maｘ（fｉｎd(mjｌ==’ ’）):ｅnd））；
　%转成时间序列数字，dａte　serial number，从０000年1月1日0时0分0秒开始计算的十进制天数
MJD_SＴART=MＪL_SＴART+6７8941。９9９999741；
　 　　n=１;i=4； 　 　
　 ｃase ’CＯMＰAＣT３＇
　 　 ％读取开始时间
　 　　 ｔ_start=char（A（２）)；
t_ｓtaｒt = deblank(ｔ_ｓｔａrt）;
　 　　 　 s = spｌiｔｓtr（ｔ_stａrｔ，＇＊＊＇，6）；

%t_sｔarｔ＿tｉme=［char（s{２})　＇年'　ｃhar（ｓ{3｝） ＇月’ chａｒ(s｛4｝) ＇日' ｃhar（s{5}） ＇时＇ chaｒ（s｛６｝） '分’ nuｍ2str（ｓtｒ2num（ｃhａｒ(s｛7}）)，＇％02d’) '秒’］ ;
　 　 　 %获取采样的开始时间，2０1３，12，7，0３，05，55
　　 　 t_s_ｔiｍe=[sｔr２nuｍ（char(s{2｝））， stｒ2num（ｃhar(ｓ{3｝））,ｓtr2num（char（s{4｝)），ｓｔr2ｎum(cｈar（s{５｝）),ｓtr2num（char(s{6｝))，str2nuｍ（ｃhar（s｛７｝）)]；
　 　 　n=１;i=3;
　 　 otｈerwiｓe
　disp（＇数据格式存在问题’）；
　 ｒetｕｒn
end
％ｎ=1;ｉ=３；
％sats＝stｒ２num（A｛i｝);
snｙggfｉｌen=sｔｒreｐ（fileｎ,’_’，’—');
%生成进度条
h=waiｔbar(0,[＇正在读取数据，请稍等……　’　chaｒ(snyｇｇｆiｌｅn)]);
ｗhilｅ　i〈lenｇth(A）;
waｉtbar（i/ｌeｎｇth（A)，ｈ）；drａwnow
　 %读取卫星编号数据行 采样时间 卫星数量 卫星编号 CＯＭAPＣT3
　　satbhｓｔr=char（A（ｉ）);
ｓ =　sｐlitsｔr（satbhstr，'*＊',３２）；
　 %如果卫星编号数据为0,说明数据缺失，退出本函数
　　if ｓtｒtｒｉm（saｔbｈstr)==’0’
　　 　dｉｓp ('数据为空!'）；
　 　　ｃloｓe(h）；
　 　 ｒｅturｎ
　eｎd
　switch ｆilelx
　 case　'COMPACT3＇
　 %记录采样时间,单位秒
　　　 　 tseｃ（ｎ)＝stｒ2nuｍ(s{1｝）;
　 ｉｆ s{2｝=＝'—1'
　 　 　　 ％如果为-1，使用上次的卫星编号行
　　 　 　　saｔbhstr＝oldsａtbｈ;
　　 　 　 　 ｓ = splitsｔr（ｓatbhsｔr,＇＊＊＇,32)；
　 elｓe
　 　olｄｓatbh=satｂｈsｔｒ；
　 eｎｄ
　 　 ％获取当前卫星数量
　 satcｏｕnt=ｓtr2num(s｛2｝）;
　 caｓｅ 　｛'ＣOMPＡCＴ’,＇COMPACT２’｝
　 ％记录采样时间,单位秒
　 　 　 ｔsec(ｎ）＝(ｎ-１）＊T_SAMP；　 　 　
　　　 　iｆ　ｓ｛1｝=＝'－1＇
　 　 ％如果为—1，使用上次的卫星编号行
　 　 sａtbｈｓtｒ=ｏｌdsａtｂh；
　 　 　 s　＝ splitstｒ（satbhstr，＇*＊’,32);
　 　 　　　 elｓe
　ｏlｄｓatｂh=sａｔbhstr；
　 　 　　enｄ
　 　　　 ％获取当前卫星数量
　　 　 　　　 satcounｔ=str２num（s｛1｝); 　 　 　 　
　　 end

　 　　 　%读取卫星数据行
　 　 satｄａtastr=char(A(i+１)）；
　sｄａta＝sｐｌitｓtr(sａｔdatastr,’＊＊’，32)；　
　 for　ｋ=3：ｓａtcount＋2
　　 sｗiｔch　filelx
　 　　 　 case　'COMPACT'
　　 　 　 ％如果是COMPACT格式，卫星编号前面没有字母，默认是GPS卫星，在前面增加字符G
　　　　　 satbh=[＇Ｇ' char（s{k—1｝）］;
　 　　 　　 　 ｃaｓｅ　’CＯMPACT2'
　 　 　 ％如果是COＭPACT2格式，卫星编号行第２个数据开始是卫星编号
　 　 　 satbh=char(s{k-１｝);
　 　 　 　 　　　　 　cａse ’COMPＡCＴ3＇
　 　 ％如果是COMＰACT３格式，卫星编号行第３个数据开始是卫星编号
　 　 　　 　 　 sａtｂh＝char（s｛k｝）；
　end

　　sｗitch satbｈ（1）;
　 　 　 cａsｅ ’Ｇ’
　 　 　 　%如果是GPS卫星，获得与编号对应的存储序号，范围1－３2
　 　 　 　 　iｎdｅｘ=ｓｔｒ２num（ｓatbｈ（２：3)）;
　 　 　 　cａse ＇R’
　 　 　 　 　 　 　 %如果是GLＯNASS卫星，获得与编号对应的存储序号，范围33-64
　 　　 　 　　 iｎｄｅx=３2+ｓtr2ｎum（ｓaｔbh（２：3))；
　 　　 　 　 cａse ’Ｃ’
　 　 　 　　　 　 　 %如果是ＢEIDOU卫星，获得与编号对应的存储序号，范围65－99
　 　 　 　 index=６4+sｔr2num（satbｈ（２：3）)；
　　　 　 　oｔhｅrｗise
　　 　 　 　　％其他情况,获得与编号对应的存储序号,范围１-３2
　　 index=ｓtｒ2num（ｓatbh(2：3));
　 eｎｄ %switch
　 　 　 　 ％获得对应的卫星数据
　 sdatasｔr=char（sdａｔa{k-２}）；
　 　 　％如果卫星数据不是数值型,用0替代
　 　 　　 if ｉsemｐtｙ（sｔr2num（sdaｔastr））==1
　 　 　　 　 ｓdａtastr＝＇0.０00';
　　 　 　ｅnd
%将卫星数据存入SＡT数组的对应位置
　 SAT（n，indｅx）=str2ｎum(ｓｄatastｒ);
　 　 ｅｎｄ　％foｒ
n=n+1；
　 %根据数据文件结构，一行是卫星编号,下一行就是对应的卫星数据
　 %程序一次读取2行数据
ｉ=i+2;
end ％whiｌe
％数据读取完毕，关闭进度条
ｃｌose（h）;
％获取已使用的卫星编号,返回给sａｔ_ｂｈ
sat_bｈ=getｓatbh(SAＴ);
％增加最后一个结束字符
sat_bh_temp＝［saｔ＿ｂｈ＇ ｛’Ｅnd’｝];
ｓａt＿bh＿eｎｄ＝sat_bh_temｐ＇；
%删去所有都为NaN值的列，返回给sａｔ_ｄata
sａt_dａtａ=ｄelnanｃol(SAT）;
%删去所有NａN值的采样时间列
t_seｃonds=deｌnancol（tsec);
％删去所有都为NaN值的行，返回saｔ_datas
ｓaｔ_datas＝delnanroｗ(sat_ｄata）；
％增加一个nａｎ列
[row,coｌ]=sizｅ（sat_datａs）；
bb（１:roｗ，1）＝ＮaN;
sat＿ｄatａs=［sat_datas bb]；
％将含有NａN值的数据替换为0,返回给s_ｄａｔa
%s_ｄata=repnaｎ20(sａｔ_data)；
%将采样时间转成顺序日期serｉａl date　ｎumber
switch　fiｌelｘ
　 　 ｃａse　｛'COＭPAＣT’,＇COMＰAＣT2’}
　 　 t_s_jd_time＝MJD_ＳＴＡRT；
case 'COMＰＡCT３＇
　 　　 　t＿s_jd_tiｍｅ=ｄatenｕｍ(ｔ＿s_ｔime）；
enｄ
for i＝1:ｌｅnｇth（t＿sｅconｄｓ）
　 ％将每次采样时间转换成相应的顺序日期
t_jd_time（i)=ｔ_s＿ｊd＿timｅ+t_secoｎds（i）/６0/60/２4；
ｅnｄ
%绘制图表+＋++++++＋+++++++++++++++++
％根据不同的文件类型，设置坐标轴的范围
［typｅ，maxy，minｙ］=get＿ｆilｅｔype（fiｌe）;
figure;ｂox on；hoｌｄ ｏn
［row，cｏl］=siｚe(sａt_daｔas）；
%绘制渐变彩色图
pcolor(t_jd_time＇,［1：col］，sat_ｄataｓ’)；
set(gcｆ，＇rｅnderｅｒ'，＇zbuffer＇）；
shadinｇ flat
set（gｃa,’xtiｃｋlａbel＇,［t＿jd_tｉme］）;
set(gcａ,'xlｉｍ',［t_jｄ_time(1）　t_jｄ_tｉｍe（end）］)
%ｔ_start_times＝[t_ｓtaｒt＿ｈ ’:’ t_start_ｍ　’：＇　t_ｓtarｔ＿s］；
dateaxis('Ｘ＇,１3)
cｂａr（＇v'，[minｙ　maxy］，tyｐe）;
seｔ(gca,＇ylｉｍ＇，［１ col＋1］）
set(gｃa，’yｔick',［1。5：1:col+０.5］）
ｓｅt（ｇca，＇ytｉｃklabel＇,sａt_bh_eｎd);
sｅt(gca,＇foｎｔsｉｚe'，７）；
colｏrｍap(fｌipud(jeｔ））；
ｃaxiｓ（［miny ｍaxy］);
％t_end_ｔime=cal2eｔ(t_ｓ_time,t＿secoｎds（ｅnd)）；
xlabeｌ（［dａtestr(t＿ｊd＿tｉme(１)） ' 　|-—－-－—-— 　采样间隔： ＇ num2ｓｔr(t_ｓecｏnds(2)—t_seconｄｓ(１)) ' 秒　　-－－—－———｜ 　 ＇ ｄaｔｅstr(t_jd_ｔｉｍe(eｎd））]）
ylabel(’卫星编号'）
％ｔimestr=seｃs2ｈms(lｅnｇｔh(saｔ));
Ｔ=ｔｉtlｅ（[’TEQC 报告文件: ' sｔｒrｅp（file，'_’，’-＇）])；
set(T，'fｏntsize',８）

%out.(ｆile(end-２:end)）＝sａt；
％ouｔ.Ｔ_sａmp=Ｔ_ＳAＭP；
out。Staｒt=daｔestr(t_jｄ_time（1））；
ｏut．Stｏp＝ｄａtestr（t_jd_ｔime（end）)；
％+＋+＋++＋＋++++++++＋+++＋+＋++++++＋+++＋＋++＋＋
%以下是本函数中所使用的一些相关函数

％+＋＋＋＋+++++＋+＋++＋+＋++＋+++++++++++＋+++++＋+++++＋+＋＋+＋
%+++++++++＋+++＋++++＋＋+＋＋++＋++++++++++＋++++++++++＋+
ｆuｎction s=ｓplｉtstr（stｒ，sｐlｉt,nuｍ);
%用指定分隔符分隔字符串函数
％str＝＇naｍe 0 　 2010-0５－06 ３3 Qt BKB'；
％'＊＊'表示特殊分隔符回车、空格、换行、制表符等
%nuｍ表示分隔字符串的个数
stｒiｎｇ=cｅll(num,1）；
inｄex=1；

teｍp＝［］；
iｆ ｓｐlｉt=＝’＊*'
%特殊分隔符
for ｉ=1：length（ｓtr）
　　　 iｆ　ｉsspacｅ(sｔr（i))
　 　 　if isempｔy（temｐ）==falｓe
　 　　 sｔring｛ｉndex｝＝temp;
　　　 indｅx＝ｉndeｘ＋１;
　 temp=[］;
　　 ｅnd
　 　 else
　 teｍｐ=［ｔemp,sｔｒ(i）］；
　　enｄ　　　
end
elｓｅ
%常规分隔符
　 foｒ i=１:length（stｒ）
if　str(i）==spｌit
　 　 if ｉｓempty(temp)==false
　 　 　 string｛indeｘ｝＝temp;
　 　 　 　inｄex=iｎdｅｘ+1;
　　 temp=［］;
　end
　 　 else
　　 temｐ=［teｍｐ，str（i）]；
　 enｄ
　 enｄ
eｎd

if iｓemｐｔy（ｔemp)==ｆalsｅ
ｓtring{index｝=tｅｍp；
eｎd
s＝ｓtrｉnｇ；
%++++++＋++＋+++＋＋++++++＋＋＋+++＋+＋+++++
functiｏn B=dｅlnancｏｌA）；
%删除矩阵中所有值都是nan的列
［x　y]=sｉze（A）；
fｏｒ i=y:-1:1
　 if alｌ（iｓnan（A（：,ｉ))）==１
　　 　A(:，i）=[];
enｄ
eｎd
B＝A;
％+++++＋＋+＋++++++＋+＋++++＋++++++++++++
funcｔion Ｂ=delnanrow(Ａ）；
％删除矩阵中所有值都是ｎan的行
[ｘ ｙ］=size（A）;
fｏr i＝ｘ：-1：1
　 if　all（isnaｎ(A（ｉ,:）））==1
　 　　Ａ（i，:)=［]；
end
end
B=Ａ;
%+++++++＋++＋+＋＋+＋++++++++＋＋＋+＋++＋+++

funcｔion　B=repｎan20(data);
%将矩阵中的NAN值用0替代
[ｄatａs，featurｅs]=size(data）;
foｒ k=1：ｆeatｕｒｅｓ
for i＝1:ｄataｓ
　 　　if　isnan（daｔa（ｉ,ｋ）)＝=1
　 　 data（i，k）=０；
eｎd
　 end
ｅnd
B=data;
％＋++＋+＋＋＋＋++＋+++++++＋++＋++++++＋++
fｕnction　satbh=getsatbｈ(satｄata);
%获取卫星数据中的卫星编号,返回一个字符串数组
s＝cｅll（32,１）;
teｍp=［］;
inｄex=１；
　 for i=1：99
ｉf ａll（isnａn(ｓatdata（:,i）））～＝1

　 switch ｉ>32
　 　 　 ｃaｓｅ ０
　　　 temp＝［＇G'，ｎum2str（i,’％02ｄ＇）]；
　 　 case　1
　 　　iｆ　i>６４
　 　 　 　　tｅmp=［’C＇，num２str(i-64，’％02d’）］；
　 else
　 　　　　 　teｍｐ=['R’,num２sｔｒ（i—３2，'％02d'）］；
　　 　ｅnｄ

　 eｎd
　 s｛ｉndｅｘ}=teｍｐ；
　　　 iｎdex=index+1；
eｎd
　 enｄ
iｆ ｉｎdeｘ〉1
　sａｔbｈ＝s（1:index－１);
eｌse
saｔｂｈ=［]；
ｅnd
％++＋＋＋++++＋+＋+++＋+＋++++++++++++＋+++++＋+++＋＋+++++++++＋＋+＋＋+++＋+++
％ FUNCTION： GＥＴ FＩLETYPＥ INFO +++++＋+++++++++＋++++++++＋+＋++＋+＋+++++

functiｏn [ｏｕt,ｍａxy，miny］=gｅt＿ｆiｌetｙpe(teｑfilｅ）；
％根据文件类型，设置坐标轴的范围
［paｔh，ｎａme,ext,veｒ］=fiｌepａrts（tｅqｆｉlｅ)；
sｗｉtｃh　exｔ
cａse ’。sｎ1'
%ｏuｔ＝'Signal　ｔo noiｓe ratｉo　Ｓ/Ｎ L1'；
　 ｏut=’L1载波上的信噪比(dBHz)'；
　 　 maｘy=60；
　 　 ｍiny＝15;　　　
　　ｃase '。ｓn2'
　 　　 %ｏut='Signal tｏ　ｎoiｓe ratiｏ Ｓ/N L2’；
　 out='Ｌ2载波上的信噪比（dBHz）’;
　 maxy=６0；
　 minｙ=１5；
cａse ｛’.mp1','。m12＇}
　　 　 ％ｏｕｔ=’Multiｐａｔh　Ｌ1’;
　 　out＝’Ｌ1多路径观测值'；
　 mａxy=１；
　 　 miny＝－１；
caｓｅ ｛’．ｍp2’，’。m21'｝
　 　 %out='Multipath Ｌ2'；
　 oｕt=’Ｌ2多路径观测值＇；
　 　 mａxy=１;
　 miny=—1；
cａse {＇.ｉｏd’，'．d12',＇d21'｝
　　 　 ％oｕｔ=＇Ｄerivative of　ionｏsｐherｉc dｅlay obｓｅrvable （m／s）'；
　 　 ouｔ=’电离层延迟观测值变化率（米/秒)’；
　 mａxy=１;
　　　　 　 mｉny=-１；
　casｅ ｛＇。iｏn’,'.i１2’,'i21＇｝
　　 　　　 maxy＝2；
　　 　miｎy=-２;
　 　％out=’Iｏnｏsｐheｒｉc　dｅlaｙ　ｏbserｖable （m）';
　 　out=’电离层延迟观测值（米)’；
case '。eｌe’
　　 　 maxy=９０；
　　　　 　 miny=0；
　 　 ％ｏut='Sａtelｌｉｔｅ　ｅlevatｉoｎ daｔa’;
　　 　 oｕt=’卫星高度角（°)'；
　 　 ｃase ’.aｚｉ'
　　 maxy=１８０；
　 ｍinｙ=-180；
　 ％out＝’Sａtｅllｉtｅ ａzｉmuthal　ｄａta';
　　 　 out='卫星方位角（°)＇;
otｈｅrwiｓe
disp（'Somethings wrｏng．．！'）
ｅnd

％ FUNCTIＯＮ: PLACE　A ＭODIＦIED COLORBAR　＋+＋＋+++++++＋+++＋++＋＋++＋++++＋++

ｆuｎction CB＝ｃｂar(ｌｏc，ranｇe，label);

% 。。.．．．。．。.。.。.。。..．．。.。..。．.．。.．。.。.。.。.。．..。..。。.．。.。。．．。...
% CＢ =　cｂaｒ(loc,rａngｅ，lａbel）
％ 　 ｐｌaceｓ a　coｌoｒbar ａｔ：
％ ｌoc = 'v＇ iｎ vertｉcal or 'h＇ in horiｚontaｌ
％ 　 　 　 poｓiｔion in cｕrrent fｉgure scaled betweｅn：
%　　 rａnge = [min max］　wｉth a:
％ labｅｌ = ’stｒiｎg'。
%
％ fonｔｓize iｓ ｒｅdｕcｅd to 10 aｎd　wｉdth oｆ　baｒ iｓ ｈalf dｅfauｌt.
% 　
％ 　Examplｅ：　　 　[Ｘ，Y，Z］＝peaks２5）；
％ 　 　 　 　　 ranｇe＝[mｉn（miｎ（Z)) maｘ（mａx(Z））]；
%　 　 　 ｐcolor(X，Y,Ｚ)；
％ 　 　 cbａr（'ｖ＇,rａngｅ，'Elevａｔｉｏn （m）＇）
% .．.．。．．.．.。。．。.。．。.．。。．..。。。.．.．。。.．....。。。..。。。．．。。．．。。。。。.。

caxiｓ（[rangｅ(1）　rangｅ（２）］)；
switch lｏc
　 　cａse 'v＇
　 　 　CＢ＝ｃolorｂar（'ｖertiｃａｌ’)；
　 ｓet(ＣＢ,’ｙlｉm',［range(1) range（2）］）;
　 　POＳ=gｅt（ＣB，'pｏｓiｔion’)；
　 　　 sｅt（CB,'poｓition＇,[ＰOS(１) POS（2) 0。03 ＰＯＳ（４)］）；
　 　　 　sｅｔ（ＣB，’fｏｎｔsize’,8）;
　　 　set（geｔ(CB,’ylabel’）,’string',laｂeｌ）；
　 　ｓet（CB，'box＇，＇on’）

　 　 cａse 'h'
　　 ＣB＝colｏrbar('horizonｔal'）；
　 　 　ｓｅｔ（CB，＇xｌｉm'，［rangｅ(1）　range(2）］)；
　 　POS＝get（CB,’positiｏn＇）；
　 set(CB，'pｏｓｉtion’,[ＰＯＳ(1) POS（2) ＰOS(３）　0。０3］);
　 set（ＣＢ，'ｆontsｉze’，８);
　 　 sｅt(get（CB，'xlａbel’），＇striｎg＇，ｌabel）
eｎd
％+＋+＋+＋++++++++++++++++++＋＋++++++++＋＋＋++++++＋++＋+++＋++++＋++
funｃｔiｏn ｊd=ｃaｌ2jd（ｃal）
% ｃａl2jd1ﻩ将公历年月日时分秒转换到儒略日。
％ jd=caｌ2jd（caｌ） 　返回儒略日
% 　cａｌ:1ｘ６矩阵，6列分别为年月日时分秒。构造cal时可以省略末尾的0
％
% 公元１５８2年10月4日24：00点之前使用儒略历,公元１5８２年10月15日00：00点之后使用公历

iｆ lengｔh（cal） 〈 6
ｃal(6)=0；
end
ｙeaｒ=ｃal（1）;
ｍoｎth=caｌ(２)；
ｄay=ｃａl(3)＋(caｌ（4)*36０0+cal(5)＊6０＋ｃａl（６））/86４00;

ｙ ＝ yｅaｒ + 4８00； ％4８０1 Ｂ．C. iｓ　a ceｎtury　yeａr and also a　ｌeap yｅar.

if( ｙeａr 〈　０　)
y ＝y+　1；　ﻩ% Ｐｌeａｓｅ ｎｏｔe that　thｅre iｓ nｏ　yeaｒ ０ A。Ｄ。
end

ｍ＝ｍoｎｔh;
ｉf( m 〈= ２ ) ％ Jａnuaｒy ａnd Februarｙ ｃome　afｔer　Decembｅr。
m ＝ m＋１2；
ﻩy =　y —　1;
end

e=fｌoor（３0.6 * (m+1））;

a＝fｌｏor(y/１0０）； ％ nｕｍｂeｒ　ｏf　ceｎturｉeｓ

％ 教皇格雷戈里十三世于15８2年2月２４日以教皇训令颁布,将1582年10月5日至14抹掉。1５８２年１０月４日过完后第二天是10月15日
ｉｆ( year　＜1582 ）｜（yｅar==1５82＆mｏｎth<10）｜(year=＝1582＆ ｍonｔh==１0　＆ｄaｙ〈15)
b ＝ -3８；
elｓe
ﻩb = fｌｏｏr(（ａ／4） - a）; % nｕmber of ｃeｎｔｕrｙ yearｓ　thaｔ aｒe not lｅap　years
ｅnd

c=fｌoor(３６5.2５*　y)； % Juｌian calendar　yeaｒｓ ａｎd ｌｅap ｙｅarｓ

jd= b + ｃ +　ｅ + day -　３2１67．5;
％＋＋+＋+++++＋+++++++++＋＋++++++＋++＋++++++++++++++++++
functｉｏn cａl=ｊd2cａl（Ｊ）
% 从儒略日计算公历年月日时分秒
％ caｌ＝ｊd2cal（J）
％ 返回的cal是１ｘ6矩阵，6列分别为年月日时分秒
％
%　公元15８2年１０月4日24:00点之前使用儒略历，公元１5８2年１０月1５日00：00点之后使用公历

if　(J 〈 1721423．5） ％ 公元1月１日0时
BＣ = 1；
elsｅ
BC = 0；
enｄ

％stａrt fｒom Julian Ｍａrcｈ　１, 4801 Ｂ。C.

if(　Ｊ <　２29９160.5　） ％ beforｅ 1５8２.１0.4。　24:00　iｓ　Juｌiａn calender
ﻩｊ0＝fｌoor(J＋０．5）;
ﻩdd＝Ｊ+０．5—j0；
elｓe % aｆｔeｒ 1582。10．1５。 00:00 is Greｇoriａn　ｃaleｎder
％number ｏｆ certurｙ yearｓ ｔhaｔ aｒe ｎot　leap yeaｒ
n１=ｆｌooｒ（（J—23420３1．5）/36524。2５/4）+1; ％17００。３。1.０
n2＝fｌoｏr((J—2３78555.5)/36524.25／4)+1; ％1８０0．3。1。０
n3=ｆloor((J—241５0７9.５)／３65２4．２5/4）+１; ％1900。３。１.０
ﻩｊ0=n1+n2＋n3+Ｊ+10;
dd=ｊ０+0。5—ｆloor（j0＋0。５）;
j0=fｌoor(j0＋０．5)；
ｅｎd

j0=j0+３２08３;
year０＝ceil（j0／365．2５)—1;
year=ｙear0-４8００;
day＝j0—flooｒ(yeａr0＊36５．25）；
ｍonth＝ｆloor（(day－０.6)/30．６)＋3;
dａｙ=dａｙ－ｒound（(mｏｎth—3)*30．６）;
if ｍｏnth〉１2
ﻩmonth=month—1２;
yeａr＝yｅａr+1;
enｄ

ｙear=yeａr—BC；

seｃ=ｒouｎd(dd＊8６400);

hour=ｆｌoｏｒ（sec／360０）；
sec＝sec—hoｕr*36０0;
min＝floor(sｅc/６0)；
sec=sec-ｍｉn＊６0;
ｃal＝[year，　moｎth, dａy, ｈoｕr， ｍin， ｓec］；
%++++＋+＋+＋++++＋++++＋+++++＋++++++＋+＋+＋++＋+++++++++＋
ｆunctiｏn mｊd＝cal2ｍjd(caｌ）
％　cａｌ2mjdﻩ将公历年月日时分秒转换到简化儒略日。
％　　ｍｊｄ＝cal2mｊd(ｃal） 返回简化儒略日
% cal：1x6矩阵，６列分别为年月日时分秒。构造cａｌ时可以省略末尾的0

jｄ　= cal2jd（cal）；
mjｄ =　jd - 2４0００００.5;
％ ＦUNCTION： SECONDＳ　TO HOURＳ， MIＮUTES ａnd SＥCONDS +＋++＋+++++＋+＋+++

function　tiｍｅstｒ＝sｅcs2hms(SECS)

HOＵRS=SEＣS/60/６０;
ｈours=flｏor（ＨOUＲS）；
MINUTES=HＯＵRS－hｏｕｒs）＊60;
minuｔｅs=flｏor（MINＵTES）;
seconｄs=（MINUTＥＳ-minuｔes）*６0;
ＨH=ｎum2str(hourｓ）；
ＭM=nuｍ２str（ｍiｎｕｔｅs）;
SＳ＝ｎum2ｓtr（seconds)；

if secｏｎds＜10;
　 SＳ=［'0'　num2str（SS)］；
else
　　 SＳ=nuｍ2ｓtｒ(SS);
eｎd
if mｉnuteｓ〈10；
　 MM=［＇０’ nｕm2ｓｔｒ(MM）]；
ｅlse
MM=num２stｒ(MＭ);
end
tｉｍestｒ=［HH　'：’ ＭM '：' ＳS］；
％++++++＋+++＋++++＋＋＋++＋+++++++＋+＋＋＋+＋++
ｆｕnctｉon ｔimｅstr＝cal2et（t＿start，secs）
%根据开始时间，时长（秒),计算结束时间
％ｔ_stａrｔ 1×6矩阵，年,月,日,时,分，秒
%secs　时长，单位秒
t_ｓ_days＝t_sｔart(4)/24+ｔ_sｔart(5)/60／24+ｔ_start(6)／60/６0/24;
t_e_dayｓ=t＿s_ｄays＋secs/６0/６0/24；
HOURＳ＝t_e_days*2４；
hours=floｏｒ（HOＵRS）；