function [pr,CO,sk]=calibrate(x,pr,T,c)
% x and pr are 1 by n vectors, T=T(x) is p by n and c is p by 1.
% tilted distribution pr so that E(T)=c
% outputs covariance matrix CO of T  and skewness
n=length(x) ;  
m=T*pr'   ;  p=length(m) ;% row vector of mean values of T
D=ones(1,n) ;
L=zeros(p,1) ; 
while sum(abs(m-c)>.01)
 m=T*pr'   ;
 V=T-repmat(m,1,n);
 CO=V.*repmat(pr,p,1)*V';    %   the covariance matrix p times p
 L=-inv(CO)*(m-c);
 pr=pr.*exp(L'*T);
 pr=pr/sum(pr);
 
end
m=T*pr'   ;
 V=T-repmat(m,1,n);
 CO=V.*repmat(pr,p,1)*V';
 sk=((V.^3)*pr')/CO^1.5;