Add support for parsing Fourier Analysis

[kgoodrick-rl]

Hi Nuno,

First of all, thank you for creating and sharing this library! I have found it very useful for my work. On a recent project I needed to parse the Fourier analysis output in LTSpice, so I added in the functionality, and this pull request contains those additions.

The Fourier results are stored in the log dataset dictionary with the following structure:

All Fourier analysis is stored in the dataset under the "fourier" key, and the base level is a dictionary keyed by the fundamental frequency as a float. Each value of the base dict is another dictionary keyed by the number of periods of the analysis. The values for this dictionary are another dictionary keyed by the name of the waveform. This is the final dictionary, and it contains a list of the results for each step in the simulation.

The results are a dictionary with three items:

1. The DC component as a float and keyed by the string 'dc'
2. The total harmonic distortion as a float and keyed by the string 'thd'
3. The harmonics stored as a pandas' dataframe and keyed by the string 'harmonics'
   a. The dataframe has columns ['Harmonic Number', 'Frequency [Hz]', 'Fourier Component', 'Normalized Component', 'Phase [degree]', 'Normalized Phase [degree]']
   b. Each row is a different harmonic

Here is an example of how to get to the result dictionary:

log.dataset["fourier"][frequency][n_periods][waveform_name][step_index]

This is functional, though I'm not crazy about how deep the nesting goes (especially if you want to access something in the dataframe), so I'm open to suggestions on alternates.

To get this working with some of the recent changes I had to change the function that detects the encoding as my version of LTSpice outputs files with cp1252 encoding. During this function I read the entire file to ensure there are no errors at any point. This is not a problem with the log files I am dealing with, but it could be if there were a large log file with many steps.

I also changed the way the log file is read such that it is in a context manager.

Let me know what you think and if there's anything in the code you'd like changed or added.

Thanks,

Kyle

[nunobrum]

Hi Kyle,

I'll only have time to make a proper revision later this week.
As a first impression, I tend to agree with you that the return format needs a bit more thought.
I've never used the Fourier Analysis myself, but, it seems to me that if you'd return the pandas DataFrame itself it would serve the user just as well, if not better. Let me know what you think.

Best Regards,
Nuno

[kgoodrick-rl]

Hi Nuno,

Sounds good. Unfortunately, a data frame alone is insufficient as the DC value and THD don't fit in the table. Newer versions of pandas have an attrs attribute that can be used to store metadata, but according to the documentation it is still experimental, so there is some risk involved with that approach. It is also less obvious where to find that information with that approach.

Here's an example of the output of a Fourier analysis to give you a better idea of what it looks like. There can be any number of rows in the table.

N-Period=3
Fourier components of V(x)
DC component:2.0013

Harmonic	Frequency	 Fourier 	Normalized	 Phase  	Normalized
 Number 	  [Hz]   	Component	 Component	[degree]	Phase [deg]
    1   	1.356e+07	2.613e+00	1.000e+00	  -73.71°	    0.00°
    2   	2.712e+07	2.613e+00	1.000e+00	  -73.71°	    0.00°
Total Harmonic Distortion: 0.000000%(16.625442%)

-Kyle

[nunobrum]

Hello Kyle,

Finally I have a moment to revise your code. Again sorry for the late reply.
Regarding the detect_encoding, it's all good, nothing to say and I'm ready to accept your changes.

On the LTSteps.py I have a few remarks.

On line 360, I think you're missing a backslash before the point if you want that the regular expression matches the point.
Also, I ask myself whether you don't want to also capture and store the value between parenthesis. In your example it seems to have more meaning than the zero.
Total Harmonic Distortion: 0.000000%(16.625442%).
I would suggest the following :

                            regex = re.search(r": (.+)%\((.+)%\)", line)
                            thd = float(regex.group(0))
                            thd1 = float(regex.group(1))

I like the pandas library, but, since it is adding a dependency, I was wondering whether it couldn't be supported by a dedicated class that could provide a better interface to the user. Able to also compute other things like relative differences of fourier components in regards to a reference signal. Simulation accuracy and data compression have an impact on the fourier analysis, so, maybe a relative distortion could make sense. It would also allow to get harmonics on a single trace, or the frequency component on a list of traces. Or to calculate SNR and SINAD when providing a noise level (different simulation).
Let me know what you think.

The library dependencies was so far an important point for me. You might note that even the usage of numpy as been so far an optional. Maybe this no longer makes sense, if considering that almost anyone making signal analysis, should have numpy and pandas installed. What are your thoughts one that?

On line 391 and a few more time below, you are using this construct var_name: type_name = value.
This is not supported in older versions of python (ex: 3.5), I would suggest not to use it, as it would mean that PyLTSpice would no longer be supported in older versions.

Thanks again for your initiative to share your code. I can already say welcome your code, irrespective of you agreeing with my suggestion regarding having dedicated class.

Best Regards,
Nuno

[kgoodrick-rl]

Hi Nuno,

Thank you for the feedback!

