function [csgz,qxx,dert0,dlj]=jjpc(bb,ll,pp)
pn=diag(pp);
nbb=bb'*pn*bb;
ww=bb'*pn*ll;
qxx=inv(nbb);
csgz=qxx*ww;
vv=bb*csgz-ll;
vtpv=ll'*pn*ll-ww'*csgz;
[rr,nn]=size(bb);
dert0=sqrt(vtpv/rr);
dlj=der0*dert0*bb*qxx*bb';