このページは機械翻訳を使用して翻訳されました。最新版の英語を参照するには、ここをクリックします。
C/Aコードを使用したGPS受信機の取得とトラッキング
この例では、複数の衛星から従来の全地球測位システム (GPS ) 中間周波数 (IF) 波形を生成し、複合信号にノイズを追加し、初期同期を実行し、初期同期操作から検出された利用可能な衛星のコード位相と搬送周波数を追跡する方法を示します。このマルチ衛星 IF 波形は、各衛星のドップラー レートをシミュレートします。トラッキングは各衛星ごとに独立して行われます。この例に示す取得とトラッキングは、coarse取得コード (C/A コード) を含むGPS信号に基づいています。
はじめに
この例の IF 波形は、 GPSベースバンド波形に基づいて生成されます。GPSベースバンド波形を生成するためにさまざまなパラメーターを設定する方法の詳細については、GPS 波形生成 の例を参照してください。GPS受信機は、取得とトラッキングに加えて、ビット同期、フレーム同期、およびデータのデコードも実行します。この図には、 GPS受信機の動作に使用されるさまざまなブロックが示されています。緑色で強調表示されたブロックがこの例の範囲内です。
初期同期およびトラッキングモジュールは、高速（38.192 MHz）でデータを処理するため、高速モジュールと呼ばれます[1]。ビット同期、フレーム同期、およびデータ デコードを実行するモジュールは、データの速度が低い (1 秒あたり 50 ビット) ため、低速モジュールと呼ばれます。これらの低速モジュールから意味のあるデータを取得するには、重要なデータを高速モジュールで処理する必要があります。この例では、MATLAB ™ の取得およびトラッキングアルゴリズムに重点を置いているため、低速モジュールから意味のある結果を得るために処理されるデータはそれほど多くありません。
この例は3つの部分に分かれています。
シミュレーションパラメーターの設定 — 波形を生成し、チャネルを設定し、 GPS受信機を設定するためのさまざまなパラメーターを設定します。
GPS IF 波形の生成 — 複数の衛星からGPS IF 波形を生成し、遅延、ドップラー、ノイズを追加します。
初期同期とトラッキング — 可視衛星を見つけて、ドップラー周波数オフセットとタイミングオフセットのcoarse値を計算します。着信信号の変化するドップラー オフセットとコード位相を追跡します。
シミュレーションパラメーターの設定
サンプル出力を制御するすべての構成パラメーターを設定します。
一般的な構成パラメーターを設定します。この例では、ShowVisualizations フラグが可視化を制御します。例内のすべてのプロットを表示するには、これを true に設定します。WriteWaveformToFile フラグは、生成された波形をファイルに書き込むことを制御します。この例で生成された IF サンプルを使用して受信機をテストする場合は、WriteWaveformToFile を true に設定してサンプルをファイルにエクスポートできます。
ShowVisualizations =true; WriteWaveformToFile =
false;
送信機の構成
波形内のGPS衛星の数を指定します。
numSat = 4;
受信機で見える衛星の疑似ランダム ノイズ (PRN) 識別番号 (ID) を指定します。PRNIDs の長さは、指定された衛星の数以上である必要があります。長さが numSat プロパティより大きい場合、最初の numSat 値のみがシミュレーションに使用されます。
PRNIDs = [7; 11; 20; 28];
GPSレガシー ナビゲーション データ全体は 25 フレームで構成され、各フレームは 1500 ビットです。このナビゲーション データには、エフェメリス、時計、アルマナック情報全体が含まれます。ナビゲーション データ全体のGPS波形を生成するには長い時間がかかる可能性があるため、この例では、ナビゲーション データの 20 ビットのみの波形を生成する方法を示します。プロパティ NumNavDataBits を使用して、指定された数のデータ ビットの波形の生成を制御できます。
% Set this value to 1 to generate the waveform from the first bit of the % navigation data. Set this to any other number to generate waveform from % the point of interest. Starting from a random point in this example. NavDataBitStartIndex = 1321; % Set this value to control the number of navigation data bits in the % generated waveform NumNavDataBits = 20;
IF 中心周波数やサンプリング レートなどの受信信号のプロパティを定義します。この例では、中心周波数が 10 MHz の信号を示しています。CenterFrequency プロパティを 0 に設定することで、このプロパティを構成して複素ベースバンド信号を生成することができます。ベースバンド データで作業する場合、SampleRate プロパティを 10.23e6 に減らすと、例の動作がはるかに速くなります。
% Set the center frequency of the transmission waveform CenterFrequency = 10e6; % In Hz % Set the sample rate of the transmission waveform SampleRate = 38.192e6; % In samples/sec
チャネル構成
各衛星について、信号対雑音比 (SNR) を指定します。ベクトル内の各要素は、対応する衛星PRN ID の SNR 値を表します。IS- GPS-200 [2]規格では、C/Aコードの場合、最小受信信号電力は-158.5dBWと規定されている（[2]の表3-Vaを参照）。以下に示す熱雑音電力の式から、雑音電力は -130dBW と計算されます。したがって、受信機が動作する必要がある最小 SNR は -28.5 dB (信号電力 (dB) - ノイズ電力 (dB)) です。この例では、干渉による劣化は無視されます。
熱雑音電力は、ボルツマン定数 、ケルビン単位の受信機温度 、および信号帯域幅 の積によって決まります。C/A コードトラッキングを備えた従来のGPS信号の場合、考慮される一般的な信号帯域幅は 24 MHz です。
この例は、-20dB という低い SNR 値でも機能します。
SNRs = [-18; -18.5; -19.5; -19]; % In dB各衛星はGPS受信機から異なる距離にあるため、各信号に生じる遅延は衛星ごとに異なります。各衛星のC/A コードチップ遅延の数を列ベクトルとして初期化します。送信時間がC/A コードチップ持続時間 0.9775 マイクロ秒の整数倍ではない可能性が高いため、この遅延値は小数になる可能性があります。
考慮される遅延値は、ライブGPS信号から取得される実際の値よりもはるかに小さくなります。遅延値が小さいほどシミュレーション時間が短縮され、取得およびトラッキング機能を詳細に示す機会が得られます。
sigdelay = [300.34; 587.21; 425.89; 312.88]; % Number of C/A-code chips delay各衛星の速度と位置は異なるため、各衛星に発生するドップラーシフトとドップラー率は異なります。衛星シナリオにおけるレイテンシとドップラー効果の計算 の例は、衛星の位置と速度に基づいてレイテンシとドップラー レートを計算する方法を示しています。衛星シナリオにおけるレイテンシとドップラー効果の計算の例ではドップラーの変化がおおよそ正弦波のプロファイルに従うため、現在の例ではドップラーは正弦波的に変化する周波数オフセットとしてモデル化されます[3]。
各衛星のピーク ドップラー シフトとドップラー レートを初期化します。この例では、-10KHz から 10KHz までのドップラーシフトを追跡できます。
peakDoppler = [3589; 4256; 8596; 9568]; % In Hz dopplerRate = [1000; 500; 700; 500]; % In Hz/sec
受信機の構成
GPS受信機では、トラッキングアルゴリズムは、それぞれ周波数ロック ループ (FLL)、位相ロック ループ (PLL)、遅延ロック ループ (DLL) を使用して、周波数、位相、遅延を追跡します。ループ帯域幅が広いほど高速トラッキングが可能になりますが、SNR が低い場合はロックが失われる可能性があります。ループ帯域幅が狭いと、SNR が低い場合にパフォーマンスは良好になりますが、位相、周波数、遅延の変化の追跡は遅くなります。これらのプロパティの解釈については、初期同期とトラッキングに示されています。
PLLNoiseBandwidth = 90; % In Hz FLLNoiseBandwidth = 4; % In Hz DLLNoiseBandwidth = 1; % In Hz
統合時間を設定します。積分時間が長くなると、位相トラッキングアルゴリズムが向上します。データの遷移は 20 ミリ秒ごとに発生するため、統合時間として最大 20 ミリ秒を設定できます。ビット開始インデックスの位置は最初は不明なので、積分時間を 20 ミリ秒にすると、ビット遷移は積分信号の中間のどこかで発生します。これにより、トラッキングループが予期しない動作をするため、最初は積分時間を 1 ミリ秒に設定し、ビット同期が完了したら適切に増やします。この例ではビット同期は対象外であるため、積分時間は 1 ミリ秒に固定されます。積分時間を増やしながら、PLL ノイズ帯域幅を減らして、受信機が非常に低い SNR 値で動作するようにします。
IntegrationTime = 1e-3; % In secondsGPS IF波形を生成する
生成された IF 波形には、同相ブランチ上の高精度コード (パルス符号) と直交位相ブランチ上のC/A コードが含まれており、 GPS衛星から L1 周波数 (1575.42 MHz) で送信されます。送信されるデータの詳細については、IS- GPS-200規格[2]に記載されています。GPS 波形生成 では、 GPSデータをベースバンド送信波形に変換する方法について説明します。このベースバンド波形は IF 波形にアップコンバートされます。各衛星に指定された SNR に基づいて、受信信号電力に合わせて IF 波形がスケーリングされます。また、各衛星信号の遅延を個別にモデル化し、すべての信号を追加して複合 IF 波形を生成します。この複合 IF 波形に、一定量の加法性ホワイト ガウス ノイズ (AWGN) を追加します。
1 つの衛星から L1 周波数 (1575.42 MHz) で従来のGPS IF 波形を生成するには、この図に示す手順を実行します。
各衛星からのスケーリングされた IF 波形を追加し、一定の白色ガウス ノイズを複合信号に追加します。各信号電力は SNR に応じて調整されるため、各衛星信号の有効 SNR は異なります。この図は複数の衛星からの波形を組み合わせたものを示しています。下の図の各波形ジェネレーターは、前の図に示されている操作で構成されています。
波形を生成する前に、 ノイズ電力を計算して、各衛星信号を適切にスケーリングします。
k = 1.38064852e-23; % Boltzmann constant in Joules/Kelvin T = 300; % Room temperature in Kelvin B = 24e6; % Bandwidth in Hz N = k*T*B; % Thermal noise power in watts
遅延モデリングのために、追加のナビゲーション データ ビットを 1 つ考慮します。これにより、信号の最大 20 ミリ秒の遅延をモデル化できるようになります。信号内のより高い遅延値をモデル化するには、numBitsForDelay を整数値で増やします。
numBitsForDelay = 1;
各衛星からの従来のGPS送信データを生成するには、HelperGPSNavigationConfig オブジェクトを初期化して使用します。このオブジェクトには、衛星によって送信された情報が含まれています。この例では、すべてのパラメーターがすべての衛星に対して一定に保たれています。これは、この例の目的はデータをデコードすることではなく、任意のランダム データで十分な取得とトラッキングを示すことであるためです。この構成オブジェクトの詳細については、GPS 波形生成 の例を参照してください。
resultsig = 0; % Generate waveform from each satellite, scale and add based on SNR for isat = 1:numSat % Initialize the legacy navigation (LNAV) configuration object lnavConfig = HelperGPSNavigationConfig("SignalType","LNAV","PRNID",PRNIDs(isat)); % Generate the navigation data bits from the configuration object lnavData = HelperGPSNAVDataEncode(lnavConfig); % Configure the GPS waveform generation properties t = lnavConfig.HOWTOW*6; % First get the initial time % HOWTOW is an indication of next subframe starting point. Because each % subframe is 300 bits long, 300 bits must be subtracted from the % initial value to get the first subframe's starting value. This value % can be negative as well. Each bit is of 20 millisecond duration and to % get time elapsed for bits, bit index must be multiplied with 20e-3. bitDuration = 20e-3; % In sec pCodeRate = 10.23e6; % In Hz numPChipsPerNavBit = bitDuration*pCodeRate; navdatalen = length(lnavData); offsetTime = mod(NavDataBitStartIndex-301,navdatalen)*bitDuration; inittime = t + offsetTime; % For modeling delay, get one extra navigation bit from previous bit navBitIndices = mod(NavDataBitStartIndex+(-1*numBitsForDelay:(NumNavDataBits-1)),navdatalen); navBitIndices(navBitIndices==0) = navdatalen; navbits = lnavData(navBitIndices); navdata = 1-2*navbits; upSampledNavData = repelem(navdata,numPChipsPerNavBit,1); % Generate P-code and C/A-code pgen = gpsPCode(PRNID = PRNIDs(isat), ... InitialTime = inittime, ... OutputCodeLength = (NumNavDataBits+numBitsForDelay)*numPChipsPerNavBit); pcode = 1-2*double(pgen()); % Reduce the power of I-branch signal by 3 dB as per IS-GPS-200 [2]. % See table 3-Va in [2]. isig = pcode/sqrt(2); cacode = 1-2*double(gnssCACode(PRNIDs(isat),"GPS")); numCACodeBlocks = (NumNavDataBits+numBitsForDelay)*bitDuration*1e3; caCodeBlocks = repmat(cacode(:),numCACodeBlocks,1); % Because C/A-code is 10 times slower than P-code, repeat each sample % of C/A-code 10 times qsig = repelem(caCodeBlocks,10,1); % Generate the baseband waveform gpsBBWaveform = (isig + 1j*qsig).*upSampledNavData; % Initialize the number of samples per bit at IF numSamplesPerBit = SampleRate*bitDuration; % Introduce Doppler and up-convert the signal to IF upconvert = HelperGPSUpConvert; upconvert.PeakDoppler = peakDoppler(isat); upconvert.DopplerRate = dopplerRate(isat); upconvert.CenterFrequency = CenterFrequency; upconvert.IFSampleRate = SampleRate; gpsIFWaveform = upconvert(gpsBBWaveform); % Get the number of samples for delay caCodeRate = 1.023e6; numDelaySamples = floor(sigdelay(isat)*upconvert.IFSampleRate/caCodeRate); % Add delay to the signal by keeping samples of previous bit at the % beginning of the signal delayedSig = gpsIFWaveform(numSamplesPerBit-numDelaySamples+1:end); % Remove the final samples to make all signals of equal length delayedSig = delayedSig(1:end-numDelaySamples); % Scale this delayed signal to appropriate power level currentSNR = 10^(SNRs(isat)/10); % Convert to linear form signalpower = currentSNR*N; scaledsig = sqrt(signalpower)*delayedSig/rms(delayedSig); % Get the composite signal by adding the current satellite signal resultsig = resultsig + scaledsig; end % Add AWGN to the resultant IF waveform numSamples = length(resultsig); % For repeatable simulations, set the random number generator to default rng default; if CenterFrequency == 0 % Generate complex noise samples noisesig = (wgn(numSamples,1,10*log10(N)) + 1j*wgn(numSamples,1,10*log10(N)))./sqrt(2); else noisesig = wgn(numSamples,1,10*log10(N)); end rxwaveform = resultsig + noisesig; % Scale the received signal for having unit power rxwaveform = rxwaveform/rms(rxwaveform); if WriteWaveformToFile == 1 bbWriter = comm.BasebandFileWriter("IFWaveform.bb", ... upconvert.IFSampleRate,upconvert.CenterFrequency); bbWriter(rxwaveform); end
GPS 信号とノイズ信号のスペクトルを可視化します。信号はノイズ電力以下です。
if ShowVisualizations == 1 ifscope = spectrumAnalyzer(SampleRate = upconvert.IFSampleRate, ... PlotAsTwoSidedSpectrum = true, ... SpectrumType = "Power", ... SpectrumUnits = "dBW", ... Title = "IF Spectrum Comparison of GPS Signal with Thermal Noise", ... ShowLegend = true, ... ChannelNames = ["GPS IF waveform spectrum" "Noise spectrum"], ... YLimits = [-190 -155]); ifscope([resultsig, noisesig]); end
受信した信号のスペクトルを可視化します。スペクトルにはノイズだけがあるかのように表示されます。受信信号は単位電力にスケーリングされるため、スペクトルは単位電力ノイズ信号レベルになります。
if ShowVisualizations == 1 rxscope = spectrumAnalyzer(SampleRate = upconvert.IFSampleRate, ... PlotAsTwoSidedSpectrum = true, ... SpectrumType = "Power", ... SpectrumUnits = "dBW", ... Title = "Received signal IF spectrum after scaling"); rxscope(rxwaveform); end
初期同期とトラッキング
通常、受信機は次のいずれかのモードで起動します[4]。
ホット スタート - 受信機は現在の位置に関する情報を持ち、利用可能な衛星のアルマナック情報を持っています。この種の起動は、受信機が目に見える衛星を一時的に追跡できなくなったときに発生します。このモードでは、特定の衛星に対して初期同期がスキップされ、トラッキングのみが実行されます。
ウォーム スタート - 受信機の電源をオフにした後、受信機はあまり移動していないため、受信機は自身のおおよその位置を認識しています。また、受信機にはおおよそのアルマナック情報があります。この種の起動は、通常、受信機が数秒から数時間オフになっており、受信機がオフになった後もあまり動かなかった場合に発生します。このモードでは、利用可能なすべての衛星の検索がスキップされます。既知の可視衛星の初期同期が実行され、検出された衛星のトラッキングが行われます。
コールド スタート - 受信機は電源がオフになる前にはその位置に関する情報を持たず、世界中のどこかへ移動している可能性があります。受信機にはアルマナックに関する情報がありません。初期同期では全ての衛星を検索し、検出された衛星のトラッキングを行います。
コールド スタートでは、既知の情報がないため、自身の位置を特定するのに最も長い時間がかかります。受信機が自身の位置を特定するのにかかる時間は、初回位置特定時間 (TTFF) と呼ばれます。この例では、受信機が最初はコールド スタート モードにあり、可能性のある 32 個のGPS衛星すべての検索を開始すると想定しています。可視衛星が検出されると、それらの衛星のトラッキングが行われます。この例では、アルマナック情報をデコードせず、初期同期によって検出された衛星で動作します。
初期同期モジュールは、可視衛星をすべて検出し、各衛星からの受信信号のcoarse ドップラー シフトとcoarse時間オフセットを推定します[5]。この例では、FFT アルゴリズムを使用してC/A コードを周波数領域に変換し、周波数領域で信号との乗算を実行し、その後 IFFT を使用して時間領域に戻します。この同期は最も高速な初期同期アルゴリズムの1つです[5]。この並列コード位相探索アルゴリズムのブロック図を図に示します。この図は受信 IF 波形のブロック図を示していますが、同じブロックがベースバンド波形でも機能します (サンプリング レートと FFT サイズが異なります)。
この並列コードフェーズ アルゴリズムを使用するには、まず取得オブジェクトを初期化し、そのプロパティを構成します。
initialsync = gnssSignalAcquirer; initialsync.SampleRate = upconvert.IFSampleRate; initialsync.IntermediateFrequency = CenterFrequency
initialsync =
gnssSignalAcquirer with properties:
GNSSSignalType: "GPS C/A"
SampleRate: 38192000
IntermediateFrequency: 10000000
FrequencyRange: [-10000 10000]
FrequencyResolution: 500
DetectionThresholdFactor: 1.9000
トラッキングモジュールは、初期同期後に残るドップラーシフトとコード位相オフセットを補正します。追跡する必要があるパラメーターは、キャリア周波数、キャリア位相、およびコード遅延の 3 つです。これらの各パラメーターはフィードバック ループを使用して追跡されます。キャリア周波数は FLL を使用して追跡され、位相は PLL を使用して追跡され、コード位相オフセットは DLL を使用して追跡されます。どのトラッキングループ (FLL、PLL、DLL) の基本構造も同じであり、次の図に示されています。コンバイナは、入力信号とフィードバック ループで生成された参照信号を結合します。弁別器は、結合器によって生成された信号を使用して、入力信号と参照信号の誤差を検出します。トラッキングループは、このエラーをゼロに減らそうとします。このエラーは、FLL の周波数、PLL の位相、DLL の遅延など、任意の物理量になる可能性があります。このエラーは、ループ フィルターと呼ばれるローパス フィルターを通じて処理されます。ループ フィルターは、ディスクリミネータ出力のノイズを低減します。ループ フィルタは、基本的には、弁別器の出力を累積してエラー (弁別器の出力) をゼロにする積分器です。トラッキングループのもう 1 つのモジュールは数値制御発振器 (NCO) です。これも単純な積分器で、ループフィルタ[1]からの出力を累積します。
ループフィルタの主な役割は、弁別器[1]の出力からノイズを除去することです。このフィルタリングにより、エラーの突然の変化は十分な収束時間の経過後にのみ反映されます。この現象は、この例の最後に示したトラッキングモジュールの結果で確認できます。ここでは、識別子の出力はトラッキングループの出力よりもノイズが多いように見えます。トラッキングループで処理できるノイズの量は、ループ ノイズ帯域幅 (BW) によって制御されます。このパラメーターはトラッキングループごとに異なります。ループ帯域幅が大きいと、ループ フィルタによってトラッキングループを通過するノイズが増え、収束が速くなるためです。同様に、ノイズ帯域幅の狭いフィルターでは、フィルターを通過するノイズが少なくなり、収束が遅くなります。この表は、さまざまなトラッキングループのノイズ帯域幅に使用される標準的な値を示しています[1]。トラッキングループにとって重要な要素はループの順序です。高次のループは信号の突然の変化に耐えることができますが、安定性は低くなります。より安定した低次のループは、動的な環境では動作できません。この例では、1 次 FLL、2 次 PLL、および 1 次 DLL を使用します。次の表は、この例で使用される 3 つのループすべてのパラメーターを示しています。これらの値を構成するには、レシーバー構成を参照してください。
識別アルゴリズム | ループフィルタ帯域幅（Hz） | ループ順序 | |
PLL | 2象限アタン | 90 | 2 |
FLL | 4象限atan2 | 4 | 1 |
DLL | (E-L)/(2*(E+L)) | 1 | 1 |
この図は、 GPS受信機のトラッキングループのさまざまな相互接続を示しています。この図では、 はプロンプト コードからの統合信号、 はハーフ チップ アーリー コードからの統合信号、 はハーフ チップ レイト コードからの統合信号です。プロンプト コードは、 C/A コードの時間同期レプリカで乗算される信号です。今回は、ループ内でロックが確立されると同期コードが取得されます。早期コードはプロンプト コードを周期的に進めることによって生成されます。遅延コードは、プロンプト コードをチップ期間の半分だけ周期的に遅延させることによって生成されます。FLL弁別器は4象限atan2弁別器[1]を使用し、指定された積分時間内に位相器が回転する角度を評価することによって周波数オフセットを計算します。 によって形成されるベクトルは位相器です。周波数オフセットがある場合、この位相器は角速度で回転します。この角速度を使用して周波数オフセットを計算します。C/A コードは直交ブランチに配置されているため、位相オフセットがある場合、この位相器は y 軸に近い位置にはなりません。この位相と y 軸のなす角度から、位相弁別器の位相オフセットを計算します。遅延識別器は、早期コードまたは遅延コードと比較した場合、プロンプト コードの相関が常に最も高くなければならないという原理に基づいて動作します。プロンプトコードの相関が最高でない場合、遅延弁別器はプロンプトコードの相関が最大になるように遅延を修正する適切なエラーを生成します[1]。
キャリアとコードのトラッキングは、可視衛星上で行われる必要があります。各衛星ごとに個別のトラッキングループが必要です。このコードは、トラッキングに必要ないくつかのパラメーターを初期化します。初期同期を実行して、可視衛星を見つけて追跡します。
% Consider the data that is 1 millisecond long. numSamples = SampleRate*IntegrationTime; [allRxInput,prevSamples] = buffer(rxwaveform,numSamples); nFrames = size(allRxInput,2); numdetectsat = 0; PRNIDsToSearch = 1:32; for iBuffer = 1:nFrames bufferWave = allRxInput(:,iBuffer); if iBuffer == 1 % This example assumes a hot start for all the satellites once % initial synchronization is complete. Hence, perform initial % synchronization only once in this example. When decoding the % almanac data too, based on the available satellites, initial % synchronization can be performed for the visible satellites only. numSamplesForInitSync = SampleRate*2e-3; % Corresponding to 2 milliseccond [y,corrval] = initialsync(allRxInput(1:numSamplesForInitSync).',1:32); % Search for all 32 GPS satellites corrval = corrval(1:numSamplesForInitSync/2,:,:); PRNIDsToSearch = y(y(:,4).IsDetected==1,1).PRNID.'; % Get row vector doppleroffsets = y(y(:,4).IsDetected==1,2).FrequencyOffset; codephoffsets = y(y(:,4).IsDetected==1,3).CodePhaseOffset; numdetectsat = length(PRNIDsToSearch); disp("The detected satellite PRN IDs: " + num2str(PRNIDsToSearch)) if ShowVisualizations == 1 figure; mesh(initialsync.FrequencyRange(1):initialsync.FrequencyResolution:initialsync.FrequencyRange(2), ... 0:size(corrval,1)-1, corrval(:,:,1)); xlabel("Doppler Offset") ylabel("Code Phase Offset") zlabel("Correlation") msg = "Correlation Plot for PRN ID: " + PRNIDsToSearch(1); title(msg) end % Initialize all the properties which must be accumulated. accuph = zeros(nFrames,numdetectsat); % Each column represents data from a satellite accufqy = zeros(nFrames,numdetectsat); accufqyerr = zeros(nFrames,numdetectsat); accupherr = zeros(nFrames,numdetectsat); accuintegwave = zeros(nFrames,numdetectsat); accudelay = zeros(nFrames,numdetectsat); accudelayerr = zeros(nFrames,numdetectsat); % Update properties for each tracking loop carrierCodeTrack = gnssSignalTracker; carrierCodeTrack.SampleRate = SampleRate; carrierCodeTrack.IntermediateFrequency = CenterFrequency; carrierCodeTrack.PLLNoiseBandwidth = PLLNoiseBandwidth; carrierCodeTrack.FLLNoiseBandwidth = FLLNoiseBandwidth; carrierCodeTrack.DLLNoiseBandwidth = DLLNoiseBandwidth; carrierCodeTrack.IntegrationTime = IntegrationTime; carrierCodeTrack.PRNID = PRNIDsToSearch; carrierCodeTrack.InitialFrequencyOffset = doppleroffsets; carrierCodeTrack.InitialCodePhaseOffset = codephoffsets; end [integwave,trackinfo] = carrierCodeTrack(bufferWave); % Accumulate the values to see the results at the end accuintegwave(iBuffer,:) = integwave; accufqyerr(iBuffer,:) = trackinfo.FrequencyError; accufqy(iBuffer,:) = trackinfo.FrequencyEstimate; accupherr(iBuffer,:) = trackinfo.PhaseError; accuph(iBuffer,:) = trackinfo.PhaseEstimate; accudelayerr(iBuffer,:) = trackinfo.DelayError; accudelay(iBuffer,:) = trackinfo.DelayEstimate; end
The detected satellite PRN IDs: 7 11 28 20
trackedSignal = accuintegwave; % Useful for further GPS receiver processing次のグラフはトラッキングループの結果を示しています。推定された位相誤差は、周波数ロック ループが収束するまで (最初のプロット) 収束しないことを示しています (2 番目のプロット)。周波数の誤差によって位相の誤差が発生するため、この動作は予想されます。また、ドップラーシフトが時間とともに変化するため、PLL と FLL の出力も時間とともに変化するように見えます。
if ShowVisualizations == 1 % for isat = 1:numdetectsat for isat = 1 % See tracking results of all the detected satellites by using above line groupTitle = "Tracking Loop Results for Satellite PRN ID:" + PRNIDsToSearch(isat); figure % Plot the frequency discriminator output subplot(2,1,1) plot(accufqyerr(:,isat)) xlabel("Milliseconds") ylabel("Frequency Error") title("Frequency Discriminator Output") % Plot the FLL output subplot(2,1,2) plot(accufqy(:,isat)) xlabel("Milliseconds") ylabel("Estimated Frequency Offset") title("FLL Output") sgtitle("FLL " + groupTitle) figure % Plot the phase discriminator output subplot(2,1,1) plot(accupherr(:,isat)) xlabel("Milliseconds") ylabel("Phase Error") title("Phase Discriminator Output") % Plot the PLL output subplot(2,1,2) plot(accuph(:,isat)) xlabel("Milliseconds") ylabel("Estimated Phase") title("PLL Output") sgtitle("PLL " + groupTitle) figure % Plot the delay discriminator output subplot(2,1,1) plot(accudelayerr(:,isat)) xlabel("Milliseconds") ylabel("Delay Error") title("Delay Discriminator Output") % Plot the DLL output subplot(2,1,2) plot(accudelay(:,isat)) xlabel("Milliseconds") ylabel("Estimated Delay") title("DLL output in the units of number of C/A-code chips") sgtitle("DLL " + groupTitle) end end
トラッキング後、信号のコンスタレーションをプロットします。
if ShowVisualizations == 1 demodsamples = accuintegwave(301:end,:)/rms(accuintegwave(:)); if ~isempty(demodsamples) scatterplot(demodsamples(:)) end end
その他の調査
この例では、4 つのGPS衛星のシミュレーションを示します。GPS対応のスマートフォンの場合、通常、任意の時点で 6 つのGPS衛星が見えます。可視GPS衛星の数を増やして結果を観察できます。ドップラー、SNR、遅延を適切に初期化します。
この例では、各衛星信号の遅延とドップラーにハードコードされた値を使用します。Satellite Communications Toolbox の シナリオの生成と可視化 を使用して、これら 2 つのパラメーターを実際のGPS衛星としてモデル化することもできます。これらの機能を使用して、送信機と受信機間の距離を計算し、光速またはその他の遅延モデルを使用して遅延をモデル化します。詳細については、End-to-End GPS Legacy Navigation Receiver Using C/A-Code の例を参照してください。
この例の レシーバー構成プロパティ は、この例で設定された条件に合わせて最適化されるように設定されています。例を動作させたい条件に合わせて、これらの受信機構成パラメーターを更新できます。たとえば、この例で使用されているものよりもはるかに低い SNR の場合、ループのノイズ帯域幅を、トラッキングアルゴリズムが適切に機能する適切な値まで減らします。ノイズ帯域幅の値を下げると、ノイズが多く移動度が低い状況でも信号を適切に追跡できるようになります。
サポート ファイル
この例では、次のヘルパー ファイルを使用します。
HelperGPSNAVDataEncode.m — 構成オブジェクト内のデータからナビゲーションデータをビットにエンコードします
HelperGPSNavigationConfig.m — GPSナビゲーションデータの構成オブジェクトを作成します
HelperGPSUpConvert.m — ドップラーシフトをモデル化し、ベースバンド信号を中間周波数にアップコンバートします。
参考文献
[1] Elliott D. KaplanとC. Hegarty編、「Understanding GPS/ GNSS」原理とアプリケーション 、第3版、 GNSSテクノロジーとアプリケーションシリーズ（ボストン、ロンドン：（アーテックハウス、2017年）。
[2] IS- GPS-200、改訂:L. NAVSTAR GPS宇宙セグメント/ナビゲーション ユーザー セグメント インターフェース。2020年5月14日; コード識別:66RP1。
[3] CCSDS -G-1.グリーンブック。問題 1。「SCCC - 定義とパフォーマンスの概要」宇宙データシステム標準に関する情報レポートワシントンD.C.:CCSDS、2019年4月。
[4] ジョン・W・ベッツ、エンジニアリング衛星ベースのナビゲーションとタイミング：全地球航法衛星システム、信号、および受信機 (John Wiley & Sons, Ltd、2015)。
[5] K. Borre編、「A ソフトウェア定義GPSおよびガリレオ受信機」単一周波数アプローチ、応用および数値調和解析（ボストン、マサチューセッツ州：(Birkhäuser、2007)。
