extract_tau1 is not used
and broken as the following code does not make much sense:
for (int i = nrMin; i < nrMax; i++ ) {
for (int h = 1; h < ni-1; h++ ) {
real deltaphi = radMed(i)*(phi(h+1)-phi(h));
real moy = 0.;
for (int j = 0; j < ns; j++ ) {
moy += dens(i,h,j)*kappa(i,h,j);
}
moy /= ns;
tau(i) += moy*deltaphi;
if(tau(i) <= 1){ // this override prev h values
ztau1(i) = radMed(i)*cos(phi(h));
phitau1(i) = phi(h);
}
}
}