thliebig / openEMS-Project

openEMS is a free and open electromagnetic field solver using the FDTD method.
418 stars 78 forks source link

How to create a custom excitation? #153

Open oberstet opened 11 months ago

oberstet commented 11 months ago

During a discussion here https://github.com/thliebig/openEMS-Project/discussions/152#discussioncomment-7730788, the question of whether it is possible, and if so, how to create and use a custom excitation in openEMS.

If that is possible or even supported, ideally from Python, but if not I'd be good with a way using C++ as well (fwiw, if so I might find time to add Python wrappers then as well).

Here is the concrete question:

In openEMS, and also in the Python API of openEMS, there are a few built-in excitation functions, like e.g. SetGaussExcite.

Is it possible to use openEMS with a custom excitation function, such as for example a continuous white noise spectrum shaped by a single Gaussian envelope over a multi-second time scale and with the source of randomness to generate the white noise coming from a PRNG instance seeded from a user supplied value?

thliebig commented 11 months ago

The matlab interface can define custom excitation's. I'm not sure if the C++ code has an API-interface exposed for it at this point (I need to check). If it does, the python wrapper for it would be a one-liner (as for all the other excitation) definitions.

But I'm always puzzled why anybody really wants this. I've been doing RF simulation with FDTD for decades now and I have never ever needed anything else than a Gaussian pulse. Any other pulse is just less efficient and will only result in longer simulation time for absolutely zero additional information. I'm sure that there are some oddball corner case were another excitation pulse/signal could be helpful/useful, but currently I would have no idea what that could be. But I'm always open to learn something new ;)

Anyways, since openEMS has this ability, we should make it available in the C++ and python API. Please feel free to add the missing pieces and open a Pull Request in the openEMS repository ;)

thliebig commented 11 months ago

I just checked, the C++ does not expose the interface, but should be easy to add... See it as a small task to get into the openEMS code? ;)

oberstet commented 11 months ago

alright;) let's see

(cpy311_7) oberstet@intel-nuci7:~/scm/3rdparty/openEMS$ find . -name "*.cpp" -exec grep -Hi "gauss" {} \;
./openems.cpp:      exc->SetupGaussianPulse(f0, fc);
./openems.cpp:void openEMS::SetGaussExcite(double f0, double fc)
./openems.cpp:  m_Exc->SetupGaussianPulse(f0, fc);
./FDTD/excitation.cpp:void Excitation::SetupGaussianPulse(double f0, double fc)
./FDTD/excitation.cpp:      CalcGaussianPulsExcitation(m_f0,m_fc,maxTS);
./FDTD/excitation.cpp:void Excitation::CalcGaussianPulsExcitation(double f0, double fc, int nTS)
./FDTD/excitation.cpp:      cerr << "Operator::CalcGaussianPulsExcitation: Requested excitation pusle would be " << Length << " timesteps or " << Length * dT << " s long. Cutting to max number of timesteps!" << endl;

https://github.com/thliebig/openEMS/blob/5f36e7f3a2367123f00999491a069aed50c6f244/FDTD/excitation.h https://github.com/thliebig/openEMS/blob/5f36e7f3a2367123f00999491a069aed50c6f244/FDTD/excitation.cpp

actually, it doesn't seem too hard. I mean, adding it right to above class - there doesn't seem to be means for custom user, after compile addition of excitation functions


I am wondering, there is a Excitation::DumpCurrentExcite ... but there is no Excitation::LoadCurrentExciteFromDump ... if there would be, plus a Python wrapper, that might actually be a lot more flexible ... than a Excitation::CalcGaussianShapedWhiteNoiseExcitation ... what do you think?

oberstet commented 11 months ago

Any other pulse is just less efficient and will only result in longer simulation time for absolutely zero additional information.

I don't think this is true. For me, and again I have zero RF experience, I just understand "numbers" to a degree, a gaussian shaped white noise excitation is the natural solution to the following.

I want an excitation that excites at all frequencies all the time (while the excitation is running), and I want the excitation start with 0 energy, reach a max half time, and disappear to 0 energy again, with the energy envelope being a smooth function of time (eg differentiable), and I want to be able to freely control the overall time (excitation energy 0 to energy 0)

the natural solution is a gaussian envelope (rgd energy by time, spectral power density) over white noise. there are some details, since the gauss curve only touches 0 for infinity .. and we don't have that much time ... so lets cutoff (step function) at -60dB ... or cutoff at 2x var(gauss) ... or something.

