Skip to content
GitLab
Commit 516e2d17 authored by Marion Bandet's avatar Marion Bandet
Browse files

update WERA calibrations

parent 5fde9a30
Hide whitespace changes
Inline Side-by-side
build/html/WERAcalibration.html
View file @ 516e2d17
<!DOCTYPE html PUBLIC "-//W3C//DTD XHTML 1.0 Transitional//EN"
"http://www.w3.org/TR/xhtml1/DTD/xhtml1-transitional.dtd">
<html xmlns="http://www.w3.org/1999/xhtml">
<head>
<meta http-equiv="Content-Type" content="text/html; charset=utf-8" />
<title>WERA calibrations &mdash; documentation_HFR_analysis 0.0.1 documentation</title>
<link rel="stylesheet" href="_static/default.css" type="text/css" />
<link rel="stylesheet" href="_static/pygments.css" type="text/css" />
<script type="text/javascript">
var DOCUMENTATION_OPTIONS = {
URL_ROOT: './',
VERSION: '0.0.1',
COLLAPSE_INDEX: false,
FILE_SUFFIX: '.html',
HAS_SOURCE: true
};
</script>
<script type="text/javascript" src="_static/jquery.js"></script>
<script type="text/javascript" src="_static/underscore.js"></script>
<script type="text/javascript" src="_static/doctools.js"></script>
<link rel="top" title="documentation_HFR_analysis 0.0.1 documentation" href="index.html" />
<link rel="next" title="Comparison WERA and ADCP radial velocities September 2013" href="Comparison_WERA-ADCP_Sept2013.html" />
<link rel="prev" title="WERA radial currents processing" href="WERAprocessing.html" />
</head>
<body>
<div class="related">
<h3>Navigation</h3>
<ul>
<li class="right" style="margin-right: 10px">
<a href="genindex.html" title="General Index"
accesskey="I">index</a></li>
<li class="right" >
<a href="Comparison_WERA-ADCP_Sept2013.html" title="Comparison WERA and ADCP radial velocities September 2013"
accesskey="N">next</a> |</li>
<li class="right" >
<a href="WERAprocessing.html" title="WERA radial currents processing"
accesskey="P">previous</a> |</li>
<li><a href="index.html">documentation_HFR_analysis 0.0.1 documentation</a> &raquo;</li>
</ul>
</div>
<div class="document">
<div class="documentwrapper">
<div class="bodywrapper">
<div class="body">
<!DOCTYPE html
PUBLIC "-//W3C//DTD HTML 4.01 Transitional//EN">
<html><head>
<meta http-equiv="Content-Type" content="text/html; charset=utf-8">
<!--
This HTML was auto-generated from MATLAB code.
To make changes, update the MATLAB code and republish this document.
--><title>WERAcalibration</title><meta name="generator" content="MATLAB 7.14"><link rel="schema.DC" href="http://purl.org/dc/elements/1.1/"><meta name="DC.date" content="2015-10-06"><meta name="DC.source" content="WERAcalibration.m"><style type="text/css">
html,body,div,span,applet,object,iframe,h1,h2,h3,h4,h5,h6,p,blockquote,pre,a,abbr,acronym,address,big,cite,code,del,dfn,em,font,img,ins,kbd,q,s,samp,small,strike,strong,sub,sup,tt,var,b,u,i,center,dl,dt,dd,ol,ul,li,fieldset,form,label,legend,table,caption,tbody,tfoot,thead,tr,th,td{margin:0;padding:0;border:0;outline:0;font-size:100%;vertical-align:baseline;background:transparent}body{line-height:1}ol,ul{list-style:none}blockquote,q{quotes:none}blockquote:before,blockquote:after,q:before,q:after{content:'';content:none}:focus{outine:0}ins{text-decoration:none}del{text-decoration:line-through}table{border-collapse:collapse;border-spacing:0}
html { min-height:100%; margin-bottom:1px; }
html body { height:100%; margin:0px; font-family:Arial, Helvetica, sans-serif; font-size:10px; color:#000; line-height:140%; background:#fff none; overflow-y:scroll; }
html body td { vertical-align:top; text-align:left; }
h1 { padding:0px; margin:0px 0px 25px; font-family:Arial, Helvetica, sans-serif; font-size:1.5em; color:#d55000; line-height:100%; font-weight:normal; }
h2 { padding:0px; margin:0px 0px 8px; font-family:Arial, Helvetica, sans-serif; font-size:1.2em; color:#000; font-weight:bold; line-height:140%; border-bottom:1px solid #d6d4d4; display:block; }
h3 { padding:0px; margin:0px 0px 5px; font-family:Arial, Helvetica, sans-serif; font-size:1.1em; color:#000; font-weight:bold; line-height:140%; }
a { color:#005fce; text-decoration:none; }
a:hover { color:#005fce; text-decoration:underline; }
a:visited { color:#004aa0; text-decoration:none; }
p { padding:0px; margin:0px 0px 20px; }
img { padding:0px; margin:0px 0px 20px; border:none; }
p img, pre img, tt img, li img { margin-bottom:0px; }
ul { padding:0px; margin:0px 0px 20px 23px; list-style:square; }
ul li { padding:0px; margin:0px 0px 7px 0px; }
ul li ul { padding:5px 0px 0px; margin:0px 0px 7px 23px; }
ul li ol li { list-style:decimal; }
ol { padding:0px; margin:0px 0px 20px 0px; list-style:decimal; }
ol li { padding:0px; margin:0px 0px 7px 23px; list-style-type:decimal; }
ol li ol { padding:5px 0px 0px; margin:0px 0px 7px 0px; }
ol li ol li { list-style-type:lower-alpha; }
ol li ul { padding-top:7px; }
ol li ul li { list-style:square; }
.content { font-size:1.2em; line-height:140%; padding: 20px; }
pre, tt, code { font-size:12px; }
pre { margin:0px 0px 20px; }
pre.error { color:red; }
pre.codeinput { padding:10px; border:1px solid #d3d3d3; background:#f7f7f7; }
pre.codeoutput { padding:10px 11px; margin:0px 0px 20px; color:#4c4c4c; }
@media print { pre.codeinput, pre.codeoutput { word-wrap:break-word; width:100%; } }
span.keyword { color:#0000FF }
span.comment { color:#228B22 }
span.string { color:#A020F0 }
span.untermstring { color:#B20000 }
span.syscmd { color:#B28C00 }
.footer { width:auto; padding:10px 0px; margin:25px 0px 0px; border-top:1px dotted #878787; font-size:0.8em; line-height:140%; font-style:italic; color:#878787; text-align:left; float:none; }
.footer p { margin:0px; }
</style></head><body><div class="content"><h2>Contents</h2><div><ul><li><a href="#3">phases and time of calibration</a></li><li><a href="#4">separate between calibrations made with Helzel's system and Pierre's system</a></li><li><a href="#5">theoretical values</a></li></ul></div><pre class="codeinput">clear <span class="string">all</span>
close <span class="string">all</span>
<span class="comment">%cd('/home/bandma01/Documents/WERAs/')</span>
cd(<span class="string">'/home/bandma01/hfr-data-processing/calibrations/'</span>)
<span class="comment">% site = 'PAO';</span>
site = <span class="string">'PAB'</span>;
cd(site)
<span class="comment">% system</span>
<span class="keyword">for</span> sys = 1:2
</pre><pre class="codeinput"> <span class="keyword">if</span> sys == 1
systm = <span class="string">'Helzel'</span>;
<span class="keyword">elseif</span> sys == 2
systm = <span class="string">'Flament'</span>;
<span class="keyword">end</span>
clear <span class="string">J</span> <span class="string">I</span> <span class="string">E</span> <span class="string">Phase</span> <span class="string">calib_file</span> <span class="string">calib_hr</span> <span class="string">calib_jd</span> <span class="string">calib_min</span> <span class="string">calib_time</span> <span class="string">calib_value</span>
clear <span class="string">calib_year</span> <span class="string">d</span> <span class="string">data</span> <span class="string">delta_phase_obs</span> <span class="string">delta_phase_th</span> <span class="string">gps_data</span> <span class="string">index_timecalib</span>
clear <span class="string">theta_pairs</span>
<span class="keyword">if</span> strcmp(site,<span class="string">'PAO'</span>) == 1
load(<span class="string">'HelzelPAO_Track_29-OCT-14 183138.mat'</span>)
<span class="keyword">elseif</span> strcmp(site,<span class="string">'PAB'</span>) == 1
<span class="keyword">if</span> strcmp(systm,<span class="string">'Helzel'</span>) == 1
load(<span class="string">'HelzelPAB_Track_30-OCT-14 175446.mat'</span>)
<span class="keyword">elseif</span> strcmp(systm,<span class="string">'Flament'</span>) == 1
load(<span class="string">'FlamentPAB_Track_30-OCT-14 182008.mat'</span>)
<span class="keyword">end</span>
<span class="keyword">end</span>
</pre><h2>phases and time of calibration<a name="3"></a></h2><pre class="codeinput"> <span class="keyword">if</span> strcmp(site,<span class="string">'PAO'</span>) == 1
cd(<span class="string">'2014302'</span>)
<span class="keyword">elseif</span> strcmp(site,<span class="string">'PAB'</span>) == 1
cd(<span class="string">'2014303'</span>)
<span class="keyword">end</span>
<span class="comment">% list files</span>
clear <span class="string">d</span>
d = dir(<span class="string">'*.mat'</span>); <span class="comment">% returns a structure with name, date, bytes, isdir, datenum</span>
Nfiles = length(d); <span class="comment">% nombre de fichiers a analyser</span>
<span class="comment">% look at the names of files and extract information</span>
<span class="keyword">for</span> c = 1:Nfiles
calib_year(c) = str2double(d(c).name(1:4));
calib_jd(c) = str2double(d(c).name(5:7));
calib_hr(c) = str2double(d(c).name(8:9));
calib_min(c) = str2double(d(c).name(10:11));
load(d(c).name,<span class="string">'phase'</span>)
Phase(:,c) = phase;
<span class="keyword">end</span>
clear <span class="string">phase</span>
delta_phase_obs = modsym(diff(Phase,1),180);
To = datenum(2014,01,01,0,0,0);
<span class="comment">% add 20 sec to center calibration time</span>
calib_time = To + calib_jd -1 + calib_hr./24 + calib_min/24/60 + 20/24/3600;
<span class="comment">% datevec(calib_time)</span>
<span class="keyword">if</span> strcmp(site,<span class="string">'PAO'</span>) == 1
<span class="comment">% calib time 2min34sec late wrt GPS</span>
calib_time = calib_time + 2/24/60 + 34/24/3600;
<span class="keyword">end</span>
</pre><h2>separate between calibrations made with Helzel's system and Pierre's system<a name="4"></a></h2><p>date of new calibration files Helzel 2014-10-30 16:52 to 17:54 Flament 2014-10-30 17:55 to 18:20</p><pre class="codeinput"> T1 = datenum(2014,10,30,17,59,00);
clear <span class="string">J</span>
<span class="keyword">if</span> strcmp(systm,<span class="string">'Helzel'</span>) == 1
J = find(calib_time &gt; T1);
<span class="keyword">elseif</span> strcmp(systm,<span class="string">'Flament'</span>) == 1
J = find(calib_time &gt; T1);
<span class="keyword">end</span>
<span class="comment">%remove unwanted data to avoid confusion!</span>
calib_year(J)=[];
calib_jd(J)=[];
calib_hr(J)=[];
calib_min(J)=[];
Phase(:,J)=[];
delta_phase_obs(:,J)=[];
calib_time(J)=[];
d(J)=[];
clear <span class="string">c</span> <span class="string">Nfiles</span>
<span class="comment">%+++++++++++++++++++++++++++++++++++++++++</span>
<span class="comment">% for each calibration file (series# 1)</span>
<span class="comment">%+++++++++++++++++++++++++++++++++++++++++</span>
<span class="keyword">for</span> calib_file = 1:length(calib_time)
time_calib = calib_time(calib_file);
<span class="comment">% find corresponding file in gps data</span>
clear <span class="string">I</span>
I = find(data(:,1) &gt; time_calib);
<span class="comment">% take closest point for now</span>
I = I(1);
<span class="comment">% datevec(data(I,1))</span>
<span class="comment">% datevec(time_calib)</span>
index_timecalib(calib_file) = I;
<span class="keyword">end</span>
<span class="comment">% check</span>
<span class="comment">% datevec(data(index_timecalib,1))</span>
<span class="comment">% isolate gps calibration data</span>
gps_data = data(index_timecalib,:); <span class="comment">% time lat lon</span>
figure(1)
hold <span class="string">on</span>
h1=plot(gps_data(:,2),gps_data(:,3),<span class="string">'.'</span>);
grid <span class="string">on</span>
<span class="keyword">if</span> strcmp(systm,<span class="string">'Helzel'</span>) == 1
set(h1,<span class="string">'marker'</span>,<span class="string">'.'</span>,<span class="string">'color'</span>,<span class="string">'b'</span>)
<span class="keyword">elseif</span> strcmp(systm,<span class="string">'Flament'</span>) == 1
set(h1,<span class="string">'marker'</span>,<span class="string">'x'</span>,<span class="string">'color'</span>,<span class="string">'k'</span>)
<span class="keyword">end</span>
<span class="comment">% load antenna coordinates</span>
cd <span class="string">..</span>
load([site <span class="string">'_Rx_coords.mat'</span>])
hold <span class="string">on</span>
plot(lon,lat,<span class="string">'x'</span>)
grid <span class="string">on</span>
title([<span class="string">'site = '</span> site ])
</pre><img vspace="5" hspace="5" src="WERAcalibration_01.png" alt=""> <img vspace="5" hspace="5" src="WERAcalibration_16.png" alt=""> <h2>theoretical values<a name="5"></a></h2><p>wavelength</p><pre class="codeinput"> lambda = 3E8/16.15E6;
<span class="comment">% initialize variables</span>
theta_pairs = nan(11,length(gps_data));
<span class="comment">% for each gps point</span>
<span class="keyword">for</span> pt = 1:length(gps_data)
<span class="comment">% compute distance between GPS and antennas</span>
<span class="keyword">for</span> ant = 1:12
clear <span class="string">E</span> <span class="string">W</span>
[E,W] = lonlat2km(gps_data(pt,2),gps_data(pt,3),lon(ant),lat(ant));
ddist = sqrt(E.^2+W.^2);
r(ant) = ddist*1000; <span class="comment">%meters</span>
<span class="keyword">end</span>
<span class="comment">%for each pair of antennas</span>
<span class="keyword">for</span> p = 1:11
<span class="comment">%delta_r</span>
delta_r(p) = r(p+1) - r(p);
<span class="comment">% theoretical phase</span>
delta_phase_th(p,pt) = modsym((delta_r(p)*360./lambda),180);
<span class="comment">% position of middle of pairs</span>
clear <span class="string">lon_mid</span> <span class="string">lat_mid</span>
lon_mid = (lon(p+1)+lon(p))/2;
lat_mid = (lat(p+1)+lat(p))/2;
<span class="comment">% distance between middle of pair and one of the antenna in the pair</span>
<span class="comment">% clear E W</span>
<span class="comment">% [E,W] = lonlat2km(lon_mid,lat_mid,lon(p),lat(p));</span>
<span class="comment">% ddist = sqrt(E.^2+W.^2);</span>
<span class="comment">% a(p) = ddist*1000; %meters</span>
<span class="comment">% distance between middle of pair and GPS pt</span>
clear <span class="string">ddsit</span> <span class="string">phasAng</span>
[ddist,phasAng] = sw_dist([gps_data(pt,3) lat_mid],[gps_data(pt,2) lon_mid],<span class="string">'km'</span>);
r_mid(p) = ddist*1000; <span class="comment">%meters</span>
zeta(p) = phasAng;
<span class="comment">% distance between pair of antennas</span>
clear <span class="string">ddsit</span> <span class="string">phasAng</span>
[ddist,phasAng] = sw_dist([lat(p) lat(p+1)],[lon(p) lon(p+1)],<span class="string">'km'</span>);
phi(p) = phasAng;
<span class="comment">% azimuth</span>
theta_pairs(p,pt) = modsym(zeta(p)-phi(p)+90,180);
<span class="comment">%rad2deg(theta_pairs)</span>
<span class="keyword">end</span>
<span class="keyword">end</span>
figure(1)
legend(<span class="string">'gps'</span>,<span class="string">'antennas'</span>,<span class="string">'location'</span>,<span class="string">'southeast'</span>)
xlabel(<span class="string">'longitude'</span>)
ylabel(<span class="string">'latitude'</span>)
<span class="comment">% for each pair of antenna</span>
<span class="keyword">for</span> p = 1:11
figure(p+1)
h2=plot(theta_pairs(p,:),delta_phase_th(p,:),<span class="string">'x'</span>);
hold <span class="string">on</span>
h3=plot(theta_pairs(p,:),-delta_phase_obs(p,:),<span class="string">'or'</span>);
grid <span class="string">on</span>
xlabel(<span class="string">'incidence (degrees)'</span>)
ylabel(<span class="string">'phase difference (degrees)'</span>)
legend(<span class="string">'theoretical'</span>,<span class="string">'observed'</span>,<span class="string">'location'</span>,<span class="string">'best'</span>)
title([site <span class="string">' - antennas '</span> num2str(p) <span class="string">'-'</span> num2str(p+1)])
<span class="keyword">if</span> strcmp(systm,<span class="string">'Helzel'</span>) == 1
set(h2,<span class="string">'color'</span>,<span class="string">'b'</span>)
set(h3,<span class="string">'color'</span>,<span class="string">'r'</span>,<span class="string">'marker'</span>,<span class="string">'o'</span>)
<span class="keyword">elseif</span> strcmp(systm,<span class="string">'Flament'</span>) == 1
set(h2,<span class="string">'color'</span>,<span class="string">'k'</span>,<span class="string">'marker'</span>,<span class="string">'*'</span>)
set(h3,<span class="string">'color'</span>,<span class="string">'r'</span>,<span class="string">'marker'</span>,<span class="string">'sq'</span>,<span class="string">'markerfacecolor'</span>,<span class="string">'r'</span>)
<span class="keyword">end</span>
<span class="comment">% compute calibration value for antennas pairs</span>
calib_value(p) = nanmean(modsym(-delta_phase_obs(p,:)-delta_phase_th(p,:),180),2);
<span class="comment">% compute calibration value for each antenna</span>
phi_ant(p+1,sys) = sum(calib_value);
<span class="keyword">end</span>
P=ones(16,3);
P(13:16,:)=0;
P(1:12,3)=phi_ant(:,1);
<span class="comment">% save file (we only keep results from Helzel</span>
<span class="keyword">if</span> strcmp(site,<span class="string">'PAO'</span>) == 1
dlmwrite(<span class="string">'calibration_pao_29OCT2014.wera'</span>,P,<span class="string">'delimiter'</span>,<span class="string">'\t'</span>,<span class="string">'precision'</span>,<span class="string">'%.7f'</span>)
<span class="keyword">elseif</span> strcmp(site,<span class="string">'PAB'</span>) == 1
dlmwrite(<span class="string">'calibration_pab_30OCT2014.wera'</span>,P,<span class="string">'delimiter'</span>,<span class="string">'\t'</span>,<span class="string">'precision'</span>,<span class="string">'%.7f'</span>)
<span class="keyword">end</span>
figure(13)
hold <span class="string">on</span>
h4=plot(calib_value,<span class="string">'.-'</span>);
grid <span class="string">on</span>
xlabel(<span class="string">'pairs of antennas'</span>)
ylabel({<span class="string">'phase difference between'</span>,<span class="string">'theoretical and observed values (degrees)'</span>})
title(site)
<span class="keyword">if</span> strcmp(systm,<span class="string">'Helzel'</span>) == 1
set(h4,<span class="string">'color'</span>,<span class="string">'b'</span>)
<span class="keyword">elseif</span> strcmp(systm,<span class="string">'Flament'</span>) == 1
set(h4,<span class="string">'color'</span>,<span class="string">'k'</span>)
legend(<span class="string">'Helzel'</span>,<span class="string">'Flament'</span>,<span class="string">'location'</span>,<span class="string">'southeast'</span>)
<span class="keyword">end</span>
figure(14)
hold <span class="string">on</span>
h5=plot([1:12],phi_ant(:,sys),<span class="string">'r'</span>);
grid <span class="string">on</span>
<span class="keyword">if</span> strcmp(systm,<span class="string">'Helzel'</span>) == 1
set(h5,<span class="string">'marker'</span>,<span class="string">'o'</span>)
<span class="keyword">elseif</span> strcmp(systm,<span class="string">'Flament'</span>) == 1
set(h5,<span class="string">'marker'</span>,<span class="string">'sq'</span>,<span class="string">'markerfacecolor'</span>,<span class="string">'r'</span>)
legend(<span class="string">'Helzel'</span>,<span class="string">'Flament'</span>,<span class="string">'location'</span>,<span class="string">'southeast'</span>)
<span class="keyword">end</span>
xlabel(<span class="string">'antenna number'</span>)
ylabel(<span class="string">'calibration phase to be applied (degrees)'</span>)
title(site)
<span class="keyword">if</span> strcmp(site,<span class="string">'PAO'</span>) == 1
<span class="keyword">break</span>
<span class="keyword">end</span>
</pre><img vspace="5" hspace="5" src="WERAcalibration_02.png" alt=""> <img vspace="5" hspace="5" src="WERAcalibration_03.png" alt=""> <img vspace="5" hspace="5" src="WERAcalibration_04.png" alt=""> <img vspace="5" hspace="5" src="WERAcalibration_05.png" alt=""> <img vspace="5" hspace="5" src="WERAcalibration_06.png" alt=""> <img vspace="5" hspace="5" src="WERAcalibration_07.png" alt=""> <img vspace="5" hspace="5" src="WERAcalibration_08.png" alt=""> <img vspace="5" hspace="5" src="WERAcalibration_09.png" alt=""> <img vspace="5" hspace="5" src="WERAcalibration_10.png" alt=""> <img vspace="5" hspace="5" src="WERAcalibration_11.png" alt=""> <img vspace="5" hspace="5" src="WERAcalibration_12.png" alt=""> <img vspace="5" hspace="5" src="WERAcalibration_13.png" alt=""> <img vspace="5" hspace="5" src="WERAcalibration_14.png" alt=""> <img vspace="5" hspace="5" src="WERAcalibration_15.png" alt=""> <img vspace="5" hspace="5" src="WERAcalibration_17.png" alt=""> <img vspace="5" hspace="5" src="WERAcalibration_18.png" alt=""> <img vspace="5" hspace="5" src="WERAcalibration_19.png" alt=""> <img vspace="5" hspace="5" src="WERAcalibration_20.png" alt=""> <img vspace="5" hspace="5" src="WERAcalibration_21.png" alt=""> <img vspace="5" hspace="5" src="WERAcalibration_22.png" alt=""> <img vspace="5" hspace="5" src="WERAcalibration_23.png" alt=""> <img vspace="5" hspace="5" src="WERAcalibration_24.png" alt=""> <img vspace="5" hspace="5" src="WERAcalibration_25.png" alt=""> <img vspace="5" hspace="5" src="WERAcalibration_26.png" alt=""> <img vspace="5" hspace="5" src="WERAcalibration_27.png" alt=""> <img vspace="5" hspace="5" src="WERAcalibration_28.png" alt=""> <img vspace="5" hspace="5" src="WERAcalibration_29.png" alt=""> <img vspace="5" hspace="5" src="WERAcalibration_30.png" alt=""> <pre class="codeinput"><span class="keyword">end</span>
<span class="comment">%add current calibration phases</span>
<span class="keyword">if</span> strcmp(site,<span class="string">'PAO'</span>) == 1
load <span class="string">calibration_pao.wera</span>
figure(14)
plot([1:12],calibration_pao(1:12,3),<span class="string">'b.-'</span>)
legend(<span class="string">'Helzel'</span>,<span class="string">'original calibration'</span>,<span class="string">'location'</span>,<span class="string">'northwest'</span>)
<span class="keyword">elseif</span> strcmp(site,<span class="string">'PAB'</span>) == 1
load <span class="string">calibration_pab.wera</span>
figure(14)
plot([1:12],calibration_pab(1:12,3),<span class="string">'b.-'</span>)
legend(<span class="string">'Helzel'</span>,<span class="string">'Flament'</span>,<span class="string">'original calibration'</span>,<span class="string">'location'</span>,<span class="string">'southeast'</span>)
<span class="keyword">end</span>
<span class="comment">% save figures</span>
save_flag = 1;
figdir = (<span class="string">'/home/bandma01/Dropbox/UQAR/MEOPAR_Marion/figures/WERAcalibration'</span>);
cd(figdir)
fig_hdl = figure(1);
<span class="keyword">if</span> save_flag
fname = [<span class="string">'configuration_'</span> site <span class="string">'_calibrations'</span>];
fig_w = 8;<span class="comment">%cm</span>
fig_h = 5;<span class="comment">%cm</span>
set(fig_hdl,<span class="string">'PaperUnits'</span>,<span class="string">'centimeters'</span>)
set(fig_hdl,<span class="string">'PaperPosition'</span>,[0 0 fig_w fig_h])
set(fig_hdl, <span class="string">'PaperSize'</span>, [fig_w fig_h]);
print(fig_hdl,<span class="string">'-dpng'</span>,<span class="string">'-r500'</span>,fname)
<span class="comment">%crop([fname '.png'])</span>
<span class="keyword">end</span>
<span class="keyword">for</span> i = 2:12
fig_hdl = figure(i);
<span class="keyword">if</span> save_flag
fname = [<span class="string">'Phases_'</span> site <span class="string">'_pair'</span> num2str(i-1) <span class="string">'-'</span> num2str(i)];
fig_w = 8;<span class="comment">%cm</span>
fig_h = 5;<span class="comment">%cm</span>
set(fig_hdl,<span class="string">'PaperUnits'</span>,<span class="string">'centimeters'</span>)
set(fig_hdl,<span class="string">'PaperPosition'</span>,[0 0 fig_w fig_h])
set(fig_hdl, <span class="string">'PaperSize'</span>, [fig_w fig_h]);
print(fig_hdl,<span class="string">'-dpng'</span>,<span class="string">'-r500'</span>,fname)
<span class="comment">% crop([fname '.png'])</span>
<span class="keyword">end</span>
<span class="keyword">end</span>
fig_hdl = figure(13);
<span class="keyword">if</span> save_flag
fname = [<span class="string">'PhaseDifference_'</span> site <span class="string">'_calibrations'</span>];
fig_w = 8;<span class="comment">%cm</span>
fig_h = 5;<span class="comment">%cm</span>
set(fig_hdl,<span class="string">'PaperUnits'</span>,<span class="string">'centimeters'</span>)
set(fig_hdl,<span class="string">'PaperPosition'</span>,[0 0 fig_w fig_h])
set(fig_hdl, <span class="string">'PaperSize'</span>, [fig_w fig_h]);
print(fig_hdl,<span class="string">'-dpng'</span>,<span class="string">'-r500'</span>,fname)
<span class="comment">% crop([fname '.png'])</span>
<span class="keyword">end</span>
fig_hdl = figure(14);
<span class="keyword">if</span> save_flag
fname = [<span class="string">'PhaseCalibration_'</span> site];
fig_w = 8;<span class="comment">%cm</span>
fig_h = 5;<span class="comment">%cm</span>
set(fig_hdl,<span class="string">'PaperUnits'</span>,<span class="string">'centimeters'</span>)
set(fig_hdl,<span class="string">'PaperPosition'</span>,[0 0 fig_w fig_h])
set(fig_hdl, <span class="string">'PaperSize'</span>, [fig_w fig_h]);
print(fig_hdl,<span class="string">'-dpng'</span>,<span class="string">'-r500'</span>,fname)
<span class="comment">% crop([fname '.png'])</span>
<span class="keyword">end</span>
</pre><p class="footer"><br>
Published with MATLAB&reg; 7.14<br></p></div><!--
##### SOURCE BEGIN #####
clear all
close all
%cd('/home/bandma01/Documents/WERAs/')
cd('/home/bandma01/hfr-data-processing/calibrations/')
% site = 'PAO';
site = 'PAB';
cd(site)
% system
for sys = 1:2
if sys == 1
systm = 'Helzel';
elseif sys == 2
systm = 'Flament';
end
clear J I E Phase calib_file calib_hr calib_jd calib_min calib_time calib_value
clear calib_year d data delta_phase_obs delta_phase_th gps_data index_timecalib
clear theta_pairs
if strcmp(site,'PAO') == 1
load('HelzelPAO_Track_29-OCT-14 183138.mat')
elseif strcmp(site,'PAB') == 1
if strcmp(systm,'Helzel') == 1
load('HelzelPAB_Track_30-OCT-14 175446.mat')
elseif strcmp(systm,'Flament') == 1
load('FlamentPAB_Track_30-OCT-14 182008.mat')
end
end
%% phases and time of calibration
if strcmp(site,'PAO') == 1
cd('2014302')
elseif strcmp(site,'PAB') == 1
cd('2014303')
end
% list files
clear d
d = dir('*.mat'); % returns a structure with name, date, bytes, isdir, datenum
Nfiles = length(d); % nombre de fichiers a analyser
% look at the names of files and extract information
for c = 1:Nfiles
calib_year(c) = str2double(d(c).name(1:4));
calib_jd(c) = str2double(d(c).name(5:7));
calib_hr(c) = str2double(d(c).name(8:9));
calib_min(c) = str2double(d(c).name(10:11));
load(d(c).name,'phase')
Phase(:,c) = phase;
end
clear phase
delta_phase_obs = modsym(diff(Phase,1),180);
To = datenum(2014,01,01,0,0,0);
% add 20 sec to center calibration time
calib_time = To + calib_jd -1 + calib_hr./24 + calib_min/24/60 + 20/24/3600;
% datevec(calib_time)
if strcmp(site,'PAO') == 1
% calib time 2min34sec late wrt GPS
calib_time = calib_time + 2/24/60 + 34/24/3600;
end
%% separate between calibrations made with Helzel's system and Pierre's system
% date of new calibration files
% Helzel 2014-10-30 16:52 to 17:54
% Flament 2014-10-30 17:55 to 18:20
T1 = datenum(2014,10,30,17,59,00);
clear J
if strcmp(systm,'Helzel') == 1
J = find(calib_time > T1);
elseif strcmp(systm,'Flament') == 1
J = find(calib_time > T1);
end
%remove unwanted data to avoid confusion!
calib_year(J)=[];
calib_jd(J)=[];
calib_hr(J)=[];
calib_min(J)=[];
Phase(:,J)=[];
delta_phase_obs(:,J)=[];
calib_time(J)=[];
d(J)=[];
clear c Nfiles
%+++++++++++++++++++++++++++++++++++++++++
% for each calibration file (series# 1)
%+++++++++++++++++++++++++++++++++++++++++
for calib_file = 1:length(calib_time)
time_calib = calib_time(calib_file);
% find corresponding file in gps data
clear I
I = find(data(:,1) > time_calib);
% take closest point for now