Skip to content
Draft
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
37 changes: 34 additions & 3 deletions Data/emeterDecode.py
Original file line number Diff line number Diff line change
@@ -1,5 +1,6 @@
from nptdms import TdmsFile
import polars as pl
import numpy as np
import matplotlib.pyplot as plt
from Data.FSLib.IntegralsAndDerivatives import *
from Data.FSLib.AnalysisFunctions import *
Expand All @@ -9,10 +10,10 @@
autoxDaniel12File = "FS-3/compEmeterData/autoxDaniel12.tdms"
accel1 = "FS-3/compEmeterData/216_univ-of-calif---santa-cruz-_250620-203307_ ACCEL-EV.tdms"
accel2 = "FS-3/compEmeterData/216_univ-of-calif---santa-cruz-_250620-205609_ ACCEL-EV.tdms"
endur1 = "FS-3/compEmeterData/216_univ-of-calif---santa-cruz-_250621-154731_ ENDUR-EV.tdms"
endur2 = "FS-3/compEmeterData/216_univ-of-calif---santa-cruz-_250621-160530_ ENDUR-EV.tdms"
endur1 = "../fs-data/FS-3/compEmeterData/216_univ-of-calif---santa-cruz-_250621-154731_ ENDUR-EV.tdms"
endur2 = "../fs-data/FS-3/compEmeterData/216_univ-of-calif---santa-cruz-_250621-160530_ ENDUR-EV.tdms"

dfLaptimes = pl.read_csv("FS-3/compLapTimes.csv")
dfLaptimes = pl.read_csv("../fs-data/FS-3/compLapTimes.csv")
firstHalf = dfLaptimes.filter(pl.col("Lap") < 12)["Time"].sum()
secondHalf = dfLaptimes.filter(pl.col("Lap") > 11)["Time"].sum()

Expand Down Expand Up @@ -93,7 +94,37 @@ def fileTodf(path):
pl.Series(arr).cast(pl.Int64).alias("Lap")
)

l = []

for i in np.unique(dfendur1["Lap"]):
# plt.plot(dfendur1.filter(pl.col("Lap") == i)[I])
l.append(dfendur1.filter(pl.col("Lap") == i)[I])
plt.show()

shortest = min([len(x) for x in l])
l2 = [x[:shortest].alias(f"Current_Lap_{i}") for i, x in enumerate(l)]
df2 = pl.DataFrame(l2)
plt.plot(df2.mean_horizontal())
plt.show()

from scipy.fft import fft, ifft

f = fft(df2.mean_horizontal().to_numpy())
# freq = np.fft.fftfreq(len(df2.mean_horizontal()), d=0.01)