oberstet commented 11 months ago

fwiw, I we had a

void Excitation::CalcGaussianShapedWhiteNoiseExcitation(... whatever param ..., int nTS)

the power spectrum would not be real white noise (infinite frequencies, straight flat spectrum), but it would be limited by

I'm not sure if the latter is limited by lumped port geometry, material or what.

FDTD_FLOAT seems likely a float32 .. and hence the swing is limited by maximum dynamic range a 32 bit float can catch ... sure, it's probably good enough for sub GHz;) I am just saying. you might not see THz or PHz or whatever frequencies in the spectrum. so it is not perfect white noise .. but probably RF_GOOD_ENOUGH =)

oberstet commented 11 months ago

btw, there are indeed

which allows to provide a string defining a custom excitation function based on

http://warp.povusers.org/FunctionParser/fparser.html


I haven't seen it in the Python API though, and it isn't going to cut it for the white noise thing, and also not for loading previously stored excitation dumps ...

oberstet commented 11 months ago

here is a sketch https://github.com/thliebig/openEMS/pull/129 ... only the skeleton so it could be discussed whether such a thing is wanted by the project anyways, if PRs are welcome, and if so, exactly how it should integrate etc

I am unsure about: code style and docs though as well, otherwise I would have added function docs with statements of what it is supposed to do before starting any implementation work ...

oberstet commented 11 months ago

rgd the energy cutoffs ... the cutoffs to the gaussian envelope, the natural way of talking about "how far out" or "how unlikely" isn't "dB", but multiples of standard deviations, see eg https://introcs.cs.princeton.edu/python/appendix_gaussian/

one can transform those into dBs .. if I haven't goofed up sth, it should boil down to:

>>> import math
>>> 
>>> for k in range(9):
...     percentage = 1 - (1 - math.erf(k / math.sqrt(2)))
...     db = int(round(math.log10(1 - percentage) * 10))
...     print("k={} standard deviations: {} or {} dB".format(k, str(round(percentage * 100, 20)) + " %", db))
... 
k=0 standard deviations: 0.0 % or 0 dB
k=1 standard deviations: 68.26894921370858 % or -5 dB
k=2 standard deviations: 95.44997361036415 % or -13 dB
k=3 standard deviations: 99.73002039367398 % or -26 dB
k=4 standard deviations: 99.99366575163337 % or -42 dB
k=5 standard deviations: 99.99994266968562 % or -62 dB
k=6 standard deviations: 99.99999980268247 % or -87 dB
k=7 standard deviations: 99.99999999974403 % or -116 dB
k=8 standard deviations: 99.99999999999987 % or -149 dB

-149 dB ... is really going too far .. maybe it makes sense (practically) to take "5 standard deviations" which correspond to -62 dB? what do you think?

thliebig commented 11 months ago

Well I was talking about a Pull request for the missing Custom Excite API. I don't think that I will add a noise excitation any time soon. I do not see any practical use case for FDTD??? And in any case, if someone wants to play with noise excitation and different variants of it, this is what the custom excite might be good for.

oberstet commented 11 months ago

Well I was talking about a Pull request for the missing Custom Excite API.

ok, understood. fwiw, if I find time, I'm going to expose it as well, added the bits in this PR (sorry, seems like the discussion is here, but the code is there, I've got confused;) anyways: https://github.com/thliebig/openEMS/pull/129

I don't think that I will add a noise excitation any time soon. I do not see any practical use case for FDTD???

right, ok. it might be a foolish idea of mine;) no worries, should I continue this, I'm just gonna run from my fork, np.

And in any case, if someone wants to play with noise excitation and different variants of it, this is what the custom excite might be good for.

are you interested in having a new excite able to read voltage/current from files previously generated ones using DumpVoltageExcite/DumpCurrentExcite? are you interested in having those functions exposed in Python as well?

unfortunately, I've started working on all of these in one PR .. but I could certainly sort that out and keep only what you want - should I find time and actually manage to make sth work that is;=

thliebig commented 7 months ago

Sorry this thread got lost on me. Having a function to read an excitation signal would be useful I guess as then everyone can easily create whatever excitation function is desired. But that function would need to be able to interpolate as the file not necessarily should have the same sampling steps?

Please remember there is no separate voltage/current excitation signal... It needs to be the same signal just sampled at different times... (off by half a timestep)

