Skip to content

Commit

Permalink
change PWM discrete integration once more
Browse files Browse the repository at this point in the history 
collect all updates until a state flip happens (or a lot of time passed), now also for lamps

then chunk up the collected time in equal timesteps, being as close to 1ms as possible (unless <1ms passed, this uses that small timestep then, there is no other way)

integrate as before, but now no more tiny leftovers need to be processed that could significantly falsify the IIR and lamp/filament temperature computations (that all rely on roughly equal timesteps)
  • Loading branch information
@toxieainc
toxieainc committed Mar 11, 2024
1 parent a4a9c22 commit 4bfa11f
Show file tree
Hide file tree
Showing 2 changed files with 48 additions and 34 deletions.
78 changes: 46 additions & 32 deletions src/wpc/core.c
View file
Expand Up @@ -2490,9 +2490,7 @@ INLINE float cube(const float x)
return x * x * x;
}

static const double BULB_INTEGRATION_PERIOD = 0.001; // do the integration in a loop of small steps

INLINE void core_eye_flicker_fusion(core_tPhysicOutput* output, const double dt, float emission)
INLINE void core_eye_flicker_fusion(core_tPhysicOutput* output, const float emission)
{
// Compute the perceived emission using a hacky simplified eye integration model
// We want to model the flicker-fusion eye/brain behavior while keeping the output latency low (filter delay) and limit the flicker (still keeping it, if it can be seen on the real machine)
Expand All @@ -2507,11 +2505,6 @@ INLINE void core_eye_flicker_fusion(core_tPhysicOutput* output, const double dt,
// By now we also use a log function to map intensities (resulting range [~0.07..~0.11], kinda magic). This is according to the Ferry-Porter law, where higher intensities need a larger frequency in order to make
// flicker fusion work, and lower intensities already work at smaller ones. (main ff frequency range is roughly 20Hz to 70Hz, BUT with a lot of outlier situations)

if(dt != BULB_INTEGRATION_PERIOD)
{
const float ratio = (float)(dt*(1./BULB_INTEGRATION_PERIOD));
emission *= /*powf(ratio,(float)(1./1.3)) approx:*/ cube(sqrtf(sqrtf(ratio))); // pow boosts the emission artificially a bit, especially for GTS3/WPC alphaseg, as the power-on cycle is somehow too short?!?
}
const float eye_integration_old[3] = { output->state.bulb.eye_integration[0], output->state.bulb.eye_integration[1], output->state.bulb.eye_integration[2] };
const float eyeIntegrationFactor = /*0.1f*/ (0.07f + 0.02f * /*logf(output->value*6.f + 1.f) approx:*/ 2.f*sqrtf(fmaxf(output->value,0.f))), revEyeIntegrationFactor = 1.0f - eyeIntegrationFactor;
output->state.bulb.eye_integration[0] = (eyeIntegrationFactor * 0.5f) * (emission +output->state.bulb.eye_emission_old) + revEyeIntegrationFactor * eye_integration_old[0];
Expand All @@ -2522,31 +2515,47 @@ INLINE void core_eye_flicker_fusion(core_tPhysicOutput* output, const double dt,

#ifdef LOG_PWM_OUT
if (output == &coreGlobals.physicOutputState[LOG_PWM_OUT])
printf("Eye flicker fusion dt=%8.15f out=%0.5f in=%0.5f\n", dt, output->value, emission);
printf("Eye flicker fusion out=%0.5f in=%0.5f\n", output->value, emission);
#endif
}

static const double BULB_INTEGRATION_PERIOD = 0.001; // do the integration in a loop of small steps

// Incandescent bulb model based on Dulli Chandra Agrawal's and others publications (Heating-times of tungsten filament incandescent lamps).
// The bulb is a varying resistor depending on filament temperature, which is heated by the current (Ohm's law)
// and cooled by radiating energy (Planck & Stefan/Boltzmann laws). The visible emission power is then evaluated from the filament temperature.
// Integration is performed after a 1ms delay and on each state flip (since the behavior is highly non linear and there isn't any driver that would trigger this at a high frequency)
// Integration is performed after (roughly) 1ms delay in (roughly) 1ms cycles
// isFlip signals a state flip, BUT we track this internally anyways now
void core_update_pwm_output_bulb(const double now, const int index, const int isFlip)
{
core_tPhysicOutput* output = &coreGlobals.physicOutputState[index];
const float U = output->state.bulb.U * (float)(((coreGlobals.binaryOutputState[index >> 3] >> (index & 7)) & 1) ^ output->state.bulb.isReversed);
while (now - output->state.bulb.integrationTimestamp > 0.) {
const float dt = (float)(now - output->state.bulb.integrationTimestamp) > (float)BULB_INTEGRATION_PERIOD ? (float)BULB_INTEGRATION_PERIOD : (float)(now - output->state.bulb.integrationTimestamp);
if (dt < (float)BULB_INTEGRATION_PERIOD && !isFlip) // Don't perform too short integration periods unless it is the initial pulse start/end
return;
// Keeps T within the range of the LUT (between room temperature and melt down point)
output->state.bulb.filament_temperature = output->state.bulb.filament_temperature < 293.0f ? 293.0f : output->state.bulb.filament_temperature > (float) BULB_T_MAX ? (float) BULB_T_MAX : output->state.bulb.filament_temperature;
const float Ut = output->state.bulb.isAC ? (1.41421356f * sinf((float)(60.0 * 2.0 * PI) * (float)(output->state.bulb.integrationTimestamp - coreGlobals.lastACZeroCrossTimeStamp)) * U) : U;
const float dT = dt * (float) bulb_heat_up_factor(output->state.bulb.bulb, output->state.bulb.filament_temperature, Ut, output->state.bulb.serial_R);
output->state.bulb.filament_temperature += dT < 1000.0f ? dT : 1000.0f; // Limit initial current surge (1ms is a bit long when emulating this part of the heating)
core_eye_flicker_fusion(output, dt, (float) bulb_filament_temperature_to_emission(output->state.bulb.filament_temperature));
output->state.bulb.integrationTimestamp += BULB_INTEGRATION_PERIOD;
const float U = (((coreGlobals.binaryOutputState[index >> 3] >> (index & 7)) & 1) ^ output->state.bulb.isReversed) ? output->state.bulb.U : 0.f;
const float dt_diff = (float)(output->state.bulb.integrationTimestamp - output->state.bulb.prevIntegrationTimestamp);

if(U != output->state.bulb.prevIntegrationValue // state flip?
|| dt_diff >= (float)(BULB_INTEGRATION_PERIOD*20.)) // or we waited long enough to get a stable discrete integration
{
// do the integration in a loop of small steps, roughly in BULB_INTEGRATION_PERIOD sized steps (but rounded up/down to have same sized cycles in here)
float countf = floorf(dt_diff*(float)(1./BULB_INTEGRATION_PERIOD) + 0.5f);
if(countf < 1.f) // single cycle/short pulses need to use the 'original' cycle time as a workaround
countf = 1.f;
const int count = (int)countf;
const float dt = dt_diff/countf;
for(int i = 0; i < count; ++i) {
// Keeps T within the range of the LUT (between room temperature and melt down point)
output->state.bulb.filament_temperature = output->state.bulb.filament_temperature < 293.0f ? 293.0f : output->state.bulb.filament_temperature > (float) BULB_T_MAX ? (float) BULB_T_MAX : output->state.bulb.filament_temperature;
const float Ut = output->state.bulb.isAC ? (1.41421356f * sinf((float)(60.0 * 2.0 * PI) * (float)(output->state.bulb.prevIntegrationTimestamp - coreGlobals.lastACZeroCrossTimeStamp)) * output->state.bulb.prevIntegrationValue) : output->state.bulb.prevIntegrationValue;
const float dT = dt * (float)bulb_heat_up_factor(output->state.bulb.bulb, output->state.bulb.filament_temperature, Ut, output->state.bulb.serial_R);
output->state.bulb.filament_temperature += dT < 1000.0f ? dT : 1000.0f; // Limit initial current surge (1ms is a bit long when emulating this part of the heating)
core_eye_flicker_fusion(output, (float)bulb_filament_temperature_to_emission(output->state.bulb.filament_temperature));
output->state.bulb.prevIntegrationTimestamp += dt;
}
output->state.bulb.prevIntegrationTimestamp = output->state.bulb.integrationTimestamp;
output->state.bulb.prevIntegrationValue = U;
}
// else: wait with the integration until some more time vanished or the state will flip
output->state.bulb.integrationTimestamp = now;

#ifdef LOG_PWM_OUT
if (index == LOG_PWM_OUT)
printf("Output #%d t=%8.5f T=%5.0f e=%0.3f V=%0.3f S=%s\n", index, now, output->state.bulb.filament_temperature, bulb_filament_temperature_to_emission(output->state.bulb.filament_temperature), output->value, ((coreGlobals.binaryOutputState[index >> 3] >> (index & 7)) & 1) ? "x" : "-");
Expand All @@ -2564,23 +2573,28 @@ void core_update_pwm_output_bulb(const double now, const int index, const int is
// - https://en.wikipedia.org/wiki/Phosphor has a large phosphor chart
// - https://www.noritake-elec.com/technology/general-technical-information/vfd-operation states a linear relation between PWM and relative brightness
// - http://www.futabahk.com.hk/FTS/-%20FTP_DataBase/Docs_05TechDocs/VFD/AN-1103A-EN%20%20VFD%20Characteristics%20and%20Operation%20Guide%20(EN).pdf states the same
// isFlip signals a state flip, BUT we track this internally anyways now
void core_update_pwm_output_led(const double now, const int index, const int isFlip)
{
core_tPhysicOutput* output = &coreGlobals.physicOutputState[index];
const int state = (coreGlobals.binaryOutputState[index >> 3] >> (index & 7)) & 1;
const float power = state ? output->state.bulb.relative_brightness : 0.f;
const float dt_diff = (float)(output->state.bulb.integrationTimestamp - output->state.bulb.prevIntegrationTimestamp);

if(power != output->state.bulb.prevIntegrationValue // state flip?
|| output->state.bulb.integrationTimestamp - output->state.bulb.prevIntegrationTimestamp >= BULB_INTEGRATION_PERIOD*20.) // or we waited long enough to get a stable discrete integration (=no 'too short' time frames)
if(power != output->state.bulb.prevIntegrationValue // state flip?
|| dt_diff >= (float)(BULB_INTEGRATION_PERIOD*20.)) // or we waited long enough to get a stable discrete integration
{
double dt;
while ((dt = output->state.bulb.integrationTimestamp - output->state.bulb.prevIntegrationTimestamp) > 0.) {
if (dt > BULB_INTEGRATION_PERIOD) // do the integration in a loop of small steps
dt = BULB_INTEGRATION_PERIOD;
if (dt == BULB_INTEGRATION_PERIOD || output->state.bulb.prevIntegrationValue != 0.f) // if there is a 'leftover' from the discrete integration, just use it if it is in the On-state, otherwise it falsifies the integration more than it helps
core_eye_flicker_fusion(output, dt, output->state.bulb.prevIntegrationValue);
output->state.bulb.prevIntegrationTimestamp += BULB_INTEGRATION_PERIOD;
}
// do the integration in a loop of small steps, roughly in BULB_INTEGRATION_PERIOD sized steps (but rounded up/down to have same sized cycles in here)
float countf = floorf(dt_diff*(float)(1./BULB_INTEGRATION_PERIOD) + 0.5f);
if(countf < 1.f) // single cycle/short pulses need to use the 'original' cycle time as a workaround
countf = 1.f;
const int count = (int)countf;
const float dt = dt_diff/countf;
const float dt_ratio = dt*(float)(1./BULB_INTEGRATION_PERIOD);
const float emission = output->state.bulb.prevIntegrationValue * /*powf(dt_ratio,(float)(1./1.3)) approx:*/ cube(sqrtf(sqrtf(dt_ratio))); // pow boosts the emission artificially a bit, especially for GTS3/WPC alphaseg, as the power-on cycle is somehow too short?!?
for(int i = 0; i < count; ++i)
core_eye_flicker_fusion(output, emission);

// No real need to do this anymore, rather let the low pass filter over/undershoot to be more 'correct'
/*if (output->state.bulb.prevIntegrationValue != 0.f && output->value > (float)(254./255.)) // Stabilize Steady-On state
output->value = 1.f;
Expand Down
4 changes: 2 additions & 2 deletions src/wpc/core.h
View file
Expand Up @@ -434,8 +434,8 @@ typedef struct {
struct
{
int bulb; /* bulb model */
int isAC; /* AC or DC ? */
int isReversed; /* is input reversed ? */
int isAC; /* AC or DC ? */ // bool
int isReversed; /* is input reversed ? */ // bool
float U; /* voltage (Volts) */
float serial_R; /* serial resistor (Ohms) */
float relative_brightness; /* Relative brightness scale */
Expand Down

0 comments on commit 4bfa11f

Please sign in to comment.