plt.plot(np.append(np.log(f[-len(f)//2:]), np.log(f[:len(f)//2])))
plt.show()

fFiltered = np.where(np.log(f) > 10, f, 0)
invF = ifft(fFiltered)
plt.plot(invF)
plt.plot(df2.mean_horizontal())
plt.show()

dfTableCurrOut = pl.DataFrame({"Current": df2.mean_horizontal().gather_every(100), "Time": np.arange(0, df2.height/100)})
dfTableCurrOut.write_csv("endur1Curr.csv")

df2.mean_horizontal()

dfendur2 = fileTodf(endur2).filter(pl.col(t) > endur2_StartTime).filter(pl.col(t) < endur2_EndTime)

Expand Down
39 changes: 4 additions & 35 deletions Data/temp.py
Original file line number Diff line number Diff line change
Expand Up @@ -3,8 +3,8 @@
import cantools.database as db

from Data.DataDecoding_N_CorrectionScripts.dataDecodingFunctions import *
from Data.AnalysisFunctions import *
from Data.integralsAndDerivatives import *
from Data.FSLib.AnalysisFunctions import *
from Data.FSLib.IntegralsAndDerivatives import *
from scipy.interpolate import CubicSpline

dbcPath = "../fs-3/CANbus.dbc"
Expand Down Expand Up @@ -70,38 +70,7 @@
etcRTDButton = "ETC_STATUS_RTD_BUTTON"
etcBrakeVoltage = "ETC_STATUS_BRAKE_SENSE_VOLTAGE"

df = read("C:/Projects/FormulaSlug/fs-data/FS-3/10112025/firstDriveMCError30.parquet")
df = df.with_columns(
df["timestamp"].alias("Time")
)
df = read("C:/Projects/FormulaSlug/fs-data/FS-3/10082025/fixed_wheels_nathaniel_inv_test_w_fault.parquet")

df = read("C:/Projects/FormulaSlug/fs-data/FS-3/10112025/firstDriveMCError30-filled-null.parquet")
df = df.with_columns(
simpleTimeCol(df)
)

fig = plt.figure()
ax = fig.add_subplot(111)

ax.plot(df[t], df[frT], label=frT, c="blue")
ax.plot(df[t], df[flT], label=flT, c="red")
ax.plot(df[t], df[brT], label=brT, c="orange")
ax.plot(df[t], df[blT], label=blT, c="cyan")
ax.set_title("Suspension Travel during First Drive with MC Fault")
ax.set_xlabel("Time")
ax.set_ylabel("Suspension Travel (mm)")
ax.legend()
plt.show()


dfNullless = df.drop_nulls(subset=[frT, flT, brT, blT])

cs = CubicSpline(dfNullless[t], dfNullless[frT])

fig = plt.figure()
ax = fig.add_subplot(111)

ax.scatter(dfNullless[t], cs(dfNullless[t]), label=frT, s=0.5)
ax.scatter(dfNullless[t], in_place_derive(cs(dfNullless[t])), label=f"Derived {frT}", s=0.5)
ax.legend()
plt.plot(df[busV])
plt.show()
82 changes: 82 additions & 0 deletions Docs/SimulationTodo.md
Original file line number Diff line number Diff line change
@@ -0,0 +1,82 @@
# Simulation Todo

1. Better drag model that takes into account aeropackage (FS-3)
1. Individual wheel models
1. Suspension
1. Travel (x)
1. Velocity (v)
1. How will this react under acceleration in any direction (steering causes lateral acceleration) (throttle/brakes causes longitudinal acceleration)
1. Wheel
1. Brake temp (more complex soon)
1. Wheel temp (more complex soon)
1. Wheel rpm/speed
1. Differential + Drivetrain losses
1. Energy loss due to chain, tripods, axle, hub/upright rubbing
1. Model rolling resistance better
1. Model limited slip differential
1. Log losses so we have an idea of energy loss in the drivetrain
1. Motor Efficiency + Heating
1. Function of temp and current draw
1. Keep track of losses
1. Efficiency loss goes into heat of motor. Need an estimate of its thermal mass and then change in temp. (Motor temp new var)
1. Cleaner logging
1. Log everything without having to add more rows constantly
1. Keep efficiency in mind
1. Tractive system heat generation (Not acc)
1. Estimate how much heat is generated in the accumulator
1. Not high priority unless we can get more data.
1. Steering model
1. Suspension Model


# New simulation architecture idea


```python
# Dynamic Vars
posX = 0
posY = 1
velX = 2
velY = 3
accelX = 4
accelY = 5

arr = np.array((simSteps, 6+1))
arr[0] = step0

def step():
newPosX = step[posX] + step[velX] * t

for i in range(simSteps):
arr[i+1] = step(arr[i])
```

```python
# Dictionary Idea #1
# Array of dictionaries where each time step gets its own dictionary
# Trivial to access a specific thing from any row

arr = np.array((simSteps))
arr[0] = step0

def step():
newPosX = arr[posX] + arr[velX] * t

for i in range(simSteps):
arr[i+1] = step(arr[i])
```

```python
# Dictionary Idea #2
# Dictionary of arrays. Each key is a column and each array is the length of the simulation
# Trivial to access columns which is typically how we access data

arr = np.array((simSteps))
arr[0] = step0

def step():
newPosX = arr[posX] + arr[velX] * t

for i in range(simSteps):
arr[i+1] = step(arr[i])
```
21 changes: 10 additions & 11 deletions FullVehicleSim/Electrical/powertrain.py
Original file line number Diff line number Diff line change
@@ -1,24 +1,23 @@
from state import VehicleState
from paramLoader import Parameters, Magic
import numpy as np
from paramLoader import *

def calcMaxMotorTorque(worldPrev:VehicleState, resistiveForces:float, maxPower:float, maxTractionTorqueAtWheel:float):
def calcMaxMotorTorque(worldArray:np.ndarray, step:int, resistiveForces:float, maxPower:float, maxTractionTorqueAtWheel:float):
'''
Motor Torque at the wheel

minimum(rpm limited torque, power limited torque, perfect traction torque)
'''
## RPM Limited Torque (Motor Controller limits it to ~ this in practice. Maybe something more like 7490ish)
if worldPrev.motorRPM > 7490:
if worldArray[step-1, varMotorRPM] > 7490:
return -1 * resistiveForces * Parameters["wheelRadius"]
if worldPrev.motorRotationsHZ != 0: ## If rolling, torque may be power limited.
maxPowerTorque = maxPower / worldPrev.motorRotationsHZ * Parameters["gearRatio"]
if worldArray[step-1, varMotorRotationsHZ] != 0: ## If rolling, torque may be power limited.
maxPowerTorque = maxPower / worldArray[step-1, varMotorRotationsHZ] * Parameters["gearRatio"]
else: ## Avoid divide by 0 error but it's just the same as the max torque that the motor can deliver (180 Nm)
maxPowerTorque = 180.0 # Nm at 0 rpm
perfectTractionTorque = Parameters["maxTorque"]
torque = min(perfectTractionTorque, maxPowerTorque, maxTractionTorqueAtWheel/Parameters["gearRatio"])
torque = min(Parameters["maxTorque"], maxPowerTorque, maxTractionTorqueAtWheel/Parameters["gearRatio"])
return torque

def calcCurrent(power, voltage):
def calcCurrent(power:float, voltage:float) -> float:
if (power / voltage) > Parameters["tractiveIMax"]:
return Parameters["tractiveIMax"]
return power / voltage
Expand All @@ -29,10 +28,10 @@ def calcMaxWheelTorque(maxMotorTorque):
'''
return maxMotorTorque * Parameters["gearRatio"]

def calcMotorForce(maxWheelTorque):
def calcMotorForce(maxWheelTorque:float) -> float:
return (maxWheelTorque / Parameters["wheelRadius"])

def calcMaxPower(voltage):
def calcMaxPower(voltage:float) -> float:
return Parameters["tractiveIMax"] * voltage

def calcVoltage():
Expand Down
2 changes: 1 addition & 1 deletion FullVehicleSim/Electrical/tractiveBattery.py
Original file line number Diff line number Diff line change
@@ -1,2 +1,2 @@
from FullVehicleSim.paramLoader import Parameters, Magic
from FullVehicleSim.paramLoader import *

9 changes: 4 additions & 5 deletions FullVehicleSim/Mech/aero.py
Original file line number Diff line number Diff line change
@@ -1,9 +1,8 @@
import numpy as np
from paramLoader import Parameters, Magic
from state import VehicleState
from paramLoader import *

def calcDrag(prevWorld:VehicleState) -> float:
return 0.5 * Parameters["airDensity"] * Parameters["dragCoeffAreaCombo"] * prevWorld.speed**2
def calcDrag(worldArray:np.ndarray, step:int) -> float:
return 0.5 * Parameters["airDensity"] * Parameters["dragCoeffAreaCombo"] * worldArray[step-1, varSpeed]**2

def calcDownForce(prevWorld:VehicleState) -> np.ndarray:
def calcDownForce(worldArray:np.ndarray, step:int) -> np.ndarray:
return np.asarray([0,0,0,0], dtype=float)
28 changes: 14 additions & 14 deletions FullVehicleSim/Mech/braking.py
Original file line number Diff line number Diff line change
@@ -1,54 +1,54 @@
from Mech import brakepadFrictionModel
from paramLoader import Parameters, Magic
from paramLoader import *
import numpy as np
from state import VehicleState
# Docs:
# https://docs.google.com/document/d/1oGsGDnY0DEKWpE3S6481A9yZ0F9qUEwWkSXJwTSz4E4/edit?tab=t.2rmbsj26c7w
# The goal of these functions are to calculate the net force on the brakes, applied reverse to heading

def brakePSI_toNewtons(psi:float) -> float:
return psi * Parameters["brakeCaliperArea"] * 4.448222 # lb force to Newtons

def calcBrakeForce(prevWorld:VehicleState, inputs) -> tuple[float,float]:
def calcBrakeForce(worldArray:np.ndarray, step:int) -> tuple[float,float]:
"""
Calculate the brake force.

FrictionCoeff(temp) * maxBrakeForce * 4 (for 4 wheels)

:param prevWorld: World State Previous
:param worldArray: World State Array
:param step: Current step index
:return: Brake Force
"""
frontBrakePSI = inputs[1]
rearBrakePSI = inputs[2]
frontBrakePSI = worldArray[step, varBrakePressureFront]
rearBrakePSI = worldArray[step, varBrakePressureRear]
frontBrakeForce = brakePSI_toNewtons(frontBrakePSI)
rearBrakeForce = brakePSI_toNewtons(rearBrakePSI)

# Calculate Brake Force
frontBrakeForce:float = brakepadFrictionModel.calcFrictionCoeff(prevWorld.frontBrakeTemperature) * frontBrakeForce * 2 * Parameters["brakeDiscRadius"] / Parameters["wheelRadius"]
rearBrakeForce:float = brakepadFrictionModel.calcFrictionCoeff(prevWorld.rearBrakeTemperature) * rearBrakeForce * 2 * Parameters["brakeDiscRadius"] / Parameters["wheelRadius"]
frontBrakeForce:float = brakepadFrictionModel.calcFrictionCoeff(worldArray[step-1, varFrontBrakeTemperature]) * frontBrakeForce * 2 * Parameters["brakeDiscRadius"] / Parameters["wheelRadius"]
rearBrakeForce:float = brakepadFrictionModel.calcFrictionCoeff(worldArray[step-1, varRearBrakeTemperature]) * rearBrakeForce * 2 * Parameters["brakeDiscRadius"] / Parameters["wheelRadius"]
return frontBrakeForce, rearBrakeForce

def calcBrakeCooling(prevWorld:VehicleState) -> tuple[float,float]:
def calcBrakeCooling(worldArray:np.ndarray, step:int) -> tuple[float,float]:
"""
Calculate the cooled brake temperature.

:param prevWorld: World State
:return: Change in Temperature
"""
frontBrakeCooling = Parameters["ambientTemperature"] + (prevWorld.frontBrakeTemperature - Parameters["ambientTemperature"]) * np.e ** (-1 / Parameters["stepsPerSecond"]/50.2)
rearBrakeCooling = Parameters["ambientTemperature"] + (prevWorld.rearBrakeTemperature - Parameters["ambientTemperature"]) * np.e ** (-1 / Parameters["stepsPerSecond"]/50.2)
frontBrakeCooling = Parameters["ambientTemperature"] + (worldArray[step-1, varFrontBrakeTemperature] - Parameters["ambientTemperature"]) * np.e ** (-1 / Parameters["stepsPerSecond"]/50.2)
rearBrakeCooling = Parameters["ambientTemperature"] + (worldArray[step-1, varRearBrakeTemperature] - Parameters["ambientTemperature"]) * np.e ** (-1 / Parameters["stepsPerSecond"]/50.2)
return frontBrakeCooling, rearBrakeCooling
#q = (initTemperature - parameters["ambientTemperature"]) * parameters["brakeMass"] * parameters["brakeSpecificHeatCapacity"]
#change = (q * parameters["brakepadThickness"])/(initTemperature * parameters["brakeThermalConductivity"] * parameters["brakeSurfaceArea"]
#return initTemperature - change

def calcBrakeHeating(prevWorld:VehicleState, inputs) -> tuple[float,float]:
def calcBrakeHeating(worldArray:np.ndarray, step:int) -> tuple[float,float]:
# Calculate Brake Force
frontBrakeForce, rearBrakeForce = calcBrakeForce(prevWorld, inputs)
frontBrakeForce, rearBrakeForce = calcBrakeForce(worldArray, step)
# Guess energy increase based on kinetic energy decrease of the vehicle.
# Assumption is 100% of kinetic energy lost goes into brake heating.
speedChange = (frontBrakeForce + rearBrakeForce) / Parameters["Mass"] / Parameters["stepsPerSecond"] # momentum impulse
energyChange = 0.5 * Parameters["Mass"] * (prevWorld.speed - (prevWorld.speed - speedChange))
energyChange = 0.5 * Parameters["Mass"] * (worldArray[step-1, varSpeed] - (worldArray[step-1, varSpeed] - speedChange))
tempChange = energyChange/(Parameters["brakeMass"] * Parameters["brakeSpecificHeatCapacity"])

# While this doesn't seem physically intuitive, it is based on the idea that the front and rear brakes share heat based on their contribution to total braking force.
Expand Down
23 changes: 11 additions & 12 deletions FullVehicleSim/Mech/general.py
Original file line number Diff line number Diff line change
@@ -1,20 +1,19 @@
from Mech.traction import calcCorneringStiffness
from paramLoader import Parameters, Magic
from state import VehicleState
from paramLoader import *
from Mech.braking import calcBrakeForce
from Mech.aero import calcDrag
from Mech.steering import calcSlipAngle, calcYawRate
from Mech.tireLoad import calcLoadTransfer
import numpy as np

def calcResistiveForces(worldPrev:VehicleState, inputs):
if worldPrev.speed <= 1e-5: # Floating point error
def calcResistiveForces(worldArray:np.ndarray, step:int):
if worldArray[step-1, varSpeed] <= 1e-5: # Floating point error
return 0
else:
frontBrakeForce, rearBrakeForce = calcBrakeForce(worldPrev, inputs)
return -1 * (calcDrag(worldPrev) + frontBrakeForce + rearBrakeForce)
frontBrakeForce, rearBrakeForce = calcBrakeForce(worldArray, step)
return -1 * (calcDrag(worldArray, step) + frontBrakeForce + rearBrakeForce)

def calculateYawRate(prevWorld:VehicleState, steerAngle:float, initAcceleration:float, heading:np.ndarray, initYawRate:float, timeSinceLastSteer:float):
def calculateYawRate(worldArray:np.ndarray, step:int, initAcceleration:float, initYawRate:float, timeSinceLastSteer:float):
"""Calculate the yaw rate of the vehicle at the current state.
This function computes the yaw rate by calculating tire loads, slip angles,
cornering stiffness, and then applying the vehicle dynamics equations.
Expand All @@ -32,9 +31,9 @@ def calculateYawRate(prevWorld:VehicleState, steerAngle:float, initAcceleration:
-----
Slip ratio is fixed at 0.15.
"""
tireLoad = calcLoadTransfer(initAcceleration * heading[0], initAcceleration * heading[1], initYawRate)
slipAngle = calcSlipAngle(initYawRate, prevWorld.velocity, steerAngle, Parameters)
tireLoad = calcLoadTransfer(initAcceleration * worldArray[step-1, varHeadingX], initAcceleration * worldArray[step-1, varHeadingY], initYawRate)
slipAngle = calcSlipAngle(worldArray, step)
slipRatio = 0.15
corneringStiffness = calcCorneringStiffness(tireLoad, slipAngle, slipRatio, prevWorld.speed, 80, 40, Parameters, Magic) # Works but unused
res = calcYawRate(initYawRate, prevWorld.speed, steerAngle, timeSinceLastSteer, corneringStiffness[0], corneringStiffness[1], Parameters)
return res
corneringStiffness = calcCorneringStiffness(tireLoad, slipAngle, slipRatio, worldArray[step-1, varSpeed], 80, 40, Parameters, Magic) # Works but unused
res = calcYawRate(initYawRate, worldArray[step-1, varSpeed], worldArray[step, varSteerAngle], timeSinceLastSteer, corneringStiffness[0], corneringStiffness[1])
return res
Loading