一、准备工作

观测值数据读取是进行数据处理的前提，通常，观测值的数据格式有rtcm、ubx、rinex，各家厂商还有自定义的格式。数据读取是最简单的工作，掌握了数据组织格式即掌握了数据读取策略，简言之，我们是对字符串进行一系列操作。我认为对o文件读取要做以下准备工作，分别是：

（1）安装Matlab的PC;

（2)  学习rinex.pdf资料；

（3）o文件准备；

（4）伪代码编写，理清数据流；

（5）代码测试

我们在使用定位方法SPP/PPP/RTK时，是按照历元来读取卫星的观测值（PR/CP），Doppler，SNR等信息，后续的处理是对所有卫星数进行循环遍历，进行不同的操作。

二、需解决的问题

问题是：要实现使用Matlab绘制不同卫星系统的卫星在观测弧段的伪距观测值时序图。

在上述问题分析中，易知卫星PRN号是唯一的，可作为key；而卫星对应着观测值（PC/CP/Doppler），可作为value。换言之，key和value可以使用map的映射来作为数据结构，其中，value可以使用结构体保存。最重要的一环已经确定了，即使用map和struct。

接下来，我们拿一颗卫星的观测弧段内的伪距值来说，可能的情况有以下几种：

（1） 在观测弧段内，卫星一直被捕获，即100%存在伪距值；

（2） 在观测弧段内，卫星在开始历元没有被捕获，而后才被捕获；

（3）在（2）中还存在未捕获，捕获，未捕获，捕获等情况；

在未捕获的历元时观测值为0，同时，我们还要注意，有伪距不一定有载波，没有载波肯定没有伪距。我们的伪代码主要是对上述三种情况进行判断即可。

三、代码实现

具体问题具体分析，解决上述提到的问题，具体代码如下：
% function: read the body of rinex obs file
clear all
close all
clc

fid=fopen(' example.obs');
if fid==-1
    disp('fail to open obs file!');
end

% 创建结构体obs，使用Map进行映射，prn为键(char)，sv为值(table)
obs=containers.Map;
record_counter=0;

% 读取body记录
while feof(fid)==0
    line=fgets(fid);
    if strcmp(line(1),'>')
        sv_sum=str2double(line(34:35));
        time=line(3:29);
        if sv_sum>0
            for i=1:sv_sum
                line=fgets(fid);
                % 键 -> prn
                key=line(1:3);
                % 值 -> sv(pr cp dop)
                sv.pr=str2double(line(6:17));
                cp=line(21:33);
                if isnan(cp)
                    sv.cp=0;
                else
                    sv.cp=str2double(line(21:33));
                end
                sv.dop=str2double(line(41:49));
                % 值逻辑判断
                if isKey(obs,key)
                    jump_sec=record_counter-length(obs(key));
                    if jump_sec~=0 % 解决掉星后重新锁定的问题
                        temp.pr=0;
                        temp.cp=0;
                        temp.dop=0;
                        tmp_strcut=repmat(temp,1,jump_sec);
                        value=[obs(key) tmp_strcut sv];
                    else % 卫星一直锁定，直接追加
                        if strcmp(key,'C18')==1
                            if length(obs(key))>5713
                                disp(time);
                            end
                        end
                        value=[repmat(obs(key),1) repmat(sv,1)]; % mark
                    end
                else % 卫星首次锁定
                    if record_counter~=0 % 第一秒未锁定
                        temp.pr=0;
                        temp.cp=0;
                        temp.dop=0;
                        tmp_strcut=repmat(temp,1,record_counter);
                        value=[tmp_strcut sv];
                    else % 第一秒锁定
                        value=sv;
                    end
                end
                obs(key)=value;
            end
        end
        record_counter=record_counter+1;
    end
end
disp('create struct obs(body) successfully!');

% 绘图
sv_sum=obs.Count; % 可视卫星数量
keySet=keys(obs); % 可视卫星名称
disp(keySet);

