function Dxs = diffsg(x,m1,m2,diffs,col)
% This is a function to calculate derivatives of the map 
% f: R^m1 -> R^m2
% defined by restricting an (m1)st level complex of delta
% circuits to an (m2)nd level one. We approximate derivatives with
% respect to perturbing each edge, and for each edge in the restriction,
% using difference quotients with a variety of sizes of
% perturbation.
%
% Support files: diffquotSeq.m
%
% INPUT: 
% x - 3^m1 by 3 matrix containing the value at which to calculate
% derivatives
% m2 - the level to restrict to
% diffs - vector of delta values (presumably approaching zero) to use
% in difference quotients
%
% OUTPUT:
% Dxs - matrix holding length(diffs)*3^(m1+1) row vectors of
% derivatives obtained from the diffquotSeq program

% initialize Dxs
Dxs = [] ;

% for each delta to use
for i = 1:length(diffs)  

  % create matrix to hold derivatives
  block = zeros(3^(m1+1),3^(m2+1)) ;  

  % for each edge in starting complex compute
  for j = 1:3^(m1+1)
    block(j,:) = diffquotSeq(x,col,j,diffs(i),m1,m2) ;
  end  
  
  % store derivatives in the matrix Dxs
  Dxs = [Dxs; block] ;
end