From mboxrd@z Thu Jan 1 00:00:00 1970 Return-Path: Received: from mailout-de.gmx.net (mailout-de.gmx.net [213.165.64.23]) by huchra.bufferbloat.net (Postfix) with SMTP id 6CDB620029B for ; Wed, 1 Aug 2012 15:58:20 -0700 (PDT) Received: (qmail invoked by alias); 01 Aug 2012 22:58:17 -0000 Received: from tsaolab-fw.caltech.edu (EHLO [192.168.50.16]) [131.215.9.89] by mail.gmx.net (mp029) with SMTP; 02 Aug 2012 00:58:17 +0200 X-Authenticated: #24211782 X-Provags-ID: V01U2FsdGVkX1+qzJYuvGiq0H012Q9VtMzxjS0HWUMFAwypC90gdV b68ve7y/jdvMru Mime-Version: 1.0 (Apple Message framework v1278) Content-Type: multipart/mixed; boundary="Apple-Mail=_D8AF6F97-A1E9-4ED9-A34F-AF2A1E299DF6" From: Sebastian Moeller In-Reply-To: Date: Wed, 1 Aug 2012 15:58:15 -0700 Message-Id: References: To: William Katsak X-Mailer: Apple Mail (2.1278) X-Y-GMX-Trusted: 0 Cc: "" Subject: Re: [Cerowrt-devel] ADSL Issue (PPPoE) X-BeenThere: cerowrt-devel@lists.bufferbloat.net X-Mailman-Version: 2.1.13 Precedence: list List-Id: Development issues regarding the cerowrt test router project List-Unsubscribe: , List-Archive: List-Post: List-Help: List-Subscribe: , X-List-Received-Date: Wed, 01 Aug 2012 22:58:21 -0000 --Apple-Mail=_D8AF6F97-A1E9-4ED9-A34F-AF2A1E299DF6 Content-Transfer-Encoding: quoted-printable Content-Type: text/plain; charset=windows-1252 Hi William, On Aug 1, 2012, at 5:34 AM, William Katsak wrote: > Thanks for the insight. I would definitely be interested in your > program for measuring the overhead. Ah great. So here is what to do: 0) IUse "traceroute 8.8.8.8" to figure out the IP address of the first = hop on the other side of your modem (but see below) 1) from my MacBooks terminal I edit and run ping_sweeper_02.sh (inlined = below) which talks around 6 hours (I just run this overnight) assign the first outside hop's IP address to TARGET #! /bin/bash # ideally one would shuffle all size*repetition to better average out = transient issues # # USEAGE # 1) use "traceroute 8.8.8.8" and find the network hops just outside of = your modem # typically for ATM the modem to DSLAM link will introduce an = additional >10ms latency=20 # so the hop directly after the first big latency jump should be = close to the DSLAM, use=20 # this hops IP address should be used as TARGET below # 2) optional confirm this IPs reliability by running "ping -c 100" = against this IP # the min and max of the ping times should stay close to the = average and the stddev # should be small as well, otherwise PINGSPERSIZE might need to = be set to >=3D100 # TECH=3DATM # CABLE, ATM, ... just to name the log = file DATESTR=3D`date +%Y%m%d_%H%M%S` # to allow multiple sequential records TARGET=3D96.34.97.78 # use traceroute something far to get the first = hop out of the home net LOG=3Dping_sweep_${TECH}_${DATESTR}.txt PINGSPERSIZE=3D50 # 480 packet sizes take 13.33 hours at = 100 repetitions, so stick to 50 SWEEPMINSIZE=3D16 # on macosx 10.7 64bit, sizes < 16 do = not have a timestamp... SWEEPMAXSIZE=3D496 # initial size plus 10 full ATM cells (of 48 = bytes payload each) n_SWEEPS=3D`expr ${SWEEPMAXSIZE} - ${SWEEPMINSIZE}` # do it echo "Sweeping ${n_SWEEPS} * ${PINGSPERSIZE} pings to ${TARGET} might = take a while..." echo "ping -c ${PINGSPERSIZE} -g ${SWEEPMINSIZE} -G ${SWEEPMAXSIZE} = ${TARGET} > ${LOG}" ping -c ${PINGSPERSIZE} -g ${SWEEPMINSIZE} -G ${SWEEPMAXSIZE} ${TARGET} = > ${LOG} & # show some progress tail -f ${LOG} echo "Done... ($0)" 3) run the following in matlab or octave (tested with octave 3.6) = (inlined and attached) as = tc_stab_parameter_guide_02('filly_qualified_path_to_the_ping_log', = your_uplink_rate_in_Kbit, your_downlink_rate_in_Kbit) and then wait a = bit (it might take a minute or two) function [ output_args ] =3D tc_stab_parameter_guide_02( sweep_fqn, = up_Kbit, down_Kbit ) %TC_STAB_PARAMETER_GUIDE Summary of this function goes here % try to read in the result from a ping sweep run % sweep_fqn: the log file of the ping sweep against the first hop = after % the DSL link % up_Kbit: the uplink rate in Kilobits per second % down_Kbit: the downlink rate in Kilobits per second % % TODO: % find whether the carrier is ATM quantized % if yes: % what is the RTT step (try to deduce the combined up and = down rates from this) % try to figure out the overhead for each packet % make sure all sizes are filled (use NaN for emty ones???) %Thoughts: % maybe require to give the nominal up and down rates, to estimate = the % RTT stepsize % ask about IPv4 or IPv6 (what about tunneling?) % the sweep should be taken directly connected to the modem to = reduce % non-ATM routing delays dbstop if error; if ~(isoctave) timestamps.(mfilename).start =3D tic; else tic(); end disp(['Starting: ', mfilename]); % control options show_mean =3D 0; % the means are noisier than the medians show_median =3D 1; % the mean is the way to go show_min =3D 1; % the min should be the best measure, but in the = ATM test sweep is too variable show_max =3D 0; % only useful for debugging show_sem =3D 0; % give some estimate of the variance show_ci =3D 1; % show the confidence interval of the mean ci_alpha =3D 0.05; % alpha for confidence interval calculation use_measure =3D 'median'; use_processed_results =3D 1; if (nargin =3D=3D 0) sweep_fqn =3D fullfile(pwd, 'ping_sweep_ATM.txt'); % was = Bridged, LLC/SNAP RFC-1483/2684 connection (overhead 32 bytes - 14 =3D = 18) % sweep_fqn =3D fullfile(pwd, = 'ping_sweep_CABLE_20120426_230227.txt'); % sweep_fqn =3D fullfile(pwd, = 'ping_sweep_CABLE_20120801_001235.txt'); up_Kbit =3D 512; down_Kbit =3D 3008; end if (nargin =3D=3D 1) up_Kbit =3D 512; down_Kbit =3D 3008; end if (nargin =3D=3D 2) down_Kbit =3D 3008; end %ATM quantum.byte =3D 48; % ATM packets are always 53 bytes, 48 thereof = payload quantum.bit =3D quantum.byte * 8; ATM_cell.byte =3D 53; ATM_cell.bit =3D ATM_cell.byte * 8; % CONFIGURE THE NEXT TWO VALUES line.down.Kbit =3D down_Kbit; line.up.Kbit =3D up_Kbit; % DONE line.down.bit =3D line.down.Kbit * 1024; line.up.bit =3D line.up.Kbit * 1024; % known packet size offsets in bytes offsets.IPv4 =3D 20; % assume no IPv4 options are used, IP 6 = would be 40bytes? offsets.IPv6 =3D 40; % not used yet... offsets.ICMP =3D 8; % ICMP header offsets.ethernet =3D 14; % ethernet header % unknown offsets is what we need to figure out to feed tc-stab... [sweep_dir, sweep_name] =3D fileparts(sweep_fqn); cur_parsed_data_mat =3D [sweep_fqn(1:end-4), '.mat']; if (use_processed_results && ~isempty(dir(cur_parsed_data_mat))) disp(['Loading processed ping data from ', = cur_parsed_data_mat]); load(cur_parsed_data_mat, 'ping'); else % read in the result from a ping sweep disp(['Processing ping data from ', sweep_fqn]); ping =3D parse_ping_output(sweep_fqn); save(cur_parsed_data_mat, 'ping'); end % analyze the data min_ping_size =3D min(ping.data(:, ping.cols.size)) - offsets.ICMP; disp(['Minimum size of ping payload used: ', num2str(min_ping_size), ' = bytes.']); known_overhead =3D offsets.ICMP + min_ping_size + offsets.IPv4; ping.data(:, ping.cols.size) =3D ping.data(:, ping.cols.size) + = known_overhead; % we know we used ping so add the 8 bytes already (since = these will not account for overhead) size_list =3D unique(ping.data(:, ping.cols.size)); per_size.header =3D {'size', 'mean', 'median', 'min', 'max', 'std', 'n', = 'sem', 'ci'}; per_size.cols =3D get_column_name_indices(per_size.header); per_size.data =3D zeros([length(size_list), length(per_size.header)]); per_size.data(:, per_size.cols.size) =3D size_list; for i_size =3D 1 : length(size_list) cur_size =3D size_list(i_size); cur_size_idx =3D find(ping.data(:, ping.cols.size) =3D=3D = cur_size); per_size.data(i_size, per_size.cols.mean) =3D = mean(ping.data(cur_size_idx, ping.cols.time)); per_size.data(i_size, per_size.cols.median) =3D = median(ping.data(cur_size_idx, ping.cols.time)); per_size.data(i_size, per_size.cols.min) =3D = min(ping.data(cur_size_idx, ping.cols.time)); per_size.data(i_size, per_size.cols.max) =3D = max(ping.data(cur_size_idx, ping.cols.time)); per_size.data(i_size, per_size.cols.std) =3D = std(ping.data(cur_size_idx, ping.cols.time), 0); per_size.data(i_size, per_size.cols.n) =3D length(cur_size_idx); per_size.data(i_size, per_size.cols.sem) =3D = per_size.data(i_size, per_size.cols.std) / sqrt(length(cur_size_idx));=09= per_size.data(i_size, per_size.cols.ci) =3D = calc_cihw(per_size.data(i_size, per_size.cols.std), = per_size.data(i_size, per_size.cols.n), ci_alpha);=09 end figure('Name', sweep_name); hold on; legend_str =3D {}; if (show_mean) % means legend_str{end + 1} =3D 'mean'; plot(per_size.data(:, per_size.cols.size), per_size.data(:, = per_size.cols.mean), 'Color', [0 1 0 ]); if (show_sem) legend_str{end + 1} =3D '+sem'; legend_str{end + 1} =3D '-sem'; plot(per_size.data(:, per_size.cols.size), = per_size.data(:, per_size.cols.mean) - per_size.data(:, = per_size.cols.sem), 'Color', [0 0.66 0]); plot(per_size.data(:, per_size.cols.size), = per_size.data(:, per_size.cols.mean) + per_size.data(:, = per_size.cols.sem), 'Color', [0 0.66 0]); end if (show_ci) legend_str{end + 1} =3D '+ci'; legend_str{end + 1} =3D '-ci'; plot(per_size.data(:, per_size.cols.size), = per_size.data(:, per_size.cols.mean) - per_size.data(:, = per_size.cols.ci), 'Color', [0 0.37 0]); plot(per_size.data(:, per_size.cols.size), = per_size.data(:, per_size.cols.mean) + per_size.data(:, = per_size.cols.ci), 'Color', [0 0.37 0]); end end if(show_median) % median +- standard error of the mean, confidence interval = would be % better legend_str{end + 1} =3D 'median'; plot(per_size.data(:, per_size.cols.size), per_size.data(:, = per_size.cols.median), 'Color', [1 0 0]); if (show_sem) legend_str{end + 1} =3D '+sem'; legend_str{end + 1} =3D '-sem'; plot(per_size.data(:, per_size.cols.size), = per_size.data(:, per_size.cols.median) - per_size.data(:, = per_size.cols.sem), 'Color', [0.66 0 0]); plot(per_size.data(:, per_size.cols.size), = per_size.data(:, per_size.cols.median) + per_size.data(:, = per_size.cols.sem), 'Color', [0.66 0 0]); end if (show_ci) legend_str{end + 1} =3D '+ci'; legend_str{end + 1} =3D '-ci'; plot(per_size.data(:, per_size.cols.size), = per_size.data(:, per_size.cols.median) - per_size.data(:, = per_size.cols.ci), 'Color', [0.37 0 0]); plot(per_size.data(:, per_size.cols.size), = per_size.data(:, per_size.cols.median) + per_size.data(:, = per_size.cols.ci), 'Color', [0.37 0 0]); end end if(show_min) % minimum, should be cleanest, but for the test data set looks = quite sad... legend_str{end + 1} =3D 'min'; plot(per_size.data(:, per_size.cols.size), per_size.data(:, = per_size.cols.min), 'Color', [0 0 1]); end if(show_max) % minimum, should be cleanest, but for the test data set looks = quite sad... legend_str{end + 1} =3D 'max'; plot(per_size.data(:, per_size.cols.size), per_size.data(:, = per_size.cols.max), 'Color', [0 0 0.66]); end title(['If this plot shows a (noisy) step function with a stepping of ', = num2str(quantum.byte), ' bytes then the data carrier is quantised, make = sure to use tc-stab']); xlabel('Approximate packet size [bytes]'); ylabel('ICMP round trip times (ping RTT) [ms]'); legend(legend_str); hold off; % looking at the different plot we see many side peaks %figure; %plot(per_size.data(1:end - 1, strmatch('size', per_size.header , = 'exact')), diff(per_size.data(:, strmatch('median', per_size.header , = 'exact'))), 'Color', [1 0 0]); % potentially clean up the data, by interpolating values with large sem % from the neighbours? % if this is ATM based the expectancy is there are RTT time quantums of = 48 % bytes (the payload of an ATM package) find this %x =3D fminsearch(@(x)sin(x^2), x0); % estimate the RTT step size % at ADSL down 3008kbit/sec up 512kbit/sec we expect, this oes not = include % processing time expected_RTT_quantum_ms =3D (ATM_cell.bit / line.down.bit + ATM_cell.bit = / line.up.bit ) * 1000; % this estimate is rather a lower bound for = fastpath , so search for best fits disp(['lower bound estimate of one ATM cell RTT is ', = num2str(expected_RTT_quantum_ms), ' ms.']); % lets search from expected_RTT_quantum_ms to 1.5 * = expected_RTT_quantum_ms % in steps of expected_RTT_quantum_ms / 100 % to allow for interleaved ATM set ups increase the search space up to = 32 % times best fastpath RTT estimate RTT_quantum_list =3D (expected_RTT_quantum_ms : expected_RTT_quantum_ms = / 100 : 32 * expected_RTT_quantum_ms); quantum_list =3D (1 : 1 : quantum.byte); differences =3D zeros([length(RTT_quantum_list) length(quantum_list)]); cumulative_differences =3D differences; all_stairs =3D zeros([length(RTT_quantum_list) length(quantum_list) = length(per_size.data(:, per_size.cols.(use_measure)))]); for i_RTT_quant =3D 1 : length(RTT_quantum_list) cur_RTT_quant =3D RTT_quantum_list(i_RTT_quant); for i_quant =3D 1 : quantum.byte [differences(i_RTT_quant, i_quant), = cumulative_differences(i_RTT_quant, i_quant), all_stairs(i_RTT_quant, = i_quant, :)] =3D get_difference_between_data_and_stair( per_size.data(:, = per_size.cols.size), per_size.data(:, per_size.cols.(use_measure)), ... = = quantum_list(i_quant), quantum.byte, 0, cur_RTT_quant ); end end % for the test DSL set the best x_offset is 21. [min_cum_diff, min_cum_diff_idx] =3D min(cumulative_differences(:)); [min_cum_diff_row_idx, min_cum_diff_col_idx] =3D = ind2sub(size(cumulative_differences),min_cum_diff_idx); best_difference =3D differences(min_cum_diff_row_idx, = min_cum_diff_col_idx); disp(['remaining ATM cell length after ICMP header is ', = num2str(quantum_list(min_cum_diff_col_idx)), ' bytes.']); disp(['ICMP RTT of a single ATM packet is ', = num2str(RTT_quantum_list(min_cum_diff_col_idx)), ' ms.']); figure('Name', 'Comparing ping data with'); hold on legend_str =3D {'ping_data', 'fitted_stair'}; plot(per_size.data(:, per_size.cols.size), per_size.data(:, = per_size.cols.(use_measure)), 'Color', [1 0 0]); plot(per_size.data(:, per_size.cols.size), = squeeze(all_stairs(min_cum_diff_row_idx, min_cum_diff_col_idx, :)) + = best_difference, 'Color', [0 1 0]); title(['Estimated RTT per quantum: ', = num2str(RTT_quantum_list(min_cum_diff_col_idx)), ' ms; ICMP header = offset in quantum ', num2str(quantum_list(min_cum_diff_col_idx)), ' = bytes']); xlabel('Approximate packet size [bytes]'); ylabel('ICMP round trip times (ping RTT) [ms]'); if (isoctave) legend(legend_str); else legend(legend_str, 'Interpreter', 'none'); end hold off % find cumulative_differences minimum, get the differences for that = matrix % and plot the corresponding stair with the real data % next extrapolate to get the y-intercept % as first approximation use the ATM cell offset and known offsets (ICMP % IPv4 min_ping_size) to estimate the number of cells used for per paket % overhead % this assumes that no ATM related overhead is >=3D ATM cell size % -1 to account for matlab 1 based indices n_bytes_overhead_2nd_cell =3D quantum.byte - = (quantum_list(min_cum_diff_col_idx) - 1); % just assume we can not = fit all overhead into one cell... pre_IP_overhead =3D quantum.byte + (n_bytes_overhead_2nd_cell - = known_overhead); % ths is the one we are after in the end disp(' '); disp(['Estimated overhead preceeding the IP header: ', = num2str(pre_IP_overhead), ' bytes']); % use http://ace-host.stuart.id.au/russell/files/tc/tc-atm/ to present = the % most likely ATM setup for a given overhead and present a = recommendation % for the stab invocation display_protocol_stack_information(pre_IP_overhead); % now turn this into tc-stab recommendations: disp(['Add the following to both the ingress and egress root qdiscs:']); disp(' '); disp(['A) Assuming the router connects over ethernet to the = DSL-modem:']); disp(['stab mtu 2048 tsize 128 overhead ', num2str(pre_IP_overhead - = offsets.ethernet), ' linklayer atm']); disp(' '); disp(['B) Assuming the router connects via PPP and non-ethernet to the = modem:']); disp(['stab mtu 2048 tsize 128 overhead ', num2str(pre_IP_overhead), ' = linklayer atm']); disp(' '); if ~(isoctave) timestamps.(mfilename).end =3D = toc(timestamps.(mfilename).start); disp([mfilename, ' took: ', num2str(timestamps.(mfilename).end), = ' seconds.']); else toc end disp('Done...'); return end function [ ping_data ] =3D parse_ping_output( ping_log_fqn ) %PARSE_PING_OUTPUT read the putput of a ping run/sweep % for further processng cur_sweep_fd =3D fopen(ping_log_fqn, 'r'); ping_data.header =3D {'size', 'icmp_seq', 'ttl', 'time'}; ping_data.cols =3D get_column_name_indices(ping_data.header); ping_data.data =3D zeros([1, length(ping_data.header)]); cur_data_lines =3D 0; % skip the first line % PING netblock-75-79-143-1.dslextreme.com (75.79.143.1): (16 ... 1000) % data bytes header_line =3D fgetl(cur_sweep_fd); while ~feof(cur_sweep_fd) cur_line =3D fgetl(cur_sweep_fd); [first_element, remainder] =3D strtok(cur_line); first_element_as_number =3D str2num(first_element); if isempty(first_element); % skip empty lines continue; end if strmatch('---', first_element) %we reached the end of sweeps break; end % now read in the data % 30 bytes from 75.79.143.1: icmp_seq=3D339 ttl=3D63 time=3D14.771= ms if ~isempty(first_element_as_number) cur_data_lines =3D cur_data_lines + 1; % size of the ICMP package ping_data.data(cur_data_lines, ping_data.cols.size) =3D = first_element_as_number; % now process the remainder while ~isempty(remainder) [next_item, remainder] =3D strtok(remainder); equality_pos =3D strfind(next_item, '=3D'); % data items are name+value pairs if ~isempty(equality_pos); cur_key =3D next_item(1: equality_pos - = 1); cur_value =3D = str2num(next_item(equality_pos + 1: end)); switch cur_key % busybox ping and macosx ping = return different key names case {'seq', 'icmp_seq'} = ping_data.data(cur_data_lines, ping_data.cols.icmp_seq) =3D cur_value; case 'ttl' = ping_data.data(cur_data_lines, ping_data.cols.ttl) =3D cur_value; case 'time' = ping_data.data(cur_data_lines, ping_data.cols.time) =3D cur_value; end end end =09 end =09 end % clean up fclose(cur_sweep_fd); return end function [ difference , cumulative_difference, stair_y ] =3D = get_difference_between_data_and_stair( data_x, data_y, x_size, = stair_x_step_size, y_offset, stair_y_step_size ) % x_size is the flat part of the first stair, that is quantum minus the % offset debug =3D 0; difference =3D []; x_start_val =3D min(data_x); x_end_val =3D max(data_x); % construct stair stair_x =3D data_x; proto_stair_y =3D zeros([data_x(end) 1]); % make sure the x_size values do not exceed the step size... if (x_size > stair_x_step_size) if mod(x_size, stair_x_step_size) =3D=3D 0 x_size =3D stair_x_step_size; else x_size =3D mod(x_size, stair_x_step_size); end end stair_y_step_idx =3D (x_start_val + x_size : stair_x_step_size : = x_end_val); proto_stair_y(stair_y_step_idx) =3D stair_y_step_size; stair_y =3D cumsum(proto_stair_y); stair_y =3D stair_y(x_start_val:x_end_val) + y_offset; if (debug) figure hold on; title(['x offset used: ', num2str(x_size), ' with quantum ', = num2str(stair_x_step_size)]); plot(data_x, data_y, 'Color', [0 1 0]); plot(stair_x, stair_y, 'Color', [1 0 0]); hold off; end difference =3D sum(abs(data_y - stair_y)) / length(data_y); cumulative_difference =3D sum(abs(data_y - (stair_y + difference))); return end % function [ stair ] =3D build_stair(x_vector, x_size, = stair_x_step_size, y_offset, stair_y_step_size ) % stair =3D []; %=20 % return % end function [columnnames_struct, n_fields] =3D = get_column_name_indices(name_list) % return a structure with each field for each member if the name_list = cell % array, giving the position in the name_list, then the = columnnames_struct % can serve as to address the columns, so the functions assitgning = values % to the columns do not have to care too much about the positions, and = it % becomes easy to add fields. n_fields =3D length(name_list); for i_col =3D 1 : length(name_list) cur_name =3D name_list{i_col}; columnnames_struct.(cur_name) =3D i_col; end return end function [ci_halfwidth_vector] =3D calc_cihw(std_vector, n, alpha) %calc_ci : calculate the half width of the confidence interval (for 1 - = alpha) % the t_value lookup depends on alpha and the samplesize n; the = relevant % calculation of the degree of freedom is performed inside = calc_t_val. % ci_halfwidth =3D t_val(alpha, n-1) * std / sqrt(n) % Each groups CI ranges from mean - ci_halfwidth to mean - = ci_halfwidth, so % the calling function has to perform this calculation... % % INPUTS: % std_vector: vector containing the standard deviations of all = requested % groups % n: number of samples in each group, if the groups have different % samplesizes, specify each groups samplesize in a vector % alpha: the desired maximal uncertainty/error in the range of [0, = 1] % OUTPUT: % ci_halfwidth_vector: vector containing the confidence intervals = half width % for each group % calc_t_val return one sided t-values, for the desired two sidedness = one has % to half the alpha for the table lookup cur_alpha =3D alpha / 2; % if n is scalar use same n for all elements of std_vec if isscalar(n) t_ci =3D calc_t_val(cur_alpha, n); ci_halfwidth_vector =3D std_vector * t_ci / sqrt(n); % if n is a vector, prepare a matching vector of t_ci values elseif isvector(n) t_ci_vector =3D n; % this is probably ugly, but calc_t_val only accepts scalars. for i_pos =3D 1 : length(n) t_ci_vector(i_pos) =3D calc_t_val(cur_alpha, n(i_pos)); end ci_halfwidth_vector =3D std_vector .* t_ci_vector ./ sqrt(n); end return end = %-------------------------------------------------------------------------= ---- function [t_val] =3D calc_t_val(alpha, n) % the t value for the given alpha and n % so call with the n of the sample, not with degres of freedom % see http://mathworld.wolfram.com/Studentst-Distribution.html for = formulas % return values follow Bortz, Statistik fuer Sozialwissenschaftler, = Springer=20 % 1999, table D page 775. That is it returns one sided t-values. % primary author S. Moeller % TODO: % sidedness of t-value??? % basic error checking if nargin < 2 error('alpha and n have to be specified...'); end % probabilty of error tmp_alpha =3D alpha ;%/ 2; if (tmp_alpha < 0) || (tmp_alpha > 1) msgbox('alpha has to be taken from [0, 1]...'); t_val =3D NaN; return end if tmp_alpha =3D=3D 0 t_val =3D -Inf; return elseif tmp_alpha =3D=3D1 t_val =3D Inf; return end % degree of freedom df =3D n - 1; if df < 1 %msgbox('The n has to be >=3D 2 (=3D> df >=3D 1)...'); % disp('The n has to be >=3D 2 (=3D> df >=3D 1)...'); t_val =3D NaN; return end % only calculate each (alpha, df) combination once, store the results persistent t_val_array; % create the t_val_array if ~iscell(t_val_array) t_val_array =3D {[NaN;NaN]}; end % search for the (alpha, df) tupel, avoid calculation if already stored if iscell(t_val_array) % cell array of 2d arrays containing alpha / t_val pairs if df <=3D length(t_val_array) % test whether the required alpha, t_val tupel exists if ~isempty(t_val_array{df}) % search for alpha tmp_array =3D t_val_array{df}; alpha_index =3D find(tmp_array(1,:) =3D=3D = tmp_alpha); if any(alpha_index) t_val =3D tmp_array(2, alpha_index); return end end else % grow t_val_array to length of n missing_cols =3D df - length(t_val_array); for i_missing_cols =3D 1: missing_cols t_val_array{end + 1} =3D [NaN;NaN]; end end end % check the sign cdf_sign =3D 1; if (1 - tmp_alpha) =3D=3D 0.5 t_val =3D t_cdf; elseif (1 - tmp_alpha) < 0.5 % the t-cdf is point symmetric around (0, = 0.5) cdf_sign =3D -1; tmp_alpha =3D 1 - tmp_alpha; % this will be undone later end % init some variables n_iterations =3D 0; delta_t =3D 1; last_alpha =3D 1; higher_t =3D 50; lower_t =3D 0; % find a t-value pair around the desired alpha value=20 while norm_students_cdf(higher_t, df) < (1 - tmp_alpha); lower_t =3D higher_t; higher_t =3D higher_t * 2; end % search the t value for the given alpha... while (n_iterations < 1000) && (abs(delta_t) >=3D 0.0001) n_iterations =3D n_iterations + 1; % get the test_t (TODO linear interpolation)=20 % higher_alpha =3D norm_students_cdf(higher_t, df); % lower_alpha =3D norm_students_cdf(lower_t, df); test_t =3D lower_t + ((higher_t - lower_t) / 2); cur_alpha =3D norm_students_cdf(test_t, df); % just in case we hit the right t spot on... if cur_alpha =3D=3D (1 - tmp_alpha) t_crit =3D test_t; break; % probably we have to search for the right t elseif cur_alpha < (1 - tmp_alpha) % test_t is the new lower_t lower_t =3D test_t; %higher_t =3D higher_t; % this stays as is... elseif cur_alpha > (1 - tmp_alpha) %=20 %lower_t =3D lower_t; % this stays as is... higher_t =3D test_t; end delta_t =3D higher_t - lower_t; last_alpha =3D cur_alpha; end t_crit =3D test_t; % set the return value, correct for negative t values t_val =3D t_crit * cdf_sign; if cdf_sign < 0 tmp_alpha =3D 1 - tmp_alpha; end % store the alpha, n, t_val tupel in t_val_array pos =3D size(t_val_array{df}, 2); t_val_array{df}(1, (pos + 1)) =3D tmp_alpha; t_val_array{df}(2, (pos + 1)) =3D t_val; return end = %-------------------------------------------------------------------------= ---- function [scaled_cdf] =3D norm_students_cdf(t, df) % calculate the cdf of students distribution for a given degree of = freedom df, % and all given values of t, then normalize the result % the extreme values depend on the values of df!!! % get min and max by calculating values for extrem t-values (e.g. = -10000000, % 10000000) extreme_cdf_vals =3D students_cdf([-10000000, 10000000], df); tmp_cdf =3D students_cdf(t, df); scaled_cdf =3D (tmp_cdf - extreme_cdf_vals(1)) /... (extreme_cdf_vals(2) - = extreme_cdf_vals(1)); return end = %-------------------------------------------------------------------------= ---- function [cdf_value_array] =3D students_cdf(t_value_array, df) %students_cdf: calc the cumulative density function for a t-distribution % Calculate the CDF value for each value t of the input array % see http://mathworld.wolfram.com/Studentst-Distribution.html for = formulas % INPUTS: t_value_array: array containing the t values for which = to % calculate the = cdf % df: degree of freedom; equals n - 1 for the = t-distribution cdf_value_array =3D 0.5 +... ((betainc(1, 0.5 * df, 0.5) / beta(0.5 * df, 0.5)) - ... (betainc((df ./ (df + t_value_array.^2)), 0.5 * df, 0.5) = /... beta(0.5 * df, 0.5))) .*... sign(t_value_array); =09 return end = %-------------------------------------------------------------------------= ---- function [t_prob_dist] =3D students_pf(df, t_arr) % calculate the probability function for students t-distribution t_prob_dist =3D (df ./ (df + t_arr.^2)).^((1 + df) / 2) /... (sqrt(df) * beta(0.5 * df, 0.5)); % % calculate and scale the cdf by hand... % cdf =3D cumsum(t_prob_dist); % discrete_t_cdf =3D (cdf - min(cdf)) / (max(cdf) - min(cdf)); % % numericaly get the t-value for the given alpha % tmp_index =3D find(discrete_t_cdf > (1 - tmp_alpha)); % t_crit =3D t(tmp_index(1)); =09 return end function in =3D isoctave () persistent inout; if isempty(inout), inout =3D exist('OCTAVE_VERSION','builtin') ~=3D 0; end; in =3D inout; return; end function [] =3D display_protocol_stack_information(pre_IP_overhead) % use http://ace-host.stuart.id.au/russell/files/tc/tc-atm/ to present = the % most likely ATM prtocol stack setup for a given overhead so the user = can % compare with his prior knowledge disp(' '); disp('According to = http://ace-host.stuart.id.au/russell/files/tc/tc-atm/'); disp(['', num2str(pre_IP_overhead), ' bytes overhead indicate']); switch pre_IP_overhead case 8 disp('Connection: IPoA, VC/Mux RFC-2684'); disp('Protocol (bytes): ATM AAL5 SAR (8) : 8'); =09 case 16 disp('Connection: IPoA, LLC/SNAP RFC-2684'); disp('Protocol (bytes): ATM LLC (3), ATM SNAP (5), ATM = AAL5 SAR (8) : 16'); =09 case 24 disp('Connection: Bridged, VC/Mux RFC-1483/2684'); disp('Protocol (bytes): Ethernet Header (14), ATM pad = (2), ATM AAL5 SAR (8) : 24'); case 28 disp('Connection: Bridged, VC/Mux+FCS RFC-1483/2684'); disp('Protocol (bytes): Ethernet Header (14), Ethernet = PAD [8] (0), Ethernet Checksum (4), ATM pad (2), ATM AAL5 SAR (8) : = 28'); =09 case 32 disp('Connection: Bridged, LLC/SNAP RFC-1483/2684'); disp('Protocol (bytes): Ethernet Header (14), ATM LLC = (3), ATM SNAP (5), ATM pad (2), ATM AAL5 SAR (8) : 32'); disp('OR'); disp('Connection: PPPoE, VC/Mux RFC-2684'); disp('Protocol (bytes): PPP (2), PPPoE (6), Ethernet = Header (14), ATM pad (2), ATM AAL5 SAR (8) : 32'); =09 case 36 disp('Connection: Bridged, LLC/SNAP+FCS RFC-1483/2684'); disp('Protocol (bytes): Ethernet Header (14), Ethernet = PAD [8] (0), Ethernet Checksum (4), ATM LLC (3), ATM SNAP (5), ATM pad = (2), ATM AAL5 SAR (8) : 36'); disp('OR'); disp('Connection: PPPoE, VC/Mux+FCS RFC-2684'); disp('Protocol (bytes): PPP (2), PPPoE (6), Ethernet = Header (14), Ethernet PAD [8] (0), Ethernet Checksum (4), ATM pad (2), = ATM AAL5 SAR (8) : 36'); =09 case 10 disp('Connection: PPPoA, VC/Mux RFC-2364'); disp('Protocol (bytes): PPP (2), ATM AAL5 SAR (8) : = 10'); =09 case 14 =20 disp('Connection: PPPoA, LLC RFC-2364'); disp('Protocol (bytes): PPP (2), ATM LLC (3), ATM = LLC-NLPID (1), ATM AAL5 SAR (8) : 14'); =09 case 40 disp('Connection: PPPoE, LLC/SNAP RFC-2684'); disp('Protocol (bytes): PPP (2), PPPoE (6), Ethernet = Header (14), ATM LLC (3), ATM SNAP (5), ATM pad (2), ATM AAL5 SAR (8) : = 40'); =09 case 44 disp('Connection: PPPoE, LLC/SNAP+FCS RFC-2684'); disp('Protocol (bytes): PPP (2), PPPoE (6), Ethernet = Header (14), Ethernet PAD [8] (0), Ethernet Checksum (4), ATM LLC (3), = ATM SNAP (5), ATM pad (2), ATM AAL5 SAR (8) : 44'); =20 =09 otherwise disp('a protocol stack this program does know = (yet)...'); end disp(' '); return; end this should give you an estimate of the required overhead.=20 Now, I only tested it with the one data set I have available so = it would be great if you could share the log file from the sweep (plus = any information about the protocol stack in use for checking of the = results=85) Oh, the script is pretty rough but since all is pretty = pedestrian it should not require to much time to understand/edit if = necessary. best sebastian --Apple-Mail=_D8AF6F97-A1E9-4ED9-A34F-AF2A1E299DF6 Content-Disposition: attachment; filename=ping_sweeper_02.sh Content-Type: application/octet-stream; name="ping_sweeper_02.sh" Content-Transfer-Encoding: 7bit #! /bin/bash # ideally one would shuffle all size*repetition to better average out transient issues # # USEAGE # 1) use "traceroute 8.8.8.8" and find the network hops just outside of your modem # typically for ATM the modem to DSLAM link will introduce an additional >10ms latency # so the hop directly after the first big latency jump should be close to the DSLAM, use # this hops IP address should be used as TARGET below # 2) optional confirm this IPs reliability by running "ping -c 100" against this IP # the min and max of the ping times should stay close to the average and the stddev # should be small as well, otherwise PINGSPERSIZE might need to be set to >=100 # TECH=CABLE # CABLE, ATM, ... just to name the log file DATESTR=`date +%Y%m%d_%H%M%S` # to allow multiple sequential records TARGET=96.34.97.78 # use traceroute something far to get the first hop out of the home net LOG=ping_sweep_${TECH}_${DATESTR}.txt PINGSPERSIZE=50 # 480 packet sizes take 13.33 hours at 100 repetitions, so stick to 50 SWEEPMINSIZE=16 # on macosx 10.7 64bit, sizes < 16 do not have a timestamp... SWEEPMAXSIZE=496 # initial size plus 10 full ATM cells (of 48 bytes payload each) n_SWEEPS=`expr ${SWEEPMAXSIZE} - ${SWEEPMINSIZE}` # do it echo "Sweeping ${n_SWEEPS} * ${PINGSPERSIZE} pings to ${TARGET} might take a while..." echo "ping -c ${PINGSPERSIZE} -g ${SWEEPMINSIZE} -G ${SWEEPMAXSIZE} ${TARGET} > ${LOG}" ping -c ${PINGSPERSIZE} -g ${SWEEPMINSIZE} -G ${SWEEPMAXSIZE} ${TARGET} > ${LOG} & # show some progress tail -f ${LOG} echo "Done... ($0) " --Apple-Mail=_D8AF6F97-A1E9-4ED9-A34F-AF2A1E299DF6 Content-Disposition: attachment; filename=tc_stab_parameter_guide_02.m Content-Type: application/octet-stream; name="tc_stab_parameter_guide_02.m" Content-Transfer-Encoding: 7bit function [ output_args ] = tc_stab_parameter_guide_02( sweep_fqn, up_Kbit, down_Kbit ) %TC_STAB_PARAMETER_GUIDE Summary of this function goes here % try to read in the result from a ping sweep run % sweep_fqn: the log file of the ping sweep against the first hop after % the DSL link % up_Kbit: the uplink rate in Kilobits per second % down_Kbit: the downlink rate in Kilobits per second % % TODO: % find whether the carrier is ATM quantized % if yes: % what is the RTT step (try to deduce the combined up and down rates from this) % try to figure out the overhead for each packet % make sure all sizes are filled (use NaN for emty ones???) %Thoughts: % maybe require to give the nominal up and down rates, to estimate the % RTT stepsize % ask about IPv4 or IPv6 (what about tunneling?) % the sweep should be taken directly connected to the modem to reduce % non-ATM routing delays dbstop if error; if ~(isoctave) timestamps.(mfilename).start = tic; else tic(); end disp(['Starting: ', mfilename]); % control options show_mean = 0; % the means are noisier than the medians show_median = 1; % the mean is the way to go show_min = 1; % the min should be the best measure, but in the ATM test sweep is too variable show_max = 0; % only useful for debugging show_sem = 0; % give some estimate of the variance show_ci = 1; % show the confidence interval of the mean ci_alpha = 0.05; % alpha for confidence interval calculation use_measure = 'median'; use_processed_results = 1; if (nargin == 0) sweep_fqn = fullfile(pwd, 'ping_sweep_ATM.txt'); % was Bridged, LLC/SNAP RFC-1483/2684 connection (overhead 32 bytes - 14 = 18) % sweep_fqn = fullfile(pwd, 'ping_sweep_CABLE_20120426_230227.txt'); % sweep_fqn = fullfile(pwd, 'ping_sweep_CABLE_20120801_001235.txt'); up_Kbit = 512; down_Kbit = 3008; end if (nargin == 1) up_Kbit = 512; down_Kbit = 3008; end if (nargin == 2) down_Kbit = 3008; end %ATM quantum.byte = 48; % ATM packets are always 53 bytes, 48 thereof payload quantum.bit = quantum.byte * 8; ATM_cell.byte = 53; ATM_cell.bit = ATM_cell.byte * 8; % CONFIGURE THE NEXT TWO VALUES line.down.Kbit = down_Kbit; line.up.Kbit = up_Kbit; % DONE line.down.bit = line.down.Kbit * 1024; line.up.bit = line.up.Kbit * 1024; % known packet size offsets in bytes offsets.IPv4 = 20; % assume no IPv4 options are used, IP 6 would be 40bytes? offsets.IPv6 = 40; % not used yet... offsets.ICMP = 8; % ICMP header offsets.ethernet = 14; % ethernet header % unknown offsets is what we need to figure out to feed tc-stab... [sweep_dir, sweep_name] = fileparts(sweep_fqn); cur_parsed_data_mat = [sweep_fqn(1:end-4), '.mat']; if (use_processed_results && ~isempty(dir(cur_parsed_data_mat))) disp(['Loading processed ping data from ', cur_parsed_data_mat]); load(cur_parsed_data_mat, 'ping'); else % read in the result from a ping sweep disp(['Processing ping data from ', sweep_fqn]); ping = parse_ping_output(sweep_fqn); save(cur_parsed_data_mat, 'ping'); end % analyze the data min_ping_size = min(ping.data(:, ping.cols.size)) - offsets.ICMP; disp(['Minimum size of ping payload used: ', num2str(min_ping_size), ' bytes.']); known_overhead = offsets.ICMP + min_ping_size + offsets.IPv4; ping.data(:, ping.cols.size) = ping.data(:, ping.cols.size) + known_overhead; % we know we used ping so add the 8 bytes already (since these will not account for overhead) size_list = unique(ping.data(:, ping.cols.size)); per_size.header = {'size', 'mean', 'median', 'min', 'max', 'std', 'n', 'sem', 'ci'}; per_size.cols = get_column_name_indices(per_size.header); per_size.data = zeros([length(size_list), length(per_size.header)]); per_size.data(:, per_size.cols.size) = size_list; for i_size = 1 : length(size_list) cur_size = size_list(i_size); cur_size_idx = find(ping.data(:, ping.cols.size) == cur_size); per_size.data(i_size, per_size.cols.mean) = mean(ping.data(cur_size_idx, ping.cols.time)); per_size.data(i_size, per_size.cols.median) = median(ping.data(cur_size_idx, ping.cols.time)); per_size.data(i_size, per_size.cols.min) = min(ping.data(cur_size_idx, ping.cols.time)); per_size.data(i_size, per_size.cols.max) = max(ping.data(cur_size_idx, ping.cols.time)); per_size.data(i_size, per_size.cols.std) = std(ping.data(cur_size_idx, ping.cols.time), 0); per_size.data(i_size, per_size.cols.n) = length(cur_size_idx); per_size.data(i_size, per_size.cols.sem) = per_size.data(i_size, per_size.cols.std) / sqrt(length(cur_size_idx)); per_size.data(i_size, per_size.cols.ci) = calc_cihw(per_size.data(i_size, per_size.cols.std), per_size.data(i_size, per_size.cols.n), ci_alpha); end figure('Name', sweep_name); hold on; legend_str = {}; if (show_mean) % means legend_str{end + 1} = 'mean'; plot(per_size.data(:, per_size.cols.size), per_size.data(:, per_size.cols.mean), 'Color', [0 1 0 ]); if (show_sem) legend_str{end + 1} = '+sem'; legend_str{end + 1} = '-sem'; plot(per_size.data(:, per_size.cols.size), per_size.data(:, per_size.cols.mean) - per_size.data(:, per_size.cols.sem), 'Color', [0 0.66 0]); plot(per_size.data(:, per_size.cols.size), per_size.data(:, per_size.cols.mean) + per_size.data(:, per_size.cols.sem), 'Color', [0 0.66 0]); end if (show_ci) legend_str{end + 1} = '+ci'; legend_str{end + 1} = '-ci'; plot(per_size.data(:, per_size.cols.size), per_size.data(:, per_size.cols.mean) - per_size.data(:, per_size.cols.ci), 'Color', [0 0.37 0]); plot(per_size.data(:, per_size.cols.size), per_size.data(:, per_size.cols.mean) + per_size.data(:, per_size.cols.ci), 'Color', [0 0.37 0]); end end if(show_median) % median +- standard error of the mean, confidence interval would be % better legend_str{end + 1} = 'median'; plot(per_size.data(:, per_size.cols.size), per_size.data(:, per_size.cols.median), 'Color', [1 0 0]); if (show_sem) legend_str{end + 1} = '+sem'; legend_str{end + 1} = '-sem'; plot(per_size.data(:, per_size.cols.size), per_size.data(:, per_size.cols.median) - per_size.data(:, per_size.cols.sem), 'Color', [0.66 0 0]); plot(per_size.data(:, per_size.cols.size), per_size.data(:, per_size.cols.median) + per_size.data(:, per_size.cols.sem), 'Color', [0.66 0 0]); end if (show_ci) legend_str{end + 1} = '+ci'; legend_str{end + 1} = '-ci'; plot(per_size.data(:, per_size.cols.size), per_size.data(:, per_size.cols.median) - per_size.data(:, per_size.cols.ci), 'Color', [0.37 0 0]); plot(per_size.data(:, per_size.cols.size), per_size.data(:, per_size.cols.median) + per_size.data(:, per_size.cols.ci), 'Color', [0.37 0 0]); end end if(show_min) % minimum, should be cleanest, but for the test data set looks quite sad... legend_str{end + 1} = 'min'; plot(per_size.data(:, per_size.cols.size), per_size.data(:, per_size.cols.min), 'Color', [0 0 1]); end if(show_max) % minimum, should be cleanest, but for the test data set looks quite sad... legend_str{end + 1} = 'max'; plot(per_size.data(:, per_size.cols.size), per_size.data(:, per_size.cols.max), 'Color', [0 0 0.66]); end title(['If this plot shows a (noisy) step function with a stepping of ', num2str(quantum.byte), ' bytes then the data carrier is quantised, make sure to use tc-stab']); xlabel('Approximate packet size [bytes]'); ylabel('ICMP round trip times (ping RTT) [ms]'); legend(legend_str); hold off; % looking at the different plot we see many side peaks %figure; %plot(per_size.data(1:end - 1, strmatch('size', per_size.header , 'exact')), diff(per_size.data(:, strmatch('median', per_size.header , 'exact'))), 'Color', [1 0 0]); % potentially clean up the data, by interpolating values with large sem % from the neighbours? % if this is ATM based the expectancy is there are RTT time quantums of 48 % bytes (the payload of an ATM package) find this %x = fminsearch(@(x)sin(x^2), x0); % estimate the RTT step size % at ADSL down 3008kbit/sec up 512kbit/sec we expect, this oes not include % processing time expected_RTT_quantum_ms = (ATM_cell.bit / line.down.bit + ATM_cell.bit / line.up.bit ) * 1000; % this estimate is rather a lower bound for fastpath , so search for best fits disp(['lower bound estimate of one ATM cell RTT is ', num2str(expected_RTT_quantum_ms), ' ms.']); % lets search from expected_RTT_quantum_ms to 1.5 * expected_RTT_quantum_ms % in steps of expected_RTT_quantum_ms / 100 % to allow for interleaved ATM set ups increase the search space up to 32 % times best fastpath RTT estimate RTT_quantum_list = (expected_RTT_quantum_ms : expected_RTT_quantum_ms / 100 : 32 * expected_RTT_quantum_ms); quantum_list = (1 : 1 : quantum.byte); differences = zeros([length(RTT_quantum_list) length(quantum_list)]); cumulative_differences = differences; all_stairs = zeros([length(RTT_quantum_list) length(quantum_list) length(per_size.data(:, per_size.cols.(use_measure)))]); for i_RTT_quant = 1 : length(RTT_quantum_list) cur_RTT_quant = RTT_quantum_list(i_RTT_quant); for i_quant = 1 : quantum.byte [differences(i_RTT_quant, i_quant), cumulative_differences(i_RTT_quant, i_quant), all_stairs(i_RTT_quant, i_quant, :)] = get_difference_between_data_and_stair( per_size.data(:, per_size.cols.size), per_size.data(:, per_size.cols.(use_measure)), ... quantum_list(i_quant), quantum.byte, 0, cur_RTT_quant ); end end % for the test DSL set the best x_offset is 21. [min_cum_diff, min_cum_diff_idx] = min(cumulative_differences(:)); [min_cum_diff_row_idx, min_cum_diff_col_idx] = ind2sub(size(cumulative_differences),min_cum_diff_idx); best_difference = differences(min_cum_diff_row_idx, min_cum_diff_col_idx); disp(['remaining ATM cell length after ICMP header is ', num2str(quantum_list(min_cum_diff_col_idx)), ' bytes.']); disp(['ICMP RTT of a single ATM packet is ', num2str(RTT_quantum_list(min_cum_diff_col_idx)), ' ms.']); figure('Name', 'Comparing ping data with'); hold on legend_str = {'ping_data', 'fitted_stair'}; plot(per_size.data(:, per_size.cols.size), per_size.data(:, per_size.cols.(use_measure)), 'Color', [1 0 0]); plot(per_size.data(:, per_size.cols.size), squeeze(all_stairs(min_cum_diff_row_idx, min_cum_diff_col_idx, :)) + best_difference, 'Color', [0 1 0]); title(['Estimated RTT per quantum: ', num2str(RTT_quantum_list(min_cum_diff_col_idx)), ' ms; ICMP header offset in quantum ', num2str(quantum_list(min_cum_diff_col_idx)), ' bytes']); xlabel('Approximate packet size [bytes]'); ylabel('ICMP round trip times (ping RTT) [ms]'); if (isoctave) legend(legend_str); else legend(legend_str, 'Interpreter', 'none'); end hold off % find cumulative_differences minimum, get the differences for that matrix % and plot the corresponding stair with the real data % next extrapolate to get the y-intercept % as first approximation use the ATM cell offset and known offsets (ICMP % IPv4 min_ping_size) to estimate the number of cells used for per paket % overhead % this assumes that no ATM related overhead is >= ATM cell size % -1 to account for matlab 1 based indices n_bytes_overhead_2nd_cell = quantum.byte - (quantum_list(min_cum_diff_col_idx) - 1); % just assume we can not fit all overhead into one cell... pre_IP_overhead = quantum.byte + (n_bytes_overhead_2nd_cell - known_overhead); % ths is the one we are after in the end disp(' '); disp(['Estimated overhead preceeding the IP header: ', num2str(pre_IP_overhead), ' bytes']); % use http://ace-host.stuart.id.au/russell/files/tc/tc-atm/ to present the % most likely ATM setup for a given overhead and present a recommendation % for the stab invocation display_protocol_stack_information(pre_IP_overhead); % now turn this into tc-stab recommendations: disp(['Add the following to both the ingress and egress root qdiscs:']); disp(' '); disp(['A) Assuming the router connects over ethernet to the DSL-modem:']); disp(['stab mtu 2048 tsize 128 overhead ', num2str(pre_IP_overhead - offsets.ethernet), ' linklayer atm']); disp(' '); disp(['B) Assuming the router connects via PPP and non-ethernet to the modem:']); disp(['stab mtu 2048 tsize 128 overhead ', num2str(pre_IP_overhead), ' linklayer atm']); disp(' '); if ~(isoctave) timestamps.(mfilename).end = toc(timestamps.(mfilename).start); disp([mfilename, ' took: ', num2str(timestamps.(mfilename).end), ' seconds.']); else toc end disp('Done...'); return end function [ ping_data ] = parse_ping_output( ping_log_fqn ) %PARSE_PING_OUTPUT read the putput of a ping run/sweep % for further processng cur_sweep_fd = fopen(ping_log_fqn, 'r'); ping_data.header = {'size', 'icmp_seq', 'ttl', 'time'}; ping_data.cols = get_column_name_indices(ping_data.header); ping_data.data = zeros([1, length(ping_data.header)]); cur_data_lines = 0; % skip the first line % PING netblock-75-79-143-1.dslextreme.com (75.79.143.1): (16 ... 1000) % data bytes header_line = fgetl(cur_sweep_fd); while ~feof(cur_sweep_fd) cur_line = fgetl(cur_sweep_fd); [first_element, remainder] = strtok(cur_line); first_element_as_number = str2num(first_element); if isempty(first_element); % skip empty lines continue; end if strmatch('---', first_element) %we reached the end of sweeps break; end % now read in the data % 30 bytes from 75.79.143.1: icmp_seq=339 ttl=63 time=14.771 ms if ~isempty(first_element_as_number) cur_data_lines = cur_data_lines + 1; % size of the ICMP package ping_data.data(cur_data_lines, ping_data.cols.size) = first_element_as_number; % now process the remainder while ~isempty(remainder) [next_item, remainder] = strtok(remainder); equality_pos = strfind(next_item, '='); % data items are name+value pairs if ~isempty(equality_pos); cur_key = next_item(1: equality_pos - 1); cur_value = str2num(next_item(equality_pos + 1: end)); switch cur_key % busybox ping and macosx ping return different key names case {'seq', 'icmp_seq'} ping_data.data(cur_data_lines, ping_data.cols.icmp_seq) = cur_value; case 'ttl' ping_data.data(cur_data_lines, ping_data.cols.ttl) = cur_value; case 'time' ping_data.data(cur_data_lines, ping_data.cols.time) = cur_value; end end end end end % clean up fclose(cur_sweep_fd); return end function [ difference , cumulative_difference, stair_y ] = get_difference_between_data_and_stair( data_x, data_y, x_size, stair_x_step_size, y_offset, stair_y_step_size ) % x_size is the flat part of the first stair, that is quantum minus the % offset debug = 0; difference = []; x_start_val = min(data_x); x_end_val = max(data_x); % construct stair stair_x = data_x; proto_stair_y = zeros([data_x(end) 1]); % make sure the x_size values do not exceed the step size... if (x_size > stair_x_step_size) if mod(x_size, stair_x_step_size) == 0 x_size = stair_x_step_size; else x_size = mod(x_size, stair_x_step_size); end end stair_y_step_idx = (x_start_val + x_size : stair_x_step_size : x_end_val); proto_stair_y(stair_y_step_idx) = stair_y_step_size; stair_y = cumsum(proto_stair_y); stair_y = stair_y(x_start_val:x_end_val) + y_offset; if (debug) figure hold on; title(['x offset used: ', num2str(x_size), ' with quantum ', num2str(stair_x_step_size)]); plot(data_x, data_y, 'Color', [0 1 0]); plot(stair_x, stair_y, 'Color', [1 0 0]); hold off; end difference = sum(abs(data_y - stair_y)) / length(data_y); cumulative_difference = sum(abs(data_y - (stair_y + difference))); return end % function [ stair ] = build_stair(x_vector, x_size, stair_x_step_size, y_offset, stair_y_step_size ) % stair = []; % % return % end function [columnnames_struct, n_fields] = get_column_name_indices(name_list) % return a structure with each field for each member if the name_list cell % array, giving the position in the name_list, then the columnnames_struct % can serve as to address the columns, so the functions assitgning values % to the columns do not have to care too much about the positions, and it % becomes easy to add fields. n_fields = length(name_list); for i_col = 1 : length(name_list) cur_name = name_list{i_col}; columnnames_struct.(cur_name) = i_col; end return end function [ci_halfwidth_vector] = calc_cihw(std_vector, n, alpha) %calc_ci : calculate the half width of the confidence interval (for 1 - alpha) % the t_value lookup depends on alpha and the samplesize n; the relevant % calculation of the degree of freedom is performed inside calc_t_val. % ci_halfwidth = t_val(alpha, n-1) * std / sqrt(n) % Each groups CI ranges from mean - ci_halfwidth to mean - ci_halfwidth, so % the calling function has to perform this calculation... % % INPUTS: % std_vector: vector containing the standard deviations of all requested % groups % n: number of samples in each group, if the groups have different % samplesizes, specify each groups samplesize in a vector % alpha: the desired maximal uncertainty/error in the range of [0, 1] % OUTPUT: % ci_halfwidth_vector: vector containing the confidence intervals half width % for each group % calc_t_val return one sided t-values, for the desired two sidedness one has % to half the alpha for the table lookup cur_alpha = alpha / 2; % if n is scalar use same n for all elements of std_vec if isscalar(n) t_ci = calc_t_val(cur_alpha, n); ci_halfwidth_vector = std_vector * t_ci / sqrt(n); % if n is a vector, prepare a matching vector of t_ci values elseif isvector(n) t_ci_vector = n; % this is probably ugly, but calc_t_val only accepts scalars. for i_pos = 1 : length(n) t_ci_vector(i_pos) = calc_t_val(cur_alpha, n(i_pos)); end ci_halfwidth_vector = std_vector .* t_ci_vector ./ sqrt(n); end return end %----------------------------------------------------------------------------- function [t_val] = calc_t_val(alpha, n) % the t value for the given alpha and n % so call with the n of the sample, not with degres of freedom % see http://mathworld.wolfram.com/Studentst-Distribution.html for formulas % return values follow Bortz, Statistik fuer Sozialwissenschaftler, Springer % 1999, table D page 775. That is it returns one sided t-values. % primary author S. Moeller % TODO: % sidedness of t-value??? % basic error checking if nargin < 2 error('alpha and n have to be specified...'); end % probabilty of error tmp_alpha = alpha ;%/ 2; if (tmp_alpha < 0) || (tmp_alpha > 1) msgbox('alpha has to be taken from [0, 1]...'); t_val = NaN; return end if tmp_alpha == 0 t_val = -Inf; return elseif tmp_alpha ==1 t_val = Inf; return end % degree of freedom df = n - 1; if df < 1 %msgbox('The n has to be >= 2 (=> df >= 1)...'); % disp('The n has to be >= 2 (=> df >= 1)...'); t_val = NaN; return end % only calculate each (alpha, df) combination once, store the results persistent t_val_array; % create the t_val_array if ~iscell(t_val_array) t_val_array = {[NaN;NaN]}; end % search for the (alpha, df) tupel, avoid calculation if already stored if iscell(t_val_array) % cell array of 2d arrays containing alpha / t_val pairs if df <= length(t_val_array) % test whether the required alpha, t_val tupel exists if ~isempty(t_val_array{df}) % search for alpha tmp_array = t_val_array{df}; alpha_index = find(tmp_array(1,:) == tmp_alpha); if any(alpha_index) t_val = tmp_array(2, alpha_index); return end end else % grow t_val_array to length of n missing_cols = df - length(t_val_array); for i_missing_cols = 1: missing_cols t_val_array{end + 1} = [NaN;NaN]; end end end % check the sign cdf_sign = 1; if (1 - tmp_alpha) == 0.5 t_val = t_cdf; elseif (1 - tmp_alpha) < 0.5 % the t-cdf is point symmetric around (0, 0.5) cdf_sign = -1; tmp_alpha = 1 - tmp_alpha; % this will be undone later end % init some variables n_iterations = 0; delta_t = 1; last_alpha = 1; higher_t = 50; lower_t = 0; % find a t-value pair around the desired alpha value while norm_students_cdf(higher_t, df) < (1 - tmp_alpha); lower_t = higher_t; higher_t = higher_t * 2; end % search the t value for the given alpha... while (n_iterations < 1000) && (abs(delta_t) >= 0.0001) n_iterations = n_iterations + 1; % get the test_t (TODO linear interpolation) % higher_alpha = norm_students_cdf(higher_t, df); % lower_alpha = norm_students_cdf(lower_t, df); test_t = lower_t + ((higher_t - lower_t) / 2); cur_alpha = norm_students_cdf(test_t, df); % just in case we hit the right t spot on... if cur_alpha == (1 - tmp_alpha) t_crit = test_t; break; % probably we have to search for the right t elseif cur_alpha < (1 - tmp_alpha) % test_t is the new lower_t lower_t = test_t; %higher_t = higher_t; % this stays as is... elseif cur_alpha > (1 - tmp_alpha) % %lower_t = lower_t; % this stays as is... higher_t = test_t; end delta_t = higher_t - lower_t; last_alpha = cur_alpha; end t_crit = test_t; % set the return value, correct for negative t values t_val = t_crit * cdf_sign; if cdf_sign < 0 tmp_alpha = 1 - tmp_alpha; end % store the alpha, n, t_val tupel in t_val_array pos = size(t_val_array{df}, 2); t_val_array{df}(1, (pos + 1)) = tmp_alpha; t_val_array{df}(2, (pos + 1)) = t_val; return end %----------------------------------------------------------------------------- function [scaled_cdf] = norm_students_cdf(t, df) % calculate the cdf of students distribution for a given degree of freedom df, % and all given values of t, then normalize the result % the extreme values depend on the values of df!!! % get min and max by calculating values for extrem t-values (e.g. -10000000, % 10000000) extreme_cdf_vals = students_cdf([-10000000, 10000000], df); tmp_cdf = students_cdf(t, df); scaled_cdf = (tmp_cdf - extreme_cdf_vals(1)) /... (extreme_cdf_vals(2) - extreme_cdf_vals(1)); return end %----------------------------------------------------------------------------- function [cdf_value_array] = students_cdf(t_value_array, df) %students_cdf: calc the cumulative density function for a t-distribution % Calculate the CDF value for each value t of the input array % see http://mathworld.wolfram.com/Studentst-Distribution.html for formulas % INPUTS: t_value_array: array containing the t values for which to % calculate the cdf % df: degree of freedom; equals n - 1 for the t-distribution cdf_value_array = 0.5 +... ((betainc(1, 0.5 * df, 0.5) / beta(0.5 * df, 0.5)) - ... (betainc((df ./ (df + t_value_array.^2)), 0.5 * df, 0.5) /... beta(0.5 * df, 0.5))) .*... sign(t_value_array); return end %----------------------------------------------------------------------------- function [t_prob_dist] = students_pf(df, t_arr) % calculate the probability function for students t-distribution t_prob_dist = (df ./ (df + t_arr.^2)).^((1 + df) / 2) /... (sqrt(df) * beta(0.5 * df, 0.5)); % % calculate and scale the cdf by hand... % cdf = cumsum(t_prob_dist); % discrete_t_cdf = (cdf - min(cdf)) / (max(cdf) - min(cdf)); % % numericaly get the t-value for the given alpha % tmp_index = find(discrete_t_cdf > (1 - tmp_alpha)); % t_crit = t(tmp_index(1)); return end function in = isoctave () persistent inout; if isempty(inout), inout = exist('OCTAVE_VERSION','builtin') ~= 0; end; in = inout; return; end function [] = display_protocol_stack_information(pre_IP_overhead) % use http://ace-host.stuart.id.au/russell/files/tc/tc-atm/ to present the % most likely ATM prtocol stack setup for a given overhead so the user can % compare with his prior knowledge disp(' '); disp('According to http://ace-host.stuart.id.au/russell/files/tc/tc-atm/'); disp(['', num2str(pre_IP_overhead), ' bytes overhead indicate']); switch pre_IP_overhead case 8 disp('Connection: IPoA, VC/Mux RFC-2684'); disp('Protocol (bytes): ATM AAL5 SAR (8) : 8'); case 16 disp('Connection: IPoA, LLC/SNAP RFC-2684'); disp('Protocol (bytes): ATM LLC (3), ATM SNAP (5), ATM AAL5 SAR (8) : 16'); case 24 disp('Connection: Bridged, VC/Mux RFC-1483/2684'); disp('Protocol (bytes): Ethernet Header (14), ATM pad (2), ATM AAL5 SAR (8) : 24'); case 28 disp('Connection: Bridged, VC/Mux+FCS RFC-1483/2684'); disp('Protocol (bytes): Ethernet Header (14), Ethernet PAD [8] (0), Ethernet Checksum (4), ATM pad (2), ATM AAL5 SAR (8) : 28'); case 32 disp('Connection: Bridged, LLC/SNAP RFC-1483/2684'); disp('Protocol (bytes): Ethernet Header (14), ATM LLC (3), ATM SNAP (5), ATM pad (2), ATM AAL5 SAR (8) : 32'); disp('OR'); disp('Connection: PPPoE, VC/Mux RFC-2684'); disp('Protocol (bytes): PPP (2), PPPoE (6), Ethernet Header (14), ATM pad (2), ATM AAL5 SAR (8) : 32'); case 36 disp('Connection: Bridged, LLC/SNAP+FCS RFC-1483/2684'); disp('Protocol (bytes): Ethernet Header (14), Ethernet PAD [8] (0), Ethernet Checksum (4), ATM LLC (3), ATM SNAP (5), ATM pad (2), ATM AAL5 SAR (8) : 36'); disp('OR'); disp('Connection: PPPoE, VC/Mux+FCS RFC-2684'); disp('Protocol (bytes): PPP (2), PPPoE (6), Ethernet Header (14), Ethernet PAD [8] (0), Ethernet Checksum (4), ATM pad (2), ATM AAL5 SAR (8) : 36'); case 10 disp('Connection: PPPoA, VC/Mux RFC-2364'); disp('Protocol (bytes): PPP (2), ATM AAL5 SAR (8) : 10'); case 14 disp('Connection: PPPoA, LLC RFC-2364'); disp('Protocol (bytes): PPP (2), ATM LLC (3), ATM LLC-NLPID (1), ATM AAL5 SAR (8) : 14'); case 40 disp('Connection: PPPoE, LLC/SNAP RFC-2684'); disp('Protocol (bytes): PPP (2), PPPoE (6), Ethernet Header (14), ATM LLC (3), ATM SNAP (5), ATM pad (2), ATM AAL5 SAR (8) : 40'); case 44 disp('Connection: PPPoE, LLC/SNAP+FCS RFC-2684'); disp('Protocol (bytes): PPP (2), PPPoE (6), Ethernet Header (14), Ethernet PAD [8] (0), Ethernet Checksum (4), ATM LLC (3), ATM SNAP (5), ATM pad (2), ATM AAL5 SAR (8) : 44'); otherwise disp('a protocol stack this program does know (yet)...'); end disp(' '); return; end --Apple-Mail=_D8AF6F97-A1E9-4ED9-A34F-AF2A1E299DF6 Content-Transfer-Encoding: quoted-printable Content-Type: text/plain; charset=iso-8859-1 >=20 > -Bill >=20 >=20 > On Tue, Jul 31, 2012 at 12:25 AM, Sebastian Moeller = wrote: >> Hi William, >>=20 >>=20 >> On Jul 30, 2012, at 11:34 AM, William Katsak wrote: >>=20 >>> Hello, >>>=20 >>> I am playing with a CeroWRT (3.3.8-6) router on my vacation in = Russia and am >>> seeing some weird behavior with simple_qos.sh that I am unsure if I >>> should attribute to a >>> bug, or to an Internet connection that is "just that bad". >>>=20 >>> Background: >>> - The router is on my wife's parents' ADSL line (according to the = modem, >>> ~3000/500). The modem is a D-Link DSL-2500U. >>>=20 >>> - Even though the link is 3000/500, and I can get speedtest.net to >>> report 2.5mbps/0.42mbps on a clean connection (direct or Cero with = no >>> QOS on), as soon as I use a host that is outside of Rostelecom = (local >>> service), it drops to 0.9/0.4mbps. This is consistent with Netalyzr >>> test: Upload 430 Kbit/sec, Download 970 Kbit/sec. This suggests that >>> even though the DSL link is higher bitrate, the ISP doesn't have the >>> outgoing bandwidth or is rate-limiting it somehow. >>>=20 >>> - I don't necessarily intend to leave the router running Cero here, >>> but I want to get a handle on the latency situation, as it makes = Skype >>> pretty messy...I am hoping to roll what I learn into a more stable >>> build of OpenWRT. >>>=20 >>> I have tried several different configurations of the modem = including: >>> 1) Default: Modem does PPPoE and hands out 192.168.1.xxx addresses. = I >>> tried just letting Cero route through that address. >>> 2) PPP-IP extension: This has the effect of the modem handling the >>> PPPoE connection and handing out the single "real" IP address over >>> DHCP. In this case Cero would see the Internet IP on ge00. >>> 3) Bridging: Allow Cero to establish the PPPoE connection and manage = it. >>>=20 >>> Right now I am in PPP-IP extension mode on the modem, and GUI QOS on >>> the router. This seems to be reliable and also keeps the latency = down, >>> although I would imagine that PPPoE on the router and the GUI QOS >>> would be fine too, but obviously I would rather use simple_qos. >>>=20 >>> The problem: >>>=20 >>> When I try simple_qos.sh, I see this: >>>=20 >>> insmod: can't insert 'cls_fw': File exists >>> insmod: can't insert 'sch_htb': File exists >>> RTNETLINK answers: No such file or directory >>> RTNETLINK answers: No such file or directory >>>=20 >>> If I run it again, the RTNETLINK errors go away...I assume this is >>> just an annoyance. >>>=20 >>> This gives me super stable ping times, etc. but a lot of websites = hang >>> loading, and the connection is unusable. If I reboot the router, the >>> connection works fine again, although the high latency comes back. >>>=20 >>> So, with all that out there, I have some questions with simple_qos: >>>=20 >>> 1) If I am using PPPoE on the router, do I need to do = IFACE=3Dpppoe-ge00 >>> or still just ge00? >>> 2) Should I set PPOE to "yes"? >>=20 >> Since your DSL connection is running PPPOE you should set PPOE = to yes in any case IF your DSL connection uses ATM as link layer (most = probable). This will just make sure that the shaper calculus the right = packet seizes to account against your link rates. But check against = http://ace-host.stuart.id.au/russell/files/tc/tc-atm/ to figure out the = prier value for overhead, as that depends on the specifics of your DSL = connection. I found that = http://www.linuxhowtos.org/manpages/8/tc-stab.htm also is quite = interesting to read to better understand the overhead parameter. But = unfortunately simple_qos does not (yet) use the generic tc-stab method = but the atm link layer adjustments specific for HTB. (Since I am on = cable right now I have no way of testing whether the tc-stab method also = works with HTB). Especially for small packets (like VoIP) if you do not = account for the the fact that ATM always sends out integer 48byte cells = and will pad if necessary, you will cause severe queueing way below = reaching the nominal link rate, as the shaper does not account for a) = the padding nor b) the 5byte ATM overhead per ATM-cell (at least I think = that is the case). >> You should do this in any case so that shaping actually has a = chance to work reliably and repeatably independent on the size = distribution of your shaped packages. >>=20 >>> 3) Is it possible that no matter what I do, the buffers at the speed >>> drop between Rostelecom and their bandwidth provider is hurting me >>> somehow? >>=20 >> If I understand correctly, yes this is going to hurt you, so = if your intended VoIP traffic leaves the Rostelecom net you might need = to specify 974/430 for the shaped rates instead of 3000/512. But that = sounds like something that is easy to test. Since your achievable uplink = (430) is still quite close to your link rate (500) I would still = recommend to look at getting link layer and overhead specified correctly = in simple_qos. >>=20 >>> 4) If 3, what to do other than yell at them? >>=20 >> As an emergency stop-gap measure shape your rates to what the = network path you are most interested in can deliver? That said I, isn't = that what codel is supposed to do automatically??? >>=20 >>>=20 >>> Overall, is anyone using Cero with a PPPoE connection with good >>> results? What kind of configuration do you have? >>=20 >> No, but I used cerowrt with a bridged ATM-based DSL connection = in the past which shares most of the issues with PPPOE over ATM. BTW = stock openwork does not account properly for ATM either so if you switch = to openwork you will need to edit some of the QOS scripts to work well = with DSL. (Last time I looked the "calculate overhead" checkbox did = something statistic I failed to fully understand). >>=20 >>>=20 >>> Sorry for the info dump, but if there is indeed a problem going on >>> with PPPoE connections, I am more than willing to be a guinea pig >>> until August 10th. I would appreciate any ideas! >>=20 >>=20 >> Oh, by the way I have some half done octave program to figure = out the actual overhead from a ping sweep, let me know if you are = interested... >>=20 >> Best Regards >> Sebastian >>=20 >>>=20 >>> Thanks! >>>=20 >>> -Bill Katsak >>> _______________________________________________ >>> Cerowrt-devel mailing list >>> Cerowrt-devel@lists.bufferbloat.net >>> https://lists.bufferbloat.net/listinfo/cerowrt-devel >>=20 --Apple-Mail=_D8AF6F97-A1E9-4ED9-A34F-AF2A1E299DF6--