import numpy as np
import matplotlib.pyplot as plt
# Time, in seconds
time = np.array(range(-10, 1000))
noiseStdDev = 3
externalTemp = np.ones(time.size) * 25 # 25 C
externalTemp[time >= 0] = 100 # Probe it put into boiling water at t=0
# Now let's simulate the thermocouple
# Assume the rise in temperature is proportional to the temperature difference
internalTemp = np.zeros(time.size) # Initialize the whole array to zero
internalTemp[0] = externalTemp[0] # Assume it's stable at room temperature
measuredTemp = np.zeros(time.size)
for i in range(len(time) - 1):
dt = externalTemp[i] - internalTemp[i]
tempRise = dt * 0.1
internalTemp[i + 1] = internalTemp[i] + tempRise
measuredTemp[i + 1] = internalTemp[i + 1] + noiseStdDev * np.random.randn(1)
# Now try to filter
filteredTemp = np.zeros(time.size)
for i in range(10, len(time)):
filteredTemp[i] = (measuredTemp[i] + measuredTemp[i-1] + measuredTemp[i-2] + measuredTemp[i-3] + measuredTemp[i-4] + \
measuredTemp[i-5] + measuredTemp[i-6] + measuredTemp[i-7] + measuredTemp[i-8] + measuredTemp[i-9]) / 10
# Make the plot
plt.plot(time, externalTemp, time, measuredTemp, time, filteredTemp)
plt.xlabel('Time [seconds]')
plt.ylabel('Temperature [C]')
plt.legend(['External', 'Measured', 'Filtered'])
plt.show()