v=299792458.0;
bf1=1.561098E9; % bds L1 频率
gf1=1.57542E9; % gps L1 频率
gps_sum=0;bds_sum=0;

for i=1:sv_sum
    prn=char(keySet(i));
    if strcmp(prn(1),'G')==1 % GPS
        gps_sum=gps_sum+1;
    else if strcmp(prn(1),'C')==1 % BDS
            bds_sum=bds_sum+1;
        end
    end
end
subfig1=ceil(sqrt(bds_sum));
subfig2=ceil(sqrt(gps_sum));

x=1:record_counter;
for i=1:sv_sum
    prn=char(keySet(i));
    info=obs(prn);
    pr=zeros(1,record_counter);cp=pr;
    for j=1:length(info)
        pr(j)=info(j).pr;
        cp(j)=info(j).cp;
    end
    if i==1 || i==bds_sum+1
        figure;
    end
    subplot(n,n,pos);
    plot(x,pr,'r',x,cp*f,'g'); 
    title(char(prn));
end

四、示例结果

我只针对工作中遇到的rinex版本数据进行处理，gps/bds双系统，单频数据，结果图如下：




示例代码比较简单，也有相应的注释。

如题，这里对数据的处理对于熟悉matlab各函数操作的同学来讲，是常规操作，但对于matlab新手或只是想借用matlab工具将存于PC端文件中的数据通过串口发送出去的工程人员来说，由于对matlab软件比较陌生，实现各种功能所应该调用的函数没有头绪，所以还是比较麻烦的。这里先简单介绍实现这项操作过程中的几个步骤，然后贴出完整代码，供大家需要时查询。
这里主要以功能分以下步骤：
matlab读取指定文件数据
不同类型数据转换为字节数据
固定时间间隔的实现
操作串口发送数据
matlab读取指定文件数据
这里介绍两种读取方式：选择读取和指定读取
选择读取
所谓选择读取即在m文件执行过程中由操作者自己选择要读取的文件，先上代码
clear all%要操作的数据行数(从首行算起)sym num ;num=693;%定义1列，num行数组变量data_lat=zeros(1,num);%选择读取文件数据[FileName PathName]=uigetfile()file=fullfile(PathName,FileName)M=csvread(file);%单独操作M矩阵的第一列1~num行数据data_lat=M(1:num,1);%画出读取的纬度数据，红色plot(data_lat,'r');
上述代码的最终实现是选择文件读取其数据到M矩阵，并将数据的第一列的1~num行以红色画出图来。
代码执行后首先会弹出选择窗口，手动选定要操作的文件
图 1 
注意这里可以自由选择不同路径下的文件，要操作的文件不一定在当前工作目录下。选择窗口的打开源于uigetfile()函数的执行，选择文件结束后函数返回被选择文件的路径和文件名
图 2
而后，fullfile()合成被选择文件的完整路径，并返回给file变量，最后通过csvread()函数将文件数据读取进来并形成M数据矩阵，如下原csv文件与M数据对比，原数据类型依次为double/double/float/float/float/float
图 3
这里采用csvread()函数读取数据是因为笔者的原数据文件为.csv格式，需根据文件不同格式选择读取函数，如xlsread()函数
此外，matlab中数组和C语言的不同在于，没有0号元素，即没有0号下标的元素，都是从1号元素开始
指定读取
指定读取即在m文件中直接指定要读取文件的完整路径，在程序执行过程中无需再选择文件
clear all%要操作的数据行数(从首行算起)sym num ;num=693;%定义1列，num行数组变量data_lat=zeros(1,num);%选择读取文件数据% [FileName PathName]=uigetfile()% file=fullfile(PathName,FileName)% M=csvread(file);%指定读取文件数据M=csvread('D:studymatlabuart3926_vehicle_global_position_0.csv');%单独操作M矩阵的第一列1~num行数据data_lat=M(1:num,1);%画出读取的纬度数据，红色plot(data_lat,'r');
不同类型数据转换为字节数据
串口常规发送8位字节数据，这里不同类型数据只有按计算机存储规则先转换为字节发送给接收端，接收端再按相同的存储规则将字节数据转换为原数据类型以完成数据的传输，这里首先存在一个数据转换的问题
%整型inti2ubuff=zeros(1,4);i2ubuff=typecast(int(12345),'uint8');%单精度浮点型floatf2ubuff=zeros(1,4);f2ubuff=typecast(single(data_lon(i)),'uint8');%双精度浮点型doubled2ubuff=zeros(1,8);d2ubuff=typecast(double(data_lat(i)),'uint8');
‘=’左侧即为不同数据类型返回的无符号字节数组数据，这里还有要注意的一项是发送端和接收端的数据存储模式，也即大端、小端问题
固定时间间隔的实现
数据准备好之后，有时需要一股脑的将数据全部发送出去，有时又需要以固定的周期，定时一部分一部分的发送出去。后者的定时作用，在单片机上一般由定时器timer或滴答定时器等实现精准的时间把控，在matlab中则需要调用函数来实现
整秒定时
对于整秒时间间隔可以调用pause()函数，如定时3秒调用 pause(3) 即可实现
压秒级定时
当要实现的时间间隔单位低于秒级时，pause()函数就无能为力了，tic toc函数可以实现
% 延时200msticwhile toc < 0.2    end
tic执行时，会获取实时系统时间t0，toc则一直获取当前系统时间t1与t0的时间差，如上，当时间间隔大于或等于0.2时，程序会跳出while循环向下执行，否则会在while中等待，注意这种延时存在一定的精度误差
操作串口发送数据
由于只执行简单的串口操作，和单片机控制串口的操作很相似
指定端口，创建端口设备
s=serial('COM18');
设置串口属性
set(s,'BaudRate',115200,'parity','none','databits',8,'stopbits',1);
属性名字，紧跟要设置成的属性
打开串口
fopen(s);
发送数据
fwrite(s,sendbuff);
关闭串口，删除内存中串口设备对象，清除工作空间中串口设备对象
fclose(s);delete(s);clear s;
有时打开串口后，程序由于某种原因未能执行完毕，中途退出后，再次运行程序会显示打开串口错误，要打开的串口被占用，这是因为上次运行出错的程序未能正确关闭串口。此时有两种选择，一是关闭matlab软件，解除占用，但是打开matlab并等待就绪会占用很长一段时间，挺让人不耐烦的。二是在命令窗口运行下列函数，解除占用
%强行关闭串口scoms=instrfindall;stopasync(scoms);fclose(scoms);delete(scoms);
完整代码
最后贴出笔者的完整代码，该代码的实现是读取.csv数据以5hz的频率将数据按固定格式发送出去
%串口发送程序%测试数据在excel表中列顺序：经度、纬度、高度、y轴(东)速度、x轴(北)速度、z轴(下)速度clear all%操作的每列数据行数，从第一行算起sym numnum=693;%中间数据global sendbuffglobal f2ubuffglobal d2ubuffglobal data_latglobal data_longlobal data_altglobal data_vyglobal data_vxglobal data_vzsendbuff=zeros(1,48);f2ubuff=zeros(1,4);d2ubuff=zeros(1,8);data_lat=zeros(1,num);data_lon=zeros(1,num);data_alt=zeros(1,num);data_vy=zeros(1,num);data_vx=zeros(1,num);data_vz=zeros(1,num);%获取日志数据[FileName PathName]=uigetfile();file=fullfile(PathName,FileName);M=csvread(file);data_lat=M(1:num,1);data_lon=M(1:num,2);data_alt=M(1:num,3);data_vy=M(1:num,4);data_vx=M(1:num,5);data_vz=M(1:num,6);%选定串口，创建串口设备对象s=serial('COM18');fclose(s);%设置串口属性set(s,'BaudRate',115200,'parity','none','databits',8,'stopbits',1); %打开串口fopen(s);%固定字节sendbuff(1)='$';sendbuff(2)='E';sendbuff(3)='X';sendbuff(4)='Y';sendbuff(5)='F';sendbuff(6)=48;sendbuff(7)=0;sendbuff(8)=20;sendbuff(9)=20;sendbuff(42)=1;sendbuff(43)=0;sendbuff(44)=0;sendbuff(45)=160;sendbuff(46)=64;sendbuff(47)=0;sendbuff(48)=63;for i=1:num%纬度d2ubuff=typecast(double(data_lat(i)),'uint8');sendbuff(10)=d2ubuff(1);sendbuff(11)=d2ubuff(2);sendbuff(12)=d2ubuff(3);sendbuff(13)=d2ubuff(4);sendbuff(14)=d2ubuff(5);sendbuff(15)=d2ubuff(6);sendbuff(16)=d2ubuff(7);sendbuff(17)=d2ubuff(8);%经度d2ubuff=typecast(double(data_lon(i)),'uint8');sendbuff(18)=d2ubuff(1);sendbuff(19)=d2ubuff(2);sendbuff(20)=d2ubuff(3);sendbuff(21)=d2ubuff(4);sendbuff(22)=d2ubuff(5);sendbuff(23)=d2ubuff(6);sendbuff(24)=d2ubuff(7);sendbuff(25)=d2ubuff(8);%高度f2ubuff=typecast(single(data_lon(i)),'uint8');sendbuff(26)=f2ubuff(1);sendbuff(27)=f2ubuff(2);sendbuff(28)=f2ubuff(3);sendbuff(29)=f2ubuff(4);%vyf2ubuff=typecast(single(data_vy(i)),'uint8');sendbuff(30)=f2ubuff(1);sendbuff(31)=f2ubuff(2);sendbuff(32)=f2ubuff(3);sendbuff(33)=f2ubuff(4);%vxf2ubuff=typecast(single(data_vx(i)),'uint8');sendbuff(34)=f2ubuff(1);sendbuff(35)=f2ubuff(2);sendbuff(36)=f2ubuff(3);sendbuff(37)=f2ubuff(4);%vzf2ubuff=typecast(single(data_vz(i)),'uint8');sendbuff(38)=f2ubuff(1);sendbuff(39)=f2ubuff(2);sendbuff(40)=f2ubuff(3);sendbuff(41)=f2ubuff(4);fwrite(s,sendbuff);% 延时200msticwhile toc < 0.2    endend%延时一段时间等待数据发送完毕再关闭串口，数据未发送完毕时，不能正常关闭串口pause(3);%关闭串口fclose(s);%删除内存中串口设备对象delete(s);%清除工作空间中串口设备对象clear s;

