Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

feat: parse parameters and units with DDParsers #390

Open
wants to merge 6 commits into
base: main
Choose a base branch
from
Open

Conversation

c-dilks
Copy link
Member

@c-dilks c-dilks commented Dec 11, 2022

Briefly, what does this PR introduce?

Use DD4hep/DDParsers to parse configuration parameters. This permits the user to include units that are automatically converted to an internal system of units. This internal unit system can be any preferred system we agree on and does not have to be DD4hep's system; e.g., it could be configured to parse parameters to EDM4hep's standard units.

Source code changes:

  • parser.{h,cpp}: new class Parser:
    • use DD4hep/DDParsers Evaluator to do the parsing and math
    • hard-coded default system of units, to enforce a standard
      • straightforward to change, if we ever want to
      • other systems are included for reference, but commented out
      • under the hood, the system is set by setSystemOfUnits.cpp
    • this Parser class can be pulled out of EICrecon, should we decide to use it elsewhere; we could also replace DD4hep Evaluator usage with that of CLHEP, if we don't want DD4hep dependence here
  • eicrecon_cli.cpp: use Parser to parse eicrecon -c user parameters
    • this is currently the only usage of Parser
    • if we move to a configuration file for configuration parameters, Parser could be used to read its numbers
  • unit test param_parser_test.cpp
    • consider adding a CI job to run this and other unit tests
    • consider relocating to src/tests/