oberstet commented 7 months ago

great, thanks for noticing and commenting!

couple of comments back from my side:

But that function would need to be able to interpolate as the file not necessarily should have the same sampling steps?

I would not mix that into the responsibilities of the load function, but put the burden completely on the user of such a function.

IOWs: the user must supply excitation data with

Ideally, the data should be read/written in an open standard format for general IQ data.

A quick search turned up a format called SigMF

https://pysdr.org/content/iq_files.html

The SigMF file itself is just the pure raw IQ in numpy (array of complex128 or what). obviously that has been recognized by the PySDR guys

SigMF and Annotating IQ Files Since the IQ file itself doesn’t have any metadata associated with it, it’s common to have a 2nd file, containing information about the signal, with the same filename but a .txt or other file extension. This should at a minimum include the sample rate used to collect the signal, and the frequency to which the SDR was tuned.

Personally, I think what would be even better, and an extension:

both the raw IQ data file, and the metadata file, and e.g. a human-readable text file (a README.txt) should then be ZIP'ed into one archive file

the 3 files contained in such a IQ ZIP should then have standard/predefined names, like

The most simple (and minimal) way to use the SigMF standard to describe a binary IQ file you have created is to rename the .iq file to .sigmf-data and create a new file with the same name but .sigmf-meta extension, and make sure the datatype field in the meta file matches the binary format of your data file.

for example, one qpsk_in_noise.zip with contents:

  1. signal.sigmf-data
  2. signal.sigmf-meta
  3. README.txt
oberstet commented 7 months ago

here is more info about SigMF:

fwiw, and here is what an AI says about SigMF - which sounds mostly cool/good/correct to me. Just that I very much prefer the native binary signal formatting that PySDR uses: numpy complex128! rather than "XML or JSON format". at least for the IQ data.

one major reason is: I can mmap a native binary with numpy complex128 and start running native code on that without any additional memory usage or pressure from parsing and recreating a run-time repr. of the IQ data

for the metadata, I'd definitely prefer JSON rather than XML, and in fact maybe YAML instead of JSON, but that's just because writing is much nicer in YAML ... anyways, details.


Q:

Can you describe (in technical detail) the open file format "SigMF" for storing/transporting general "IQ sample" data?

A:

SigMF (Signal Metadata Format) is an open file format designed for storing and transporting general "IQ sample" data, commonly used in software-defined radio (SDR) applications and other signal processing tasks. The format aims to provide a standardized way to describe and exchange signal datasets along with their associated metadata. Here's a technical description of the SigMF format:

  1. File Structure:

    • SigMF files are typically stored in either XML or JSON format, although XML is more commonly used.
    • The file consists of two main sections: metadata and data.
    • Metadata includes information about the signal such as sample rate, center frequency, data type, and any additional annotations.
    • Data section contains the actual IQ samples.
  2. Metadata:

    • SigMF metadata provides descriptive information about the signal dataset.
    • Common metadata fields include:
      • global metadata: Information applicable to the entire dataset.
      • captures metadata: Information specific to individual captures within the dataset.
    • Metadata fields may include:
      • sample_rate: Sampling rate of the IQ samples.
      • center_frequency: Center frequency of the signal.
      • data_type: Data type of the IQ samples (e.g., complex64, complex128).
      • description: Brief description of the signal.
      • annotations: Additional notes or annotations related to the signal.
      • author: Name of the person or entity who captured the data.
      • license: Licensing information for the data.
      • timestamp: Timestamp indicating when the data was captured.
    • Metadata fields can be extended to accommodate specific requirements of different applications.
  3. Data:

    • The data section contains the IQ samples.
    • IQ samples are typically stored as interleaved pairs of real (I) and imaginary (Q) components.
    • Samples can be stored in binary format or as raw text data.
    • The data format may vary based on the specified data_type in the metadata.
    • SigMF supports both contiguous and non-contiguous data storage.
  4. Usage:

    • SigMF files can be used for a variety of purposes, including storing captured RF data, exchanging datasets between different SDR platforms, and archiving signal recordings.
    • They provide a standardized way to describe and share signal datasets, making it easier for researchers and developers to work with diverse datasets across different platforms.
  5. Extensibility:

    • SigMF is designed to be extensible, allowing users to define custom metadata fields and adapt the format to specific use cases or applications.
    • Extensions can be added to accommodate additional metadata requirements or specialized signal processing tasks.

