From 93724b4a3afe887d9f2e6afdb7a5cf896515c620 Mon Sep 17 00:00:00 2001 From: TheFGFSEagle Date: Mon, 27 Jun 2022 22:12:55 +0200 Subject: [PATCH] Added script for VSPAERO .history file -> JSBSim FDM table conversion, added sinusoidal method to utils.interpolator --- .gitignore | 0 LICENSE | 0 README.md | 0 __init__.py | 0 aircraft/vsphist2jsbtable.py | 123 +++++++++++++++++++++++++++++++++++ utils/__init__.py | 0 utils/interpolator.py | 57 +++++++++------- 7 files changed, 156 insertions(+), 24 deletions(-) mode change 100644 => 100755 .gitignore mode change 100644 => 100755 LICENSE mode change 100644 => 100755 README.md mode change 100644 => 100755 __init__.py create mode 100755 aircraft/vsphist2jsbtable.py mode change 100644 => 100755 utils/__init__.py diff --git a/.gitignore b/.gitignore old mode 100644 new mode 100755 diff --git a/LICENSE b/LICENSE old mode 100644 new mode 100755 diff --git a/README.md b/README.md old mode 100644 new mode 100755 diff --git a/__init__.py b/__init__.py old mode 100644 new mode 100755 diff --git a/aircraft/vsphist2jsbtable.py b/aircraft/vsphist2jsbtable.py new file mode 100755 index 0000000..0aed015 --- /dev/null +++ b/aircraft/vsphist2jsbtable.py @@ -0,0 +1,123 @@ +#!/usr/bin/env python +#-*- coding:utf-8 -*- + +import argparse +import os, sys +from fgtools.utils.interpolator import Interpolator + +class Case: + def __init__(self, textcase): + textcase.pop(0) # get rid of the Name Value Unit legend line + self.Sref = float(textcase.pop(0).split()[1]) + self.Bref = float(textcase.pop(0).split()[1]) + self.Cref = float(textcase.pop(0).split()[1]) + self.Xcg = float(textcase.pop(0).split()[1]) + self.Ycg = float(textcase.pop(0).split()[1]) + self.Zcg = float(textcase.pop(0).split()[1]) + self.Mach = float(textcase.pop(0).split()[1]) + self.AoA = float(textcase.pop(0).split()[1]) + self.Beta = float(textcase.pop(0).split()[1]) + self.Rho = float(textcase.pop(0).split()[1]) + self.Vinf = float(textcase.pop(0).split()[1]) + self.RollRate = float(textcase.pop(0).split()[1]) + self.PitchRate = float(textcase.pop(0).split()[1]) + self.YawRate = float(textcase.pop(0).split()[1]) + textcase.pop(0) # get rid of the Solver case line + textcase.pop(0) # get rid of another legend line + lastiterindex = textcase.index("Skin Friction Drag Break Out:") - 1 + self.CL, self.CDo, self.CDi, self.CDtot, self.CS, self.LD, self.E, self.CFx, self.CFy, self.CFz, self.CMx, self.CMy, self.CMz, self.CDtrefftz, self.TQS = map(float, textcase.pop(lastiterindex).split()[4:]) # we don't need Mach, AoA, Beta again + + def __str__(self): + return "Case(" + ", ".join(str(k) + " = " + str(v) for k, v in vars(self).items()) + ")" + + def __repr__(self): + return self.__str__() + +def get_cases(path): + path = os.path.abspath(path) + + if not path.endswith(".history"): + print("The specified file", path, "is not a VSPAERO .history file - exiting") + sys.exit(1) + + with open(path, "r") as histf: + textcases = [] + lines = histf.readlines()[1:] + textcase = [] + cases = {} + + for line in lines: + if line.startswith("****"): + textcases.append(list(filter(None, textcase))) + textcase = [] + continue + textcase.append(line.strip()) + textcases.append(list(filter(None, textcase))) + + for textcase in textcases: + case = Case(textcase) + if not case.AoA in cases: + cases[case.AoA] = {} + cases[case.AoA][case.Beta] = case + + return cases + +def print_table(cases, coeff, indent, precision): + coeffs = {} + for AoA in cases: + for Beta in cases[AoA]: + if not hasattr(cases[AoA][Beta], coeff): + print("The coefficient", coeff, "does not exist in the specified .history file - exiting") + sys.exit(1) + + if not AoA in coeffs: + coeffs[AoA] = {} + coeffs[AoA][Beta] = getattr(cases[AoA][Beta], coeff) + + print("") + print(indent + 'aero/alpha-deg') + print(indent + 'aero/beta-deg') + #print(indent + 'velocities/mach') + print(indent + "") + print(indent + indent + indent + indent.join(map(str, coeffs[list(coeffs.keys())[0]].keys()))) + for AoA in coeffs: + print(indent + indent + str(AoA), end="") + for Beta in coeffs[AoA]: + print(indent + ("%." + str(precision) + "f") % coeffs[AoA][Beta], end="") + print() + print(indent + "") + print("
") + + +if __name__ == "__main__": + argp = argparse.ArgumentParser(description="vsphist2jsbtable.py - Takes a VSPAERO .history file and a coefficient as input and outputs a JSBSim interpolation table for that coefficient from the VSPAERO cases") + + argp.add_argument( + "-c", "--coeff", + help="Name of the coefficient to produce a table for", + required=True, + ) + + argp.add_argument( + "--indentation", + help="The argument of this option will be used as one level of indentation, defaults to a tab", + default="\t" + ) + + argp.add_argument( + "-p", "--precision", + help="How many decimal places the numbers in the table should have, defaults to 6", + type=int, + default=6 + ) + + argp.add_argument( + "input_file", + help="VSPAERO .history file", + ) + + args = argp.parse_args() + + cases = get_cases(args.input_file) + print_table(cases, args.coeff, args.indentation, args.precision) + diff --git a/utils/__init__.py b/utils/__init__.py old mode 100644 new mode 100755 diff --git a/utils/interpolator.py b/utils/interpolator.py index e9304cb..dd8fa39 100755 --- a/utils/interpolator.py +++ b/utils/interpolator.py @@ -1,13 +1,13 @@ #!/usr/bin/env python #-*- coding:utf-8 -*- +import math + class Interpolator: def __init__(self): self._indexes = [] self._values = [] self._sorted = False - - self.methods = {"linear": self._interpolate_linear} def add_value(self, index, value): if type(index) not in (int, float) or type(value) not in (int, float): @@ -26,7 +26,7 @@ class Interpolator: self.add_value(i, v) def interpolate(self, index, extrapolate=True, method="linear", sort=True): - if not method in self.methods: + if not hasattr(self, "_interpolate_" + method): raise NotImplementedError(f"Interpolator.interpolate: interpolation method '{method}' not yet supported") if len(self._indexes) < 2: @@ -38,7 +38,22 @@ class Interpolator: self._values.sort() self._sorted = True - return self.methods[method](index, extrapolate) + if index in self._indexes: + return self._values[self._indexes.index(index)] + + if self._indexes[0] < index < self._indexes[-1]: + return getattr(self, "_interpolate_" + method)(index, extrapolate) + else: + if not extrapolate: + if index < self._indexes[0]: + return self._values[0] + else: + return self._values[-1] + else: + if index < self._indexes[0]: + return self._values[1] + (index - self._indexes[1]) / (self._indexes[0] - self._indexes[1]) * (self._values[0] - self._values[1]) + else: + return self._values[-2] + (index - self._indexes[-2]) / (self._indexes[-1] - self._indexes[-2]) * (self._values[-1] - self._values[-2]) def _find_neighbours(self, index): lower = upper = 0 @@ -53,30 +68,24 @@ class Interpolator: return lower, upper def _interpolate_linear(self, index, extrapolate=True): - if index in self._indexes: - return self._values[self._indexes.index(index)] - - if self._indexes[0] < index < self._indexes[-1]: - lower, upper = self._find_neighbours(index) - return self._values[lower] + (self._values[upper] - self._values[lower]) * (index - self._indexes[lower]) / (self._indexes[upper] - self._indexes[lower]) - else: - if not extrapolate: - if index < self._indexes[0]: - return self._values[0] - else: - return self._values[-1] - else: - if index < self._indexes[0]: - return self._values[1] + (index - self._indexes[1]) / (self._indexes[0] - self._indexes[1]) * (self._values[0] - self._values[1]) - else: - return self._values[-2] + (index - self._indexes[-2]) / (self._indexes[-1] - self._indexes[-2]) * (self._values[-1] - self._values[-2]) + lower, upper = self._find_neighbours(index) + return self._values[lower] + (self._values[upper] - self._values[lower]) * (index - self._indexes[lower]) / (self._indexes[upper] - self._indexes[lower]) + + def _interpolate_sinusoidal(self, index, extrapolate=True): + lower, upper = self._find_neighbours(index) + return self._values[lower] + (self._values[upper] - self._values[lower]) * math.sin(math.radians((index - self._indexes[lower]) / (self._indexes[upper] - self._indexes[lower]) * 90)) # run test if run directly if __name__ == "__main__": - print("Test results") + print("Test results:") + print() i = Interpolator() i.add_values((0, 10, 20), (0, 20, 30)) - for test_val in (-5, 0, 1, 2, 3.5, 5.55555, 9, 10, 15, 100): - print(test_val, i.interpolate(test_val)) + test_vals = (-5, 0, 1, 2, 3.5, 5.55555, 10, 15, 35, 100) + for method in ("linear", "sinusoidal"): + print(method.capitalize() + ":") + for test_val in test_vals: + print(test_val, "=>", i.interpolate(test_val, method=method)) + print()