automaticSouthColdLoad

% When do you want to read from, till and which schedules do you want to pull out
arcDir = '/data/arc'

0001 %% When do you want to read from, till and which schedules do you want to pull out
0002 %arcDir = '/data/arc'
0003 arcDir = '/elephant/CBASS_SOUTHARC/arc'
0004 date_start = '01-Jan-2013';
0005 date_end = '03-Jun-2014';
0006 path_to_redscript='reduc/redScript.red';
0007 sched_to_reduce = 'coldLoadStep';
0008 %sched_to_reduce = 'source_raster_';
0009 
0010 start_date_num = datenum(date_start);
0011 end_date_num = datenum(date_end);
0012 
0013 %% Select schedules which fall on the correct date
0014 [home,installeddir]=where_am_i();
0015 location = [home,'/',installeddir];
0016 %file = [location,'/','log/obs_log.html'];
0017 %file = '/home/cbassuser/cbass_analysis/log/obs_logSa.html'; %CBASS DELL
0018 file = '/elephant/cbass_analysisSa/log/obs_logSa.html'
0019 fid = fopen(file,'r');
0020 
0021 if fid < 0
0022     disp(['File failed to open: ',file]);
0023 end
0024 
0025 start_dates = {};
0026 end_dates = {};
0027 schedule = {};
0028 source_name = {};
0029 temp = 0;
0030 
0031 for j=0:(end_date_num - start_date_num)
0032     
0033     d_string = datestr(start_date_num + j);
0034     date = d_string;
0035     disp( date )
0036     fid = fopen(file,'r');
0037     
0038     i=0;
0039     
0040     while 1
0041         
0042         tline = fgetl(fid);
0043         if(~ischar(tline))
0044             %disp(['Closing file']);
0045             ST = fclose(fid);
0046             break;
0047         end;
0048         
0049         % Modified to ensure that the date appears in the START time.
0050         dloc = regexp(tline,date);
0051         if isempty(dloc)
0052             continue
0053         end
0054         
0055         % Well, is it?
0056         if (dloc(1) < 10)
0057             
0058             if ( length(regexp(tline,date)) ~=0)
0059                 i=i+1;
0060                 % disp(tline);
0061                 splitline = regexp(tline, '\s+', 'split');
0062                 % If the second element contains a / it is a pathname rather than a time
0063                 % i.e. no end time was written becuase the control system crashed.
0064                 if (~isempty(strfind(char(splitline(2)), '/')) )
0065                     start_dates(end+1) =  {'Error'};
0066                     end_dates(end+1) =  {'Error'};
0067                     schedule(end+1)  =  {'Error'};
0068                     source_name(end+1)={'Error'};
0069                 else
0070                     start_dates(end+1) =  splitline(1);
0071                     end_dates(end+1)  =  splitline(2);
0072                     schedule(end+1)  =  splitline(3);
0073                     % Grab source name if it is a raster or source scan
0074 %                     if((strfind(char(schedule(end)),'source')))
0075 %                         source_name(end+1) = splitline(6);
0076 %                     else
0077 %                         source_name(end+1) = schedule(end);% Just use sched name
0078 %                     end
0079                     temp = i+1;
0080                 end
0081             end
0082             
0083             
0084         end
0085     end
0086 end
0087 %% Sort out empty cells (If there are any)
0088 % find empty cells
0089 emptyCells = cellfun(@isempty,schedule);
0090 % remove empty cells
0091 schedule(emptyCells) = [];
0092 start_dates(emptyCells) = [];
0093 end_dates(emptyCells) = [];
0094 
0095 %% Select which of these schedules are the correct ones
0096 %sched_to_reduce = 'coldLoadStep';
0097 
0098 start_string = {};
0099 end_string = {};
0100 sched_string = {};
0101 
0102 for i=1:length(start_dates)
0103     
0104     sched_bits = regexp(char(schedule(i)), '/', 'split');
0105     sched_name = char(sched_bits(length(sched_bits)));
0106     % strip out possible open bracket if there is one.
0107     sched_name = strrep(sched_name, '(','');
0108     
0109     % Check for error state.
0110     if strcmp(char(start_dates(i)),'Error')
0111         continue
0112     end
0113     
0114     if (exist('sched_to_reduce','var'))
0115         stripped_sched_to_reduce=sched_to_reduce;
0116         if (~isempty(regexp(sched_to_reduce,'^~')))
0117             stripped_sched_to_reduce = regexprep(sched_to_reduce,'^~','');
0118         end
0119         
0120         does_sched_match = ~isempty(regexp(sched_name,stripped_sched_to_reduce,'match'));
0121         if (does_sched_match) % This is a case sensitive regex style match
0122             start_string{end+1} =  strrep(char(start_dates(i)),'/','');
0123             end_string{end+1} =  strrep(char(end_dates(i)),'/','');
0124             sched_string{end+1} = char(schedule(i));
0125             continue;
0126         end;
0127     end;
0128     
0129 end
0130 
0131 start_string
0132 
0133 
0134 %% Begin finding the balance temp...
0135 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
0136 %% Cold Load Step Balance
0137 % Get params we want
0138 imax = length( start_string );
0139 nchan = 128;
0140 TtempL = zeros( imax, nchan );
0141 TtempR = zeros( imax, nchan );
0142 for i=1:imax
0143     % Read in the data
0144     if ( tstr2mjd(char(end_string(i))) - tstr2mjd(char(start_string(i))) )<0.125
0145         
0146         a = char( [start_string(i),end_string(i)] )
0147         if tstr2mjd(a(1,:))<tstr2mjd('25-Oct-2013:12:00:00')
0148             d = read_arcSouthColdLoad1( a(1,:),a(2,:),'',arcDir )
0149             LLt = [ d.antenna0.roach1.LL fliplr(d.antenna0.roach2.LL) ];
0150             RRt = [ d.antenna0.roach1.RR fliplr(d.antenna0.roach2.RR) ];
0151             load1t = [ d.antenna0.roach1.load1 fliplr(d.antenna0.roach2.load1) ];
0152             load2t = [ d.antenna0.roach1.load2 fliplr(d.antenna0.roach2.load2) ];
0153             loadT = d.antenna0.thermal.ccTemperatureLoad;
0154             interp_ratio = length(LLt)/length(loadT);
0155             LL(1:1:length(loadT),1:128) = LLt( 1:interp_ratio:length(LLt(:,1)),1:128 );
0156             RR(1:1:length(loadT),1:128) = RRt( 1:interp_ratio:length(RRt(:,1)),1:128 );
0157             load1(1:1:length(loadT),1:128) = load1t( 1:interp_ratio:length(load1t(:,1)),1:128 );
0158             load2(1:1:length(loadT),1:128) = load2t( 1:interp_ratio:length(load2t(:,1)),1:128 );
0159         else
0160             d = read_arcSouthColdLoad( a(1,:),a(2,:),'',arcDir )
0161             LL = [ d.antenna0.roach1.LL fliplr(d.antenna0.roach2.LL) ];
0162             RR = [ d.antenna0.roach1.RR fliplr(d.antenna0.roach2.RR) ];
0163             load1 = [ d.antenna0.roach1.load1 fliplr(d.antenna0.roach2.load1) ];
0164             load2 = [ d.antenna0.roach1.load2 fliplr(d.antenna0.roach2.load2) ];
0165             loadT = d.antenna0.thermal.ccTemperatureLoadFast;
0166         end
0167         indMax = length( loadT );
0168         
0169         % LL - Fit a straight line to the data
0170         x = [ 1:1:length(LL(:,1)) ]';
0171         for chan=1:nchan
0172             difference = LL(:,chan)-load2(:,chan);
0173             p = polyfit( x,difference,1 );
0174             ind = -p(2)/p(1);
0175             ind = round( ind );
0176             if ind<indMax & ind>0
0177                 tempL(chan) = loadT(ind);
0178             else
0179                 tempL(chan) = NaN;
0180             end
0181         end
0182         
0183         % RR - Fit a straight line to the data
0184         x = [ 1:1:length(RR(:,1)) ]';
0185         for chan=1:nchan
0186             difference = RR(:,chan)-load1(:,chan);
0187             p = polyfit( x,difference,1 );
0188             ind = -p(2)/p(1);
0189             ind = round( ind );
0190             if ind<indMax & ind>0
0191                 tempR(chan) = loadT(ind);
0192             else
0193                 tempR(chan) = NaN;
0194             end
0195         end
0196        
0197         % TtempL/R how they change with mjd
0198         TtempL(i,:) = tempL;
0199         TtempR(i,:) = tempR;
0200         mjd(i) = tstr2mjd( a(1,:) );
0201         
0202     else
0203         disp('Rejected schedule, too long')
0204         TtempL(i,:) = NaN;
0205         TtempR(i,:) = NaN;
0206         mjd(i) = tstr2mjd( a(1,:) );
0207     end
0208 end
0209 
0210 %% Prepare To Plot The Data
0211 TtempL(TtempL==0)=NaN;
0212 TtempR(TtempR==0)=NaN;
0213 [year1 month1 day1] = mjd2date(mjd);
0214 dateno = datenum( year1, month1, day1 );
0215 chpl = [1:1:128];
0216 
0217 X = repmat( chpl, 1, length(dateno) );
0218 tempY = repmat( dateno', 1, length(chpl) );
0219 Y = reshape( tempY',1, numel(tempY) );
0220 ZL = reshape(TtempL',1,numel(TtempL) );
0221 ZR = reshape(TtempR',1,numel(TtempR) );
0222 
0223 % Old Plotting
0224 % figure
0225 % %subplot(1,2,1)
0226 % scatter3(X, Y, isnan(ZL), 5, 'fill','sk')
0227 % hold
0228 % scatter3(X, Y, ZL, 20, ZL,'fill','s')
0229 % xlabel( 'Channel No' )
0230 % ylabel( 'Date' )
0231 % zlabel( 'Temp' )
0232 % title( 'LL and Load 2 Balance' )
0233 % xlim( [0 128] )
0234 % colorbar('location','westoutside')
0235 % caxis([8 30])
0236 % view(2)
0237 % datetick('y',29)
0238 % camroll(90)
0239 %
0240 % hold off
0241 %
0242 % figure
0243 % %subplot(1,2,2)
0244 % scatter3(X, Y, isnan(ZR), 5, 'fill','sk')
0245 % hold
0246 % scatter3(X, Y, ZR, 20, ZR,'fill','s')
0247 % xlabel( 'Channel No' )
0248 % ylabel( 'Date' )
0249 % zlabel( 'Temp' )
0250 % title( 'RR and Load 1 Balance' )
0251 % xlim( [0 128] )
0252 % colorbar('location','westoutside')
0253 % caxis([8 30])
0254 % view(2)
0255 % datetick('y',29)
0256 % camroll(90)
0257 
0258 
0259 
0260 %% New Plotting
0261 
0262 figure
0263 scatter3(Y, X, isnan(ZL'), 5, 'fill','sk')
0264 ylabel( 'Channel No' )
0265 xlabel( 'Date' )
0266 zlabel( 'Temp' )
0267 title( 'LL and Load 2 Balance' )
0268 datetick('x',29,'keepticks','keeplimits')
0269 ylim( [0 128] )
0270 view(2)
0271 hold
0272 scatter3(Y, X, ZL', 20, ZL','fill','s')
0273 colorbar
0274 caxis([8 30])
0275 xlim( [min(Y) max(Y)] )
0276 
0277 figure
0278 scatter3(Y, X, isnan(ZR'), 5, 'fill','sk')
0279 ylabel( 'Channel No' )
0280 xlabel( 'Date' )
0281 zlabel( 'Temp' )
0282 title( 'RR and Load 1 Balance' )
0283 datetick('x',29,'keepticks','keeplimits')
0284 ylim( [0 128] )
0285 view(2)
0286 hold
0287 scatter3(Y, X, ZR', 20, ZL','fill','s')
0288 colorbar
0289 caxis([8 30])
0290 xlim( [min(Y) max(Y)] )
0291 
0292 
0293 save(['loadBalTemp_' date_start '_' date_end '.mat'],'X','Y','ZL','ZR')