ntpTimeAdjustSouth

function to correct the 1PPS timing on the roach and the servo to proper
time as referenced to NTP atomic clock.

0001 function d=ntpTimeAdjustSouth(d)
0002 
0003 %function to correct the 1PPS timing on the roach and the servo to proper
0004 %time as referenced to NTP atomic clock.
0005 
0006 
0007 
0008 
0009 
0010 %Find the pps edge and correct for fact that at first 100Hz sample
0011 ppsIndex = find(d.antenna0.roach1.fpgaClockStamp==2499976); %these are the indexes of the data smaples where the PPS edge has triggered and one 0.01s accumulation has also passed by. i.e at what we think is 0.01s after the last seconds
0012 ppsOffset=d.antenna0.roach1.ntpUSeconds(ppsIndex)/1e6-0.01; %we calculate the offset (in fractions of a second) between NTP and 1PPS by looking at the NTP microsecond time at the indexes located above.
0013 
0014  %%polyfit to the ntp time- straight line
0015 ppsPoly = polyfit(ppsIndex,ppsOffset,1);
0016 ppsOffsetFit=polyval(ppsPoly,[1:length(d.antenna0.roach1.utc)]); %fit a straight line to the date using the polynomial above
0017  
0018 
0019 %  plot(ppsIndex,ppsOffset,'r')
0020 %   hold all
0021 %  plot([1:length(d.antenna0.roach1.utc)],offsetPoly)
0022 
0023 %Correct the UTC values using the offset polyvalues.
0024 offsetCorrectedUTC = d.antenna0.roach1.utc - (ppsOffsetFit/86400)';
0025 
0026 offsetCorrectedUTCPoly=polyfit([1:length(offsetCorrectedUTC)],offsetCorrectedUTC',1); %
0027 offsetCorrectedUTCFit=polyval(offsetCorrectedUTCPoly,[1:length(offsetCorrectedUTC)]);
0028 
0029 
0030 data=offsetCorrectedUTC; %rename the data for use in the while loop
0031 dataFit=offsetCorrectedUTCFit; %rename the fitted data for use in the hwile loop
0032 fittingAccuracy=0.1;
0033 diffDecrement=2; %while loop increment/decrement of startDiff. Will divide the startDiff
0034 startDiff=10; %ad hoc starting fitted offset. Anything that is greater that startDiff will be removed from the fit
0035 fittedDiff=86400*(dataFit-data');
0036 diffCalc=fittedDiff;
0037 
0038 counter=1;
0039 close all
0040 dataStart=data;
0041 dataStartIndex=[1:length(offsetCorrectedUTC)];
0042 %plot([1:length(offsetCorrectedUTC)],data,'.')
0043 
0044 while((mean(fittedDiff.^2)>fittingAccuracy) && counter<=10)
0045    indexNew=find(diffCalc<=startDiff);
0046    newPoly=polyfit(indexNew,data(indexNew)',1); %
0047    dataFit=polyval(newPoly,[1:length(offsetCorrectedUTC)]);
0048    fittedDiff=86400*(dataFit(indexNew)-data(indexNew)');
0049    diffCalc=86400*(dataFit-data');
0050    startDiff=startDiff./diffDecrement;
0051    counter=counter+1;
0052  
0053 end
0054 %plot([1:length(offsetCorrectedUTC)],dataFit)
0055 
0056 %close all
0057 
0058 %plot(86400.*(dataFit-data'))
0059 
0060 d.antenna0.roach1.utc=dataFit;
0061 
0062 %%%%%%%%%ROACH2%%%%
0063 %Find the pps edge and correct for fact that at first 100Hz sample
0064 ppsIndex = find(d.antenna0.roach2.fpgaClockStamp==2499976); %these are the indexes of the data smaples where the PPS edge has triggered and one 0.01s accumulation has also passed by. i.e at what we think is 0.01s after the last seconds
0065 ppsOffset=d.antenna0.roach2.ntpUSeconds(ppsIndex)/1e6-0.01; %we calculate the offset (in fractions of a second) between NTP and 1PPS by looking at the NTP microsecond time at the indexes located above.
0066 
0067  %%polyfit to the ntp time- straight line
0068 ppsPoly = polyfit(ppsIndex,ppsOffset,1);
0069 ppsOffsetFit=polyval(ppsPoly,[1:length(d.antenna0.roach2.utc)]); %fit a straight line to the date using the polynomial above
0070  
0071 
0072 %  plot(ppsIndex,ppsOffset,'r')
0073 %   hold all
0074 %  plot([1:length(d.antenna0.roach2.utc)],offsetPoly)
0075 
0076 %Correct the UTC values using the offset polyvalues.
0077 offsetCorrectedUTC = d.antenna0.roach2.utc - (ppsOffsetFit/86400)';
0078 
0079 offsetCorrectedUTCPoly=polyfit([1:length(offsetCorrectedUTC)],offsetCorrectedUTC',1); %
0080 offsetCorrectedUTCFit=polyval(offsetCorrectedUTCPoly,[1:length(offsetCorrectedUTC)]);
0081 
0082 
0083 data=offsetCorrectedUTC; %rename the data for use in the while loop
0084 dataFit=offsetCorrectedUTCFit; %rename the fitted data for use in the hwile loop
0085 fittingAccuracy=0.1;
0086 diffDecrement=2; %while loop increment/decrement of startDiff. Will divide the startDiff
0087 startDiff=10; %ad hoc starting fitted offset. Anything that is greater that startDiff will be removed from the fit
0088 fittedDiff=86400*(dataFit-data');
0089 diffCalc=fittedDiff;
0090 
0091 counter=1;
0092 close all
0093 dataStart=data;
0094 dataStartIndex=[1:length(offsetCorrectedUTC)];
0095 %plot([1:length(offsetCorrectedUTC)],data,'.')
0096 
0097 while((mean(fittedDiff.^2)>fittingAccuracy) && counter<=10)
0098    indexNew=find(diffCalc<=startDiff);
0099    newPoly=polyfit(indexNew,data(indexNew)',1); %
0100    dataFit=polyval(newPoly,[1:length(offsetCorrectedUTC)]);
0101    fittedDiff=86400*(dataFit(indexNew)-data(indexNew)');
0102    diffCalc=86400*(dataFit-data');
0103    startDiff=startDiff./diffDecrement;
0104    counter=counter+1;
0105  
0106 end
0107 %plot([1:length(offsetCorrectedUTC)],dataFit)
0108 
0109 %close all
0110 
0111 %plot(86400.*(dataFit-data'))
0112 
0113 d.antenna0.roach2.utc=dataFit;
0114 
0115 %%%%%%%%%% SERVO Timing
0116 UTC=d.antenna0.servo.utc;
0117 UTCPoly=polyfit([1:length(UTC)],UTC',1); %
0118 UTCFit=polyval(UTCPoly,[1:length(UTC)]);
0119 
0120 
0121 data=UTC; %rename the data for use in the while loop
0122 dataFit=UTCFit; %rename the fitted data for use in the hwile loop
0123 fittingAccuracy=0.1;
0124 diffDecrement=2; %while loop increment/decrement of startDiff. Will divide the startDiff
0125 startDiff=10; %ad hoc starting fitted offset. Anything that is greater that startDiff will be removed from the fit
0126 fittedDiff=86400*(dataFit-data');
0127 diffCalc=fittedDiff;
0128 
0129 counter=1;
0130 close all
0131 dataStart=data;
0132 dataStartIndex=[1:length(UTC)];
0133 %plot([1:length(offsetCorrectedUTC)],data,'.')
0134 
0135 while((mean(fittedDiff.^2)>fittingAccuracy) && counter<=10)
0136    indexNew=find(diffCalc<=startDiff);
0137    newPoly=polyfit(indexNew,data(indexNew)',1); %
0138    dataFit=polyval(newPoly,[1:length(UTC)]);
0139    fittedDiff=86400*(dataFit(indexNew)-data(indexNew)');
0140    diffCalc=86400*(dataFit-data');
0141    startDiff=startDiff./diffDecrement;
0142    counter=counter+1;
0143  
0144 end
0145 %plot([1:length(offsetCorrectedUTC)],dataFit)
0146 
0147 %close all
0148 
0149 %plot(86400.*(dataFit-data'))
0150 
0151 d.antenna0.servo.utc=dataFit;
0152 
0153 
0154