trackTest

 d = read_arc('24-sep-2009:22:05:00', '25-sep-2009:02:00:00');
 d = read_arc('28-sep-2009:22:10:00', '29-sep-2009:06:24:27');
d =read_arc('16-oct-2009:19:35:00', '16-oct-2009:20:00:00');

0001 function [azErr, elErr, az, el] = trackTest(d, secOff)
0002 
0003   % d = read_arc('24-sep-2009:22:05:00', '25-sep-2009:02:00:00');
0004 % d = read_arc('28-sep-2009:22:10:00', '29-sep-2009:06:24:27');
0005 %d =read_arc('16-oct-2009:19:35:00', '16-oct-2009:20:00:00');
0006 
0007 
0008 % the locations in servo.utc are actually from teh previous second
0009 d.antenna0.servo.utc = d.antenna0.servo.utc+(secOff/(24*60*60));
0010 azRec = d.antenna0.servo.fast_az_pos;
0011 elRec = d.antenna0.servo.fast_el_pos;
0012 azErr = d.antenna0.servo.fast_az_err;
0013 elErr = d.antenna0.servo.fast_el_err;
0014 
0015 
0016 % need to re-interpolate the utc for the refraction corrections.
0017 Avals = d.antenna0.tracker.refraction(:,1);
0018 Bvals = d.antenna0.tracker.refraction(:,2);
0019 A = interp1(d.array.frame.utc, Avals, d.antenna0.servo.utc);
0020 B = interp1(d.array.frame.utc, Bvals, d.antenna0.servo.utc);
0021 % apply refraction corrections
0022 Pvals = d.array.weather.pressure;
0023 Tvals = d.array.weather.airTemperature;
0024 P = interp1(d.array.frame.utc, Pvals, d.antenna0.servo.utc);
0025 T = interp1(d.array.frame.utc, Tvals, d.antenna0.servo.utc);
0026 
0027 
0028 % now we calculate the lst
0029 jd = d.antenna0.servo.utc + 2400000.5;
0030 jd0 = floor(d.antenna0.servo.utc) + 2400000.5;
0031 H = 24*(jd - jd0);
0032 D = jd - 2451545.0;
0033 D0 = jd0 - 2451545.0;
0034 
0035 long=-118.2822;
0036 lat=37.2339;
0037 el=1222.00;
0038 
0039 % now we get the LST:
0040 GMT = 6.697374558 + 0.06570982441908*D0 + 1.00273790935*H;
0041 GMT = rem(GMT, 24);
0042 eps = 23.4393 - 0.0000004*D;
0043 L = 280.47 + 0.98565*D;
0044 omega = 125.04 - 0.052954*D;
0045 deltaPsi = -0.000319*sin(omega) - 0.000024*sin(2*L);
0046 eqeq = deltaPsi.*cos(eps);
0047 GAST = GMT + eqeq;
0048 lst2 = GAST + long/15;
0049 lst2(lst2<0) = lst2(lst2<0)+24;
0050 lst = lst2*15;
0051 
0052 % now for hour angle:
0053 ra = median(d.antenna0.tracker.equat_geoc(:,1))*15; % in degrees
0054 dec = median(d.antenna0.tracker.equat_geoc(:,2));  % in degrees
0055 ha = lst - ra;  % also in degrees
0056 
0057 % convert all to radians:
0058 HA = ha*pi/180;
0059 de = dec*pi/180;
0060 La = lat*pi/180;
0061 
0062 [az, el] = hdl2ae(HA,de,La);
0063 
0064 % apply refraction correction.
0065 R = 0.00452.*P./( (273.15+T ).*tan(el))*pi/180;
0066 el = el+R;
0067 
0068 az = az*180/pi;
0069 el = el*180/pi;
0070 az(az<0) = az(az<0)+360;
0071 % we need to add the online pointing model applied
0072 model = [-0.1155, -0.0976, -0.0245, -0.0277, -0.3234, 0.1846, 0, ...
0073        -0.5755, -0.1497];
0074 [az, el] = pointing_model(model, az, el);
0075 
0076 % apply the refraction corrections to the elevation.
0077 %el = el.*pi/180;
0078 %sinel = sin(el);
0079 %cosel = cos(el);
0080 
0081 %el2 = (A.*cosel.*(sinel).^3 + B.*sinel.*cosel.^3)./( sinel.^4 + A.*sinel.^2 + 3*B.*cosel.^2 );
0082 %el1 = el - (A.*cot(el) + B.*(cot(el)).^3);
0083 %el2 = el2.*180/pi;
0084 
0085 % naed to apply the encoder zeros:
0086 %az = az + d.antenna0.tracker.encoder_off(1,1);
0087 %el = el + d.antenna0.tracker.encoder_off(1,2);
0088 
0089 % now for plots:
0090 
0091 time = (lst-lst(1))/15;
0092 clf;
0093 figure(1)
0094 plot( time,(elRec-el)*60, 'k');
0095 xlabel('time');
0096 ylabel('Error, arcminutes');
0097 title('Elevation Error');
0098 
0099 figure(2)
0100 plot( time,(azRec-az)*60, 'r');
0101 %legend('Elevation', 'Azimuth');
0102 xlabel('time');
0103 ylabel('Error, arcminutes');
0104 title('Azimuth Error');
0105 
0106 
0107 nanmean(elRec-el)*60*60
0108 nanmean(azRec-az)*60*60
0109 
0110