1. I agree that we should save the second THD. I wasn't able to find documentation on what it is so I didn't bother saving it, but someone may want it, so it'd be good to have, and I will update to capture this information. If you know what the difference between the two is let me know and I'll give them meaningful names otherwise, I'll call the one outside of parenthesis thd and the other thd_alt.
2. I think a dedicated class to store this information makes sense. I wanted results quickly for my project and already had pandas, so it wasn't an issue, but you're right that it doesn't make sense to include it for just this one case. Especially since it isn't able to hold all the information elegantly. The ability to get at the data from different directions will also be nice, the big, nested dictionary is not the most user friendly. This might take a little more time, but I will try to work on it in the next few weeks.
3. I agree that probably everyone using this library will probably have numpy and pandas installed anyways, but unless it's really needed, I agree that it would be nice to not require it. In this case I think a good (or maybe even better) solution without pandas is possible, so probably good to leave it out.
4. I really like the type hinting in Python, I think it makes code much more readable and useable. How far back do you want to support? It looks like all of the officially supported python versions do fully support type hinting. It wouldn't be a problem to take these out though if you want to support older versions of Python.

I'm going to be pretty busy the next few weeks, but I will try to find time to work on the class for storing harmonic information. Let me know if you have any more thoughts!

-Kyle

[nunobrum]

Hi Kyle,

I'm happy that we are converging. I was afraid that you wouldn't.
Defining a good interface and document it can be quite time consuming and perhaps some iterations, so, I'm happy that you are willing to go a bit beyond the simple "code that does the it".

I agree with you on point 4. I really like type hinting, as it brings some formalism to a language that otherwise has none. Type hinting has very often prevented some bugs, so, I'm definitely in favor of it.
I see in this PEP that introduced it it was first made available in the python 3.6, and, python.org the 3.5 version is no longer active.
So, I reckon it makes sense to adopt it.

On the second value of THD, I've spent some time searching for an explanation. The best I could was an answer directly from Mike Engelhardt which says:
The 2nd figure an alternative method of computing THD. It nominally includes
all harmonics but numerically difficulting in identifying the fundamental
limits it actually use. It's just bonus information for the guy that asked for
the feature. You can ignore it. But is the number is higher the previous
number, it is an indication that you didn't include enough harmonics in
your THD calculation.

Generally, most people just ignore the second value.
If you're up to it, you can follow the discussion in these two threads.
https://www.diyaudio.com/community/threads/ltspice-four-thd-results-xxxx-xxxx.309887/
https://www.diyaudio.com/community/threads/whats-ltspice-doing-there.263590/#post4097547

Before moving forward and start implementing, I think it makes sense to take a moment to define an interface. Having it defined, will help define an implementation.

Here are some ideas to seed a possible interface. Note that none of this is binding in any way.
Just launching some ideas so to start the discussion.

So, having said this disclaimer, here it goes:
The THD analysis data could be accessed by

log = LTSpiceLogReader(filename)
THD_analysis = log.THD

In this object I could imagine it returning all the data you are now storing.
At this moment, I haven't yet the feeling to what makes sense in terms of priority on the first access.
By looking how the data is structured on the log file and following the same approach as the rest of the code, I'd suggest:
a) Net Name,
b) Frequency
c) Step [if applicable.]

So, it would translate to the following:

n = THD_analysis['net_name'].n_harmonics  # Returns the number of harmonics
dc = THD_analysis['net_name'].dc_Component # Returns the DC component
thd = THD_analysis['net_name'].thd   # Calculated LTSpice calculated distortion
freq_info= THD_analysis['net_name'].harmonics(1) # Returns an object with data relative to the fundamental. 
freq_info= THD_analysis['net_name'].fundamental  # same as previous
freq3_info= THD_analysis['net_name'].harmonics(3) # Returns the 3rd harmonic frequency information.
freq3000 = freq3.frequency  # Supposing the fundamental is 1000, this should return 3000
freq3_info_alt = THD_analysis['net_name'].frequency(3000) # alternative way of accessing frequency information
fcomp3 = freq3_info.fourier   # Fourier component accessed by harmonic, if steps exist, returns a list
norm3 = freq3_info.norm  # normalized fourier component
phase3 = freq3_info.phase  # phase information
phnorm3 = freq3_info.phasenorm  # normalized phase
freq3step5comp = freq3_info.step(step_number)

As said... any of this is binding... just brainstorming.

Have a great Sunday ahead,

Nuno

[nunobrum]

Hi again,

I just remembered that more than one .FOUR command can be issued for the same net and with different frequencies.
In the case, my approach will not work.
This definitely needs some refinement.

Cheers,
Nuno

[nunobrum]

Accepting the pull request.
I'll just move the pandas dependency to be local and only required when Fourier analysis is required.
Later this dependency will be dropped in favor of a dedicated class.

[nunobrum]

Hi @kgoodrick-rl ,

I hope you had a good summer and matching holidays.
I would like to move forward creating a class for representing the THD Analysis data. Did you had any chance to make any progress on that?
Best Regards,
Nuno

[kgoodrick-rl]

Hi @nunobrum,