嵌入式系统通过串口向PC发送的数据往往是16进制数，但是我们直观的比较数据则是用10进制的形式，这就需要将16进制数转化成10进制数。
嵌入式系统通过串口调试助手将采集的数据发送至PC后，一般先保存成TXT文档。此时的第一个问题就是如何通过Matlab读取TXT文件的值。
Matlab的textread函数就可以实现该功能，当然用fopen也可以，但textread 比较简单。其中textread的使用可以参考help文档。
读取txt文件后，就需要将16进制的数据转化成10进制了，此时需要用到Matlab的另外一个函数hex2dec(a)，这个函数的具体用法也见help。
下面就是读取二进制txt文件的例子：
假如txt文档中内容为0001 0010 0011 0100 1000，保存在pin.txt文档中
使用a=textread('pin.txt','%s')'
a ='0001' '0010' '0011' '0100' '1000'
可以看到数据保存为了char格式。
使用bin2dec
b=bin2dec(a)'
b =1 2 3 4 8
可以看到成功地转换成了十进制文件。
十六进制文件：0001 0010 0011 0100 1000 A B C AA
a=textread('pin.txt','%s')'
a ='0001' '0010' '0011' '0100' '1000' 'A' 'B' 'C' 'AA'
可以看到成功读取了文件。
b=hex2dec(a)'
b =1 16 17 256 4096 10 11 12 170
读取完毕。
将16进制数据转换成10进制后还需要将数据进行保存。此时要用到Matlab的dlmwrite函数，输入dlmwrite('c:b.txt',b,'delimiter', ' ')，就可以将得出数据存入新的b.txt文档。
上面就是如何通过Matlab读取十六进制数据并处理了。你学会了没有？

