The canonical way to approach this problem is I-Q demodulation.
In RF terms, we're using a synchronous detector/demodulator. A real AC signal has in-phase (I) and out-of-phase (Q) components, or sin and cos if you prefer, and if we simply multiplied a signal (point by point) by a sine wave, we'd lose half that information. So we multiply by cosine as well, and this preserves all information in the signal.
$$
\begin{eqnarray}
I & = \int_0^{2\pi} x(t) \sin t \, dt \\
Q & = \int_0^{2\pi} x(t) \cos t \, dt
\end{eqnarray}
$$
Written here as integral over a single cycle, giving the correlations for that cycle. In practice, this might be done by taking a weighted average of the product over some time window (that is to say, a lowpass filter), as in an analog receiver; or by changing \$d\$s to \$\Delta\$s and integral to summation to give the DSP case (in other words, what you'd implement on an MCU). (A receiver might then take the magnitude (\$\sqrt{I^2 + Q^2}\$) for AM, or the angle change per cycle for FM, for example; but we're doing a much simpler kind of detection here.)
Of course, \$x(t)\$ is one signal; you have two (current and voltage), so you can either assume rather than measure one (when the circuit permits such -- using an op-amp to force voltage, or read current with minimal voltage drop, for example), or repeat the analysis for both. (Your resources seem limited, suggesting the double measurement would be undesirable?)
In case both are measured, you need to divide one by the other to get the normalized phase and magnitude:
$$\hat{Z} = R + jX = \frac{I_V + jQ_V}{I_I + jQ_I}$$
(Denoting complex impedance Z with a hat, and voltage and current readings with respective subscripts.)
Measuring just one, you're still doing the same thing (you can't get impedance from units of voltage and current other than by taking their ratio!), but using an assumed value for one of the complex numbers: this will be a calibration parameter, which ideally should be set on initial power-up testing. It may depend on frequency, in which case a lookup table would be helpful (calibrating over a range of frequencies).
And that's it. The \$X\$ can be converted to equivalent L or C, and then we have the R-C or R-L series equivalent circuit at the measured frequency, or we can apply the suitable transformations to get the parallel equivalents.
Which, as long as it's fresh in my mind -- I recently wrote that for a calculator on my website:
Vector Impedance by Resistor Divider | Calculators | Seven Transistor Labs, LLC
Opening dev. console shows the code quoted below (plain old JS):
function vec_Update() {
var R = parseFloat(document.getElementById("vec_R").value);
var V1 = parseFloat(document.getElementById("vec_V1").value);
var V2 = parseFloat(document.getElementById("vec_V2").value);
var th = parseFloat(document.getElementById("vec_th").value);
var F = parseFloat(document.getElementById("vec_F").value);
th = th * Math.PI / 180;
F = 2 * Math.PI * F;
var H = V2 / V1;
var HR = H * Math.cos(th), HI = H * Math.sin(th);
var ZR = (R * (HR - H * H)) / (1 - 2 * HR + H * H);
var ZX = (R * HI) / (1 - 2 * HR + H * H);
document.getElementById("vec_ZR").innerHTML = ZR.toPrecision(9);
document.getElementById("vec_ZX").innerHTML = ZX.toPrecision(9);
document.getElementById("vec_ZMag").innerHTML = Math.sqrt(ZX * ZX + ZR * ZR).toPrecision(9);
document.getElementById("vec_ZArg").innerHTML = (Math.atan2(ZX, ZR) * 180 / Math.PI).toPrecision(5);
document.getElementById("vec_ESR").innerHTML = ZR.toPrecision(9);
document.getElementById("vec_ESL").innerHTML = (ZX / F).toPrecision(9);
document.getElementById("vec_ESC").innerHTML = (-1 / (F * ZX)).toPrecision(9);
var ZPR = (ZR * ZR + ZX * ZX) / ZR;
var ZPX = (ZR * ZR + ZX * ZX) / ZX;
document.getElementById("vec_EPR").innerHTML = ZPR.toPrecision(9);
document.getElementById("vec_EPL").innerHTML = (ZPX / F).toPrecision(9);
document.getElementById("vec_EPC").innerHTML = (-1 / (F * ZPX)).toPrecision(9);
}
This takes in a magnitude and angle as measured by an oscilloscope, but it's immediately decomposed into I and Q parts (HR and HX). A more accurate measurement would be had by taking in the raw waveforms and doing the correlations, but, well, that's rather hard to do from a browser(!).
This also uses a slightly more complicated expression, since the voltage divider has to be accounted for: it's not measuring current and voltage, it's measuring a ratio of voltages. The reference resistance (assumption: its phase is very close to 0°) is a calibration parameter here. In your case, the ZR (\$R\$) and ZX (\$X\$) are given as above.
Incidentally, this method (phase angle by zero-crossing) kind of has the same problem as your proposal; they're quite different in method of course, but in the sense that, they aren't measuring exactly the intended aspect (fundamental sin/cos components), but aspects only related to it (real/apparent power there, zero-crossing here). As such, both are sensitive to components other than the fundamental -- harmonics and noise.
And, if you are in such a situation, where assumptions can be taken, perhaps fairly strong ones -- maybe taking ratios of magnitudes (as yours), or zero-crossing points, or assuming a scalar is wholly real or imaginary as the case may be -- as long as you're comfortable with those assumptions, that can be effective too. Common example: capacitor ESR testers generally measure magnitude voltage, assuming it to be a resistive component; but ESL, and C itself if small, also participate; ESL especially, since the waveforms are often sharp (e.g. square wave). This makes the readings diverge rather dramatically for long lead lengths, or at small values, but such a device is still useful enough as both of those can be controlled. Using one on in-circuit electrolytics of modest or greater value, meets both of those assumptions.