Some usage examples (defaulting to DD4hep's unit system, 1=cm=GeV):

2022-12-07-011830_925x488_scrot

It can also handle lists, parsing each element (this example uses EDM4hep's unit system, 1=ns=GeV):

2023-04-24-202609_1055x39_scrot

Example unit systems:

EDM4hep units:

1 = millimeter = GeV = nanosecond = radian

Geant4 units:

1 = millimeter = MeV = nanosecond = nanoampere (via e+ charge) = kelvin = mole = candela = radian = steradian

This was assumed by reco_flags.py prior to #365

DD4hep units:

Our current eic-shell build of DD4hep uses the following system:

1 = centimeter = GeV = second = nanoampere = kelvin = mole = candela = radian = steradian

NOTE: see #365 for more discussion

The parser was fully compatible with the legacy reco_flags.py, and its results have been cross checked. Consistency with previous versions of EICrecon was verified by comparing run_eicrecon_reco_flags.py output, with the -c option used for eicrecon. The only differences are minor, typically dropping trailing zeros (e.g., "4.0" becomes "4").

What kind of change does this PR introduce?

  • Bug fix (issue #__)
  • New feature (issue #__)
  • Documentation update
  • Other: __

Please check if this PR fulfills the following:

  • Tests for the changes have been added
  • Documentation has been added / updated
  • Changes have been communicated to collaborators

Does this PR introduce breaking changes? What changes might users need to make to their code?

No

Does this PR change default behavior?

Yes, the specification of parameters and units is enhanced

@c-dilks c-dilks changed the title feat!: parse parameters and units with DD4hep (CLHEP) Evaluator feat!: parse parameters and units with DDParsers Dec 11, 2022
@c-dilks
Copy link
Member Author

c-dilks commented Dec 15, 2022

@sly2j the EDM4hep units are

1 = millimeter = GeV = nanosecond = radian

correct? What about the others (see PR description above).. in particular nanoampere is very slightly different between Geant4 and DD4hep.

Also, I don't think we necessarily need this, but do we have a file similar to DD4hepUnits.h for EDM4hep/eic? It could be convenient someday...

@veprbl
Copy link
Member

veprbl commented Apr 24, 2023

So it is an issue that EDM4hep units don't match up with DD4hep. I was hoping we could adopt a single unit convention across the codebase. We are, obviously, tied to the the units in the EDM, so it would make sense to use its convention by default. While we do, currently, rely on DD4hepUnits.h, that is not a dependency. The only values in DD4hep units that we will see, will probably come from geometry (is TGeo unit-agnostic?).

@c-dilks
Copy link
Member Author

c-dilks commented Apr 24, 2023

We are, obviously, tied to the the units in the EDM

I think that's the general preference, but not sure if that's official yet.

While we do, currently, rely on DD4hepUnits.h, that is not a dependency.

Yeah, this #include can be dropped if we choose EDM, but we still need DD4hep headers in #include <Evaluator/*>, unless we switch to CLHEP's Evaluator, from which DD4hep's is based.

@veprbl
Copy link
Member

veprbl commented Apr 24, 2023

We are, obviously, tied to the the units in the EDM

I think that's the general preference, but not sure if that's official yet.

Well, we'd probably need at least alternative constants for our system of units.

While we do, currently, rely on DD4hepUnits.h, that is not a dependency.

Yeah, this #include can be dropped if we choose EDM, but we still need DD4hep headers in #include <Evaluator/*>, unless we switch to CLHEP's Evaluator, from which DD4hep's is based.

That is not a problem if we can leverage this:
https://github.com/AIDASoft/DD4hep/blob/master/DDParsers/src/Evaluator/setSystemOfUnits.cpp
Note that there is a wrapper that may not expose the full evaluator functionality. In my experience, dropping it is fine.

@c-dilks c-dilks marked this pull request as ready for review April 24, 2023 23:33
@c-dilks c-dilks changed the title feat!: parse parameters and units with DDParsers feat: parse parameters and units with DDParsers Apr 24, 2023
@c-dilks
Copy link
Member Author

c-dilks commented Apr 25, 2023

That is not a problem if we can leverage this:
https://github.com/AIDASoft/DD4hep/blob/master/DDParsers/src/Evaluator/setSystemOfUnits.cpp
Note that there is a wrapper that may not expose the full evaluator functionality. In my experience, dropping it is fine.

This PR already uses its wrapper, effectively calling

setSystemOfUnits(/*EDM4hep system*/);

You make a good point, we can leverage this for anything related to units, even unit conversion:

REQUIRE(std::abs( 1.0 / m_parser->units("GeV") - 1.0 ) < DIFF_THRESHOLD ); // [GeV] -> [GeV]
REQUIRE(std::abs( 1.0 / m_parser->units("keV") - 1.0e+6 ) < DIFF_THRESHOLD ); // [GeV] -> [keV]
REQUIRE(std::abs( 1.0 / m_parser->units("s") - 1.0e-9 ) < DIFF_THRESHOLD ); // [ns] -> [s]
REQUIRE(std::abs( 1.0 / m_parser->units("nm") - 1.0e+6 ) < DIFF_THRESHOLD ); // [mm] -> [nm]
REQUIRE(std::abs( 1240/1e9/1e6 / m_parser->units("eV*nm") - 1240.0 ) < DIFF_THRESHOLD ); // [GeV*mm] -> [eV*nm]

It could use some sugar, but this demonstrates that

  1. the parser can read a user-specified parameter string, in whatever units, evaluate any math operations, and automatically convert it to a number in the internal unit system (e.g., EDM4hep)
  2. given some number, assuming it's in the internal unit system, convert it to some preferred units (well, just SI prefix conversions in the most common case, but the main point of all this is to globally set which SI prefix corresponds to "1" for each unit type)

@veprbl
Copy link
Member

veprbl commented May 3, 2023

I suggest we use something like eic/EDM4eic#35

@c-dilks
Copy link
Member Author

c-dilks commented May 3, 2023

I suggest we use something like eic/EDM4eic#35

Good idea! We can certainly then just use these here to set the parser's scales:

// EDM4hep default system
// 1 = millimeter = GeV = nanosecond = radian
m_eval = std::make_unique<dd4hep::tools::Evaluator::Object>(
1.e+3, // millimeter (per meter)
1./1.60217733e-25/1.e+3, // GeV (per kilogram); cf. Geant4 unit system below
1.e+9, // nanosecond (per second)
1./1.60217733e-10, // nanoampere, via e+ charge (per ampere); cf. Geant4 unit system below
1.0, // kelvin
1.0, // mole
1.0, // candela
1.0 // radian
);

something like

    m_eval = std::make_unique<dd4hep::tools::Evaluator::Object>(
        edm4eic::unit::m
        // ...
   )

@c-dilks
Copy link
Member Author

c-dilks commented May 3, 2023

We could then consider adding some converter functionality, e.g. something to do

number_in_edm_units = number_from_DD4hep * edm4hep::unit::mm / dd4hep::mm;

Here's some pseudocode to sketch the general idea

  • this is beyond the scope of this PR but uses two dd4hep::tools::Evaluator::Object instances to do the work
  • users can still multiply and divide units to do the conversions themselves, but IMO it's easy to make mistakes doing this (for example, edm4hep::unit::mm / dd4hep::m)
  • dictionary lookups in Evaluator::Object are likely slower than just using the constants as is
  • maybe overkill, but then we get everything in both systems
/* dd4hep::Evaluator objects: contains dictionary for a system of units
 * - this is what DD4hep uses to parse the compact file strings
 * - they can be initialized with any custom unit system: you just need to say what is "1"
 *   for each unit type (see DD4hep source code for examples)
 * - we may be able to use something simpler; there is a version from CLHEP too
 */
enum unit_systems { kDD, kEDM, nSys };
std::shared_ptr<dd4hep::tools::Evaluator::Object> eval[nSys];
eval[kDD]  = std::make_shared<dd4hep::tools::Evaluator::Object>(/* DD4hep units specification */);
eval[kEDM] = std::make_shared<dd4hep::tools::Evaluator::Object>(/* EDM4hep units specification */);

/* general converter of `val` from one unit system `from` to another `to`
 * - calls `Evaluator::Object::eval` to do the math and the dictionary lookup of `units_expr`
 * - `units_expr` can be one unit, or an expression of units, e.g., "eV*nm/s"
 * - `from` and `to` are `enum unit_systems` numbers, to specify the unit systems
 *   - if set to `-1`, use "natural unit system", where lack of an SI prefix really means "1"
 *     (in contrast, EDM4hep's length unit is "mm", so "m" does not mean "1" in that system)
 */
double ConvertUnits(
    double      val,
    std::string units_expr,
    int         from,
    int         to
    )
{
  auto conv = [&eval, &units_expr, &nSys] (int u) {
    if(u==-1)
      return 1.0;
    else if(u>=0 && u<nSys) {
      auto result = eval[u]->eval(units_expr);
      if(result.status() == dd4hep::tools::Evaluator::OK)
        return result.result();
    }
    m_log->error("error");
    return val;
  };
  return val * conv(to) / conv(from);
}

/* specialized converters (sugar)
 * - uses ConvertUnits to do the conversion DD4hep <-> EDM4hep
 * - example usage:
 *    
 *     auto sensor_size = ConvertUnitsDD2EDM(
 *       DD4hep_service->get_number("sensor_size"), // pseudocode; has DD4hep units
 *       "mm" // the units
 *       );
 *     // result: `sensor_size` will be in EDM4hep units
 *     // NOTE: "mm" could be any length unit, e.g. "m", or "nm"; this is because
 *     // "edm4hep::mm" / dd4hep::mm == "edm4hep::m" / dd4hep::m == etc.
 *
 */
double ConvertUnitsDD2EDM(double val, std::string units_expr) {
  return ConvertUnits(val, units_expr, kDD, kEDM);
}
double ConvertUnitsEDM2DD(double val, std::string units_expr) {
  return ConvertUnits(val, units_expr, kEDM, kDD);
}

/* resolve units (sugar)
 * - by default, any number in an algorithm can be assumed to have units in some unit system
 * - printing this number out may not be clear, since the result will be that
 *   number together with the power of 10 representing the SI prefix(es), or
 *   other factors such as PI
 * - use these "resolvers" to get the raw number in "natural" units,  where
 *   lack of an SI prefix really means "1"
 * - example usage
 *
 *     // get a sensor size, in EDM4hep units (same as above example)
 *     auto sensor_size = ConvertUnitsDD2EDM(
 *       DD4hep_service->get_number("sensor_size"), // pseudocode; has DD4hep units
 *       "mm" // the units
 *       );
 *     // print it in [nm]; we know the unit system is EDM, so use `ResolveUnitsEDM`
 *     m_log->trace("sensor size is {} nm", ResolveUnitsEDM(sensor_size, "nm"));
 *
 */
double ResolveUnitsDD(double val, std::string units_expr) {
  return ConvertUnits(val, units_expr, kDD, -1);
}
double ResolveUnitsEDM(double val, std::string units_expr) {
  return ConvertUnits(val, units_expr, kEDM, -1);
}

@veprbl
Copy link
Member

veprbl commented May 3, 2023

I was hoping we could somehow dodge the need for conversion. Perhaps, wrap it into the DD4hep service API? Or, perhaps, override the global TGeoSystemOfUnits plus fix (potential) DD4hep bugs?

@c-dilks
Copy link
Member Author

c-dilks commented May 4, 2023

Avoiding (automating) conversions would be ideal. What about Acts::UnitConstants? Both ACTS and DD4hep have geometry services, so that would be a good place for automatic conversion. Do we have any other unit systems?

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

2 participants