将从键盘输入的任何数据，通过dec2bin函数将其转换为二进制流之后通过2ASK调制，再通过sound函数将其通过pc自带的蜂鸣器发声，最后保存为声音文件。
clear all
clc
%***************由键盘输入文本文件并转化为bit信号********
% 获取字符进行二进制编码
user_input = input('>', 's');
x = user_input;
% 将字符转换为二进制编码
mess = dec2bin(user_input, 7);
bitstream = [];

% 将字符二进制编码转换为比特流
for i = 1:length(user_input())
    for k = 1:7
        bitstream = [bitstream double(mess(i, k))-48];
    end
end

a = bitstream;
len = length(bitstream);
t = 0:0.01:(len-0.01);
m = a(ceil(t+0.01));

%**************2ASK调制（OOK-通断键控）***************
%2ASK调制：S2ASK = m(t)*cos(wc*t)=m(t)*cos(2*pi*f*t)
%m(t)为数字信号，后者为载波
% 调制基带信号
fc=4000;          %载波频率
fs=80000;         %采样频率
ts=0:1/fs:(100*len-1)/fs;
tzxh1=cos(2*pi*fc*ts);%载波

ask = m.*tzxh1;%得到2ASK调制后的信号，采用点乘是因为矩阵元素相乘
startbit = [1 1 1 1];%开始标志码
stopbit = [0 1 1 0];%结束标志码
ask = [startbit ask stopbit];