In summary, SigMF is a flexible and extensible file format for storing and transporting IQ sample data, providing a standardized way to describe signal datasets and their associated metadata. Its adoption helps promote interoperability and facilitates collaboration in the field of software-defined radio and signal processing.

oberstet commented 7 months ago

one final thought: what about multi-channel IQ data?

say I have 4 excitation ports defined in my openEMS model (and no, I haven't done that yet even though I am very much interested .. for phased arrays and such ..)

can I load IQ excitation data for all 4 ports from 1 file and run the excitation on all ports synchronized / in parallel?


Q:

How would I use SigML for multi-channel IQ data e.g. in the context of Phased Arrays?

A:

In the context of phased arrays or multi-channel IQ data, SigMF can be extended to accommodate the specific requirements of handling multiple channels. Here's how you could use SigMF for multi-channel IQ data:

  1. Metadata Extension:

    • Extend the SigMF metadata to include information specific to each channel.
    • Add metadata fields to describe characteristics of each channel, such as antenna position, gain, phase, polarization, or any other relevant parameters.
    • Define a structure within the metadata to organize channel-specific information, such as an array of channel descriptions.
  2. Channelized Data:

    • Store IQ samples for each channel separately.
    • Organize the data section of the SigMF file to include separate sections or blocks for each channel.
    • In the data section, ensure that IQ samples for each channel are stored in the appropriate order and format, maintaining synchronization between channels.
  3. Channel Enumeration:

    • Include metadata fields to specify the number of channels and their ordering.
    • Use a consistent numbering scheme or identifier for each channel to ensure clarity and consistency.
  4. Timestamp Synchronization:

    • If the IQ samples across channels are captured synchronously, include timestamp information in the metadata to indicate the timing relationship between samples from different channels.
    • Ensure that timestamps are accurately recorded and aligned across channels to facilitate coherent processing.
  5. Array Geometry:

    • If applicable, include metadata fields to describe the physical geometry of the phased array, such as element spacing, array layout, and beamforming configuration.
    • Provide annotations or additional metadata to describe the beamforming weights or algorithms used to process the multi-channel data.
  6. Data Representation:

    • Depending on the application, consider the appropriate representation of IQ samples for each channel (e.g., interleaved or separate real and imaginary components).
    • Choose a suitable data type and format for storing multi-channel IQ samples, considering factors such as precision, storage efficiency, and compatibility with downstream processing algorithms.

By extending SigMF in this way, you can effectively represent multi-channel IQ data in the context of phased arrays or similar applications. This enables standardized storage, exchange, and analysis of multi-channel RF datasets, facilitating interoperability and collaboration among researchers and developers working in the field.


Right. I hadn't thought of storing phased array geometry metadata ..

thliebig commented 7 months ago

Wow this is going very wild here? IQ data? The FDTD method is a time-domain method. It has no concept of frequencies or complex numbers what so ever.

The only thing that would make any kind of sense to read in as en excitation signal (from a file) would be a simple ASCII file with two columns, one is time, the second is the time-signal, e.g. a Gaussian pulse or whatever the user desires and that the FDTD method can handle, e.g. the signal must be band-limitied...

And I repeat, only a single, simple, scalar time-signal is possible, nothing multi-channel, nothing complex or anything of that sort.

I think that there is still a major misconception about what FDTD is, what it can do and what it is supposed to do....

thliebig commented 7 months ago

I think it would be better we do this in a github Discussion not as an issue? And then we sort out what the feature request actually is. If it is something related to simulation setup or more to a higher level data post-processing, e.g. in a circuit simulation where IQ and multi-channel and all that is another story...

oberstet commented 7 months ago

IQ data? The FDTD method is a time-domain method. It has no concept of frequencies or complex numbers what so ever.

Not sure I am following, isn't

It needs to be the same signal just sampled at different times... (off by half a timestep)

hinting at 90° (the "half timestep"), and isn't

there is no separate voltage/current excitation signal.

just "a detail" .. as in: it's 2 vectors of real numbers?

but in any case, you are right, probably better in a discussion format .. lemme know / pls link here once/if you start one ..

oberstet commented 7 months ago

btw: as repo owner you should also be able to convert this whole issue in aggregate into a discussion preserving the contents above .. using one button click (top right usually .. "convert to discussion" or sth)