Peak deconvolution knowing the weight function in python

Question

I am trying to calculate the time constants of the following theoretical Foster model using PySpice. I got the output response and its derivative respect to z [z=ln(t)]. Foster model

Nevertheless, I am not capable of performing the deconvolution of this curve to obtain the peaks. I know in this case the peaks can be detected without convolution, but this is just a test. The weight function is known, which is w(z) = exp^(z-exp^(z)) and has the following form: weight function

I have to perform the deconvolution of da/dz (derivative of otuput response respect to z) to get the peaks and the time constants of the Foster model. I tried several things such as Understanding scipy deconvolve or Python deconvolution giving unexpected result, but nothing is working. I am not very familiar with deconvolution maths, therefore I do not know if I am doing something wrong. Both output response and its derivative are correct, but the deconvolution is not. Result

How should I perform this deconvolution in python? The following code is the actual one. Thank you ver much for your help, I tried to summarize the problem as much as I could.

# ---MODULE IMPORT---
import numpy as np
import matplotlib.pyplot as plt
import matplotlib.ticker as ticker
import scipy.fft
from scipy.fft import fft, ifft
import PySpice.Logging.Logging as Logging
from PySpice.Spice.Netlist import Circuit
from PySpice.Spice.Netlist import SubCircuit
from PySpice.Unit import *

logger = Logging.setup_logging()

####### FUNCTION DEFINITION ######
def conv(input):
    return np.exp(input-np.exp(input))

####### CREATION OF RC BLOCK ######

class Foster_RC(SubCircuit):
    __nodes__ = ('t_1', 't_2')
    def __init__(self, name, r=1@u_kOhm, c=10@u_uF):
        SubCircuit.__init__(self, name, *self.__nodes__)
        self.C(1, 't_1', 't_2', c)
        self.R(1, 't_1', 't_2', r)

class Cauer_RC(SubCircuit):
    __nodes__ = ('t_1', 't_2')
    def __init__(self, name, r=1@u_kOhm, c=10@u_uF):
        SubCircuit.__init__(self, name, *self.__nodes__)
        self.C(1, 't_1', circuit.gnd, c)
        self.R(1, 't_1', 't_2', r)
        
###### PARAMETERS ######

R1, R2, R3 = 2@u_kOhm, 1@u_kOhm, 5@u_kOhm

C1, C2, C3 = 0.1@u_uF, 50@u_uF, 0.5@u_mF

###### MODEL CONSTRUCTION ######

circuit = Circuit('Foster model')
circuit.subcircuit(Cauer_RC(name='Foster_BLOCK_1 ', r=R1, c=C1))
circuit.subcircuit(Cauer_RC(name='Foster_BLOCK_2 ', r=R2, c=C2))
circuit.subcircuit(Cauer_RC(name='Foster_BLOCK_3 ', r=R3, c=C3))

circuit.PulseCurrentSource('input', circuit.gnd, 'n_test', initial_value=0@u_A,
                           pulsed_value=1@u_A, pulse_width=30@u_s, period=30@u_s, delay_time=0@u_s)

circuit.V('ammeter', 'n_test', 'n_in', 0@u_V)
circuit.X(1, 'Foster_BLOCK_1', 'n_in', 'n_1')
circuit.X(2, 'Foster_BLOCK_2', 'n_1', 'n_2')
circuit.X(3, 'Foster_BLOCK_3', 'n_2', circuit.gnd)


###### SIMULATION RUN ######

simulator = circuit.simulator()

analysis = simulator.transient(step_time=0.1@u_ms, end_time=30@u_s)

###### CALCULATION ######

t = np.array(analysis.time)
response = np.array(analysis['n_in'])[1:]
z = np.log(t[1:])

derivative = []
for i in range(len(response)):
    if i+1 == len(response): continue
    diff = (response[i+1]-response[i])/abs(z[i+1]-z[i])
    derivative.append(diff)

derivative = np.array(derivative)
conv_function = conv(np.linspace(-4, 1.5, len(derivative)))

recovered = ifft(fft(derivative)/fft(conv_function))

###### PLOTTING ######

fig, axs =plt.subplots(nrows=3, ncols=1, sharex=True)

axs[0].plot(t[1:], response)
axs[1].plot(t[1:-1], derivative)
axs[2].plot(np.linspace(t[1], t[-1], len(recovered)), recovered)

axs[0].grid()
axs[0].grid(which='minor', color='#CCCCCC', linestyle='--')
axs[1].grid()
axs[1].grid(which='minor', color='#CCCCCC', linestyle='--')
axs[2].grid()
axs[2].grid(which='minor', color='#CCCCCC', linestyle='--')

axs[2].set_xlabel('Time[s]', fontsize=12)
axs[2].set_xscale('log')

axs[0].set_ylabel('a(t)', fontsize=12)
axs[1].set_ylabel('da/dz', fontsize=12)
axs[2].set_ylabel('Deconv(da/dz)', fontsize=12)

axs[2].set_xlim(min(t[1:]), max(t))

axs[0].yaxis.set_major_locator(ticker.MultipleLocator(1000))
axs[0].yaxis.set_minor_locator(ticker.MultipleLocator(500))
axs[1].yaxis.set_major_locator(ticker.MultipleLocator(500))
axs[1].yaxis.set_minor_locator(ticker.MultipleLocator(250))


fig.suptitle('Time constant calculation', fontweight='bold', fontsize=16)

plt.show()

I tried different ways of deconvolution, but I am not getting a proper peak detection.