基于Matlab的TEQC绘图程序代码

%在matlab下新建一个m文件，将以下代码直接拷贝进去，即可执行。%需要一个TEQC生成的plot文件作为参数functionout=teqcplot3(files);%读取TEQC生成的Plot文件，绘制数据图表，支持Copmact、Compact2、Compact3格式%选取一个TEQC的Plot文件%格式说明%*.sn1载波L1的信噪比Signaltonoiseratio(S/N)%*.sn2载波L2的信噪比Signaltonoiseratio(S/N)CarrierL2%*.iod*.d12*.d21电离层延迟观测值变化率(米/秒)Derivativeofionosphericdelayobservable(m/s)%*.ion*.i12*.i21电离层延迟观测值(米)Ionosphericdelayobservable(m)%*.mp1*.m12载波L1的多路径误差MultipathCarrierL1%*.mp2*.m21载波L2的多路径误差MultipathCarrierL2%*.azi卫星方位角（°）Satelliteazimuthaldata(degrees)%*.ele卫星高度角（°）Satelliteelevationdata(degrees)ifnargin==0[filen,path]=uigetfile('*.sn1;*.sn2;*.iod;*.ion;*.mp1;*.mp2;*.azi;*.ele;*.i12;*i21;*.m12;m21;*.d12;*d21',...'请选择TEQC报告文件：');else[path,filen,ext]=fileparts(files);path=[path'\'];filen={[filenext]};end%读取这个文件file=char(filen);%按行读取文件至数组A[A]=importdata([pathfile],'\t');%定义SAT，存放卫星数据%GPS有32颗卫星，存放序号1-32，%GLONASS有32颗卫星，存放序号33-64，%BEIDOU有35颗卫星，存放序号65-99SAT(1:length(A),1:99)=NaN;%sats(1:length(A),1:99)=NaN;%存放采样时间，单位秒tsec(1:length(A))=NaN;%读取文件的第一行filelx=A{1};%判断是哪种格式switchfilelxcase'COMPACT'%读取数据采样间隔t_samp=char(A(3));%读取开始时间mjl=char(A(4));%读取数据采样间隔T_SAMP=str2num(t_samp(max(find(t_samp=='')):end));%读取数据采样开始时间MJL_START=str2num(mjl(max(find(mjl=='')):end));%转成时间序列数字，dateserialnumber，从0000年1月1日0时0分0秒开始计算的十进制天数MJD_START=MJL_START+678941.999999741;%i为行号n=1;i=5;case'COMPACT2'%读取数据采用间隔t_samp=char(A(2));%读取开始时间mjl=char(A(3));%读取数据采样间隔T_SAMP=str2num(t_samp(max(find(t_samp=='')):end));%读取数据采样开始时间MJL_START=str2num(mjl(max(find(mjl=='')):end));%转成时间序列数字，dateserialnumber，从0000年1月1日0时0分0秒开始计算的十进制天数MJD_START=MJL_START+678941.999999741;n=1;i=4;case'COMPACT3'%读取开始时间t_start=char(A(2));t_start=deblank(t_start);s=splitstr(t_start,'**',6);%t_start_time=[char(s{2})'年'char(s{3})'月'char(s{4})'日'char(s{5})'时'char(s{6})'分'num2str(str2num(char(s{7})),'%02d')'秒'];%获取采样的开始时间，2013,12，7,03,05,55t_s_time=[str2num(char(s{2})),str2num(char(s{3})),str2num(char(s{4})),str2num(char(s{5})),str2num(char(s{6})),str2num(char(s{7}))];n=1;i=3;otherwisedisp('数据格式存在问题');returnend%n=1;i=3;%sats=str2num(A{i});snyggfilen=strrep(filen,'_','-');%生成进度条h=waitbar(0,['正在读取数据，请稍等……'char(snyggfilen)]);whileilength(A);waitbar(i/length(A),h);drawnow%读取卫星编号数据行采样时间卫星数量卫星编号COMAPCT3satbhstr=char(A(i));s=splitstr(satbhstr,'**',32);%如果卫星编号数据为0，说明数据缺失，退出本函数ifstrtrim(satbhstr)=='0'disp('数据为空！');close(h);returnendswitchfilelxcase'COMPACT3'%记录采样时间，单位秒tsec(n)=str2num(s{1});ifs{2}=='-1'%如果为-1，使用上次的卫星编号行satbhstr=oldsatbh;s=splitstr(satbhstr,'**',32);elseoldsatbh=satbhstr;end%获取当前卫星数量satcount=str2num(s{2});case{'COMPACT','COMPACT2'}%记录采样时间，单位秒tsec(n)=(n-1)*T_SAMP;ifs{1}=='-1'%如果为-1，使用上次的卫星编号行satbhstr=oldsatbh;s=splitstr(satbhstr,'**',32);elseoldsatbh=satbhstr;end%获取当前卫星数量satcount=str2num(s{1});end%读取卫星数据行satdatastr=char(A(i+1));sdata=splitstr(satdatastr,'**',32);fork=3:satcount+2switchfilelxcase'COMPACT'%如果是COMPACT格式，卫星编号前面没有字母，默认是GPS卫星，在前面增加字符Gsatbh=['G'char(s{k-1})];case'COMPACT2'%如果是COMPACT2格式，卫星编号行第2个数据开始是卫星编号satbh=char(s{k-1});case'COMPACT3'%如果是COMPACT3格式，卫星编号行第3个数据开始是卫星编号satbh=char(s{k});endswitchsatbh(1);case'G'%如果是GPS卫星，获得与编号对应的存储序号，范围1-32index=str2num(satbh(2:3));case'R'%如果是GLONASS卫星，获得与编号对应的存储序号，范围33-64index=32+str2num(satbh(2:3));case'C'%如果是BEIDOU卫星，获得与编号对应的存储序号，范围65-99index=64+str2num(satbh(2:3));otherwise%其他情况，获得与编号对应的存储序号，范围1-32index=str2num(satbh(2:3));end%switch%获得对应的卫星数据sdatastr=char(sdata{k-2});%如果卫星数据不是数值型，用0替代ifisempty(str2num(sdatastr))==1sdatastr='0.000';end%将卫星数据存入SAT数组的对应位置SAT(n,index)=str2num(sdatastr);end%forn=n+1;%根据数据文件结构，一行是卫星编号，下一行就是对应的卫星数据%程序一次读取2行数据i=i+2;end%while%数据读取完毕，关闭进度条close(h);%获取已使用的卫星编号，返回给sat_bhsat_bh=getsatbh(SAT);%增加最后一个结束字符sat_bh_temp=[sat_bh'{'End'}];sat_bh_end=sat_bh_temp';%删去所有都为NaN值的列，返回给sat_datasat_data=delnancol(SAT);%删去所有NaN值的采样时间列t_seconds=delnancol(tsec);%删去所有都为NaN值的行，返回sat_datassat_datas=delnanrow(sat_data);%增加一个nan列[row,col]=size(sat_datas);bb(1:row,1)=NaN;sat_datas=[sat_datasbb];%将含有NaN值的数据替换为0，返回给s_data%s_data=repnan20(sat_data);%将采样时间转成顺序日期serialdatenumberswitchfilelxcase{'COMPACT','COMPACT2'}t_s_jd_time=MJD_START;case'COMPACT3't_s_jd_time=datenum(t_s_time);endfori=1:length(t_seconds)%将每次采样时间转换成相应的顺序日期t_jd_time(i)=t_s_jd_time+t_seconds(i)/60/60/24;end%绘制图表++++++++++++++++++++++++++%根据不同的文件类型，设置坐标轴的范围[type,maxy,miny]=get_filetype(file);figure;boxon;holdon[row,col]=size(sat_datas);%绘制渐变彩色图pcolor(t_jd_time',[1:col],sat_datas');set(gcf,'renderer','zbuffer');shadingflatset(gca,'xticklabel',[t_jd_time]);set(gca,'xlim',[t_jd_time(1)t_jd_time(end)])%t_start_times=[t_start_h':'t_start_m':'t_start_s];dateaxis('X',13)cbar('v',[minymaxy],type);set(gca,'ylim',[1col+1])set(gca,'ytick',[1.5:1:col+0.5])set(gca,'yticklabel',sat_bh_end);set(gca,'fontsize',7);colormap(flipud(jet));caxis([minymaxy]);%t_end_time=cal2et(t_s_time,t_seconds(end));xlabel([datestr(t_jd_time(1))'|--------采样间隔:'num2str(t_seconds(2)-t_seconds(1))'秒--------|'datestr(t_jd_time(end))])ylabel('卫星编号')%timestr=secs2hms(length(sat));T=title(['TEQC报告文件:'strrep(file,'_','-')]);set(T,'fontsize',8)%out.(file(end-2:end))=sat;%out.T_samp=T_SAM