for i = 1:1:100
    sound(ask, fs); %发声
end
audiowrite('ask.wav',ask,fs); %保存为声音文件

也可以读取声音文件：
[sig,fs]=audioread('ask.wav');

将声音文件读取出来，并通过相关解调，抽样判决，还原原本的文本输入：
sig1=sig';
sig2=sig1(1,5:(len*100+4));

% 2ASK相干解调
[b11,a11]=ellip(5,0.5,60,[2000,6000]*2/80000);%带通椭圆滤波器设计
e_ask1=filter(b11,a11,sig2);
e_ask2=2*e_ask1.*tzxh1*2;%相乘器

[b12,a12]=ellip(5,0.5,60,1000*2/80000);       %低通滤波器设计
ask_xgjt=filter(b12,a12,e_ask2);       %相干解调后，抽样判决前的结果

% 2ASK的相干解调法的抽样判决结果与原数据比较
for i=0:(len-1)
    if ask_xgjt(1,(i+1)*100)>0.5
        ask_hyjt(i*100+1:i*100+100)=ones(1,100);
    else
        ask_hyjt(i*100+1:i*100+100)=zeros(1,100); 
    end
end   

ask_hyjt=[ask_hyjt ask_hyjt(1,length(ask_hyjt)-1)];
output=[];
for q=0:1:((length(ask_hyjt)-1)/7/100-1)
    row = [];
    for k = 0:6
        % 将7位比特转化为一个字符串
        row = [row int2str(ask_hyjt(1,q*700+100*k+1))];
        % 转化为ascii编码
        incoming = bin2dec(row);
        % 转为字符
        received_message = char(incoming);
    end
 output=[output received_message];
end

PS：这是参考后理解的代码，并非全部原创编写代码。