Sorry for the delay in getting back to you! I just moved, so things have been a little hectic the last several weeks. Hope things are going well for you too! I have finally found time to come back to this (thanks for the reminder) and I think I have a solution that I think should work quite well. It's currently implemented with dataclasses, but if we want python <3.6 compatibility I think we could implement the dataclass features we need on our own. The first pass at this can be seen here: https://github.com/kgoodrick-rl/PyLTSpice/tree/fourier_interface

As a high level overview, the data is stored in the log's fourier attribute. The attribute is a list with elements corresponding to the steps in the simulation. Each element is a dictionary (I subclassed UserDict in case we want to add some helper methods / properties) keyed by a FourierKey object with FourierAnalysis objects as the values. FourierKey objects contain the waveform string, the fundamental frequency, the number of harmonics, and the number of periods in the analysis. FourierAnalysis objects contain the DC component, both THDs, and a HarmonicsList of FourierHarmonic. HarmonicsList is a normal list, but indexing starts at 1. This means that if you want the information for the third harmonic you do not have to index it with 2. The FourierHarmonic objects contain the harmonic number, the frequency, the Fourier component, the phase, and the normalized versions of magnitude and phase.

As an example here is how you would use it to get the fourier component of the third harmonic of the V(a) waveform:

key = FourierKey('V(a)', frequency, n_harmonics, n_periods) # These are all parameters specified in the Fourier spice command
third_harmonic = log.fourier[0][key][3].fourier_component # 0 is the step

The reason I am using a key object instead of just the net name is that someone could theoretically run multiple Fourier analysis on the same waveform, but with different fundamentals, periods, or number of harmonics. We could maybe have a helper function though that figures everything out if there is only one Fourier analysis for a given waveform?

The way it's implemented now the step indexing happens at the beginning, but I think it would be possible to have a redefined __getitem__ function that would let the user do it either way?

Let me know what you think.

Best,

Kyle

[nunobrum]

Hi Kyle,

sorry for the late reply. This was a very busy week and only today have some time to check your solution.
I have to say that it looks really good. So good that I think that I don't have to test it. ;-)

As for the avoiding the dataclasses for 3.6 (Note 1) compatibility, I think we could eventually avoid dataclasses, but, it would make the code much more verbose. I confess to be a bit divided regarding this. I have no real perception of what is the python version adopting across the community, but, I'm going to settle with the gut feeling that most people have versions higher then 3.7.

Back to implementation. I like the idea of keying all data in a single structure, which avoids the long sequences of [ ][ ] of the previous implementation. I was worried about the safety of letting the compiler hash all that information, but, as far as I understood of the dataclasses, it is.

I was a bit surprised by the implementation of appending the empty dictionary to a list and then filling that dictionary using the last position of the list.
self.fourier[-1][fourier_key] = fourier_analysis
I reckon you are doing this to simplify the implementation for processing new step information, right?

Finally, the __getitem__ can be redefined to have one way or the other. This was actually what was done on the main LTSpiceLogReader class, where the user can pass the step number or directly a measure name.
The only thing that is slicing was forbidden exactly for that reason. But, now that I think about it, this could be also solved.
Anyway... the class FourierResults can always be expanded to include other functionality. For the time being the current structure looks quite elegant and can hold the data safely. So, I'm OK to go with it.

(Note 1) I think it's actually 3.7 (https://peps.python.org/pep-0557/)

All the best and wishing you a good rest of the weekend,

Nuno

[kgoodrick-rl]

Hi Nuno,

I'm glad you like the solution! You are correct that dataclasses are python 3.7 and above, my bad. My thought is that anyone writing new code will probably be on python 3.9 or 3.10, for anyone with existing code on older versions they should be ok by just continuing to use the older version of PyLTSpice.

That is my understanding of dataclasses as well. I have used them as keys many times and have not had any problems doing so. You are also right that implementing them ourselves will be much more verbose. That was my initial plan, but it ended up being considerably more effort, so I reverted to dataclasses.

I'm not particularly set on this way of indexing the results while parsing them, it just seemed the most obvious to me when I was writing it. We could potentially use self.step_count to index the list, but I'm not sure how this would work if step_set is not empty when __init__ is called. I'd welcome your thoughts on this part.

Should we modify the indexing and think about a fancier _getitem__ first or should I open a pull request now?

Thanks,

Kyle

[nunobrum]

Hi Kyle,

I'm OK with you just opening a pull request.
As I see it, the current implementation already does a decent job of parsing the Fourier analysis.
Of course things can always be improved, but, as long as we are not moving backwards by breaking things that were working or introducing fragilities or poor implementations, then there is nothing to lose.
The thing that worries me the most is introducing bugs that would lead to false simulation results.

So far I've been testing all the code I write manually, but I'm thinking about using automating the code tests with unittests.
This is huge cleanup pending for a long time already. The plan for the upcoming months is to automate tests, so, can you provide a testbench for your code? No urgency. I have other things to deal with in the short term.
Besides my daytime job, I saw this post in a repo I follow asweigart/pyperclip#229. I need to migrate the setup.py into the new protocol and I have requests from people using NGSpice for a code that can read both LTSpice and NGSpice.

Cheers,
Nuno