Skip to content

Retrieval

The following functions are methods from the class ProfilesData

Aerosol Extinction

This module is used to calculate extinction profiles from single attenuated backscatter profiles using an a priori.

backward_inversion(data, iref, apriori, rayleigh)

Backward (Klett [Klett1985]_ ) inversion method.

.. [Klett1985] Klett, J. D. (1985). Lidar inversion with variable backscatter/extinction ratios. Applied optics, 24(11), 1638-1643.

Parameters:

Name Type Description Default
data array_like

1D array of a single profile of the attenuated backscatter coefficient.

 required
iref float

Index of the reference altitude returned by :func:aprofiles.retrieval.ref_altitude.get_iref().

 required
apriori dict

A priori values used to constrain the inversion. Valid keys include: - lr (float): Lidar Ratio (in sr). - aod (float): AOD (unitless).

 required
rayleigh

class:aprofiles.rayleigh.RayleighData): Instance of the RayleighData class.

 required

Raises:

Type Description
NotImplementedError

AOD apriori is not implemented yet.

Returns:

Type Description
array_like

Extinction coefficient (in m⁻¹).
Source code in aprofiles/retrieval/extinction.py 
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
def backward_inversion(data, iref, apriori, rayleigh):
    """
    Backward (Klett [Klett1985]_ ) inversion method.

    .. [Klett1985] Klett, J. D. (1985). Lidar inversion with variable backscatter/extinction ratios. Applied optics, 24(11), 1638-1643.

    Args:
        data (array_like): 1D array of a single profile of the attenuated backscatter coefficient.
        iref (float): Index of the reference altitude returned by :func:`aprofiles.retrieval.ref_altitude.get_iref()`.
        apriori (dict): A priori values used to constrain the inversion. Valid keys include:
            - **lr** (float): Lidar Ratio (in sr).
            - **aod** (float): AOD (unitless).
        rayleigh (:class:`aprofiles.rayleigh.RayleighData`): Instance of the `RayleighData` class.

    Raises:
        NotImplementedError: AOD apriori is not implemented yet.

    Returns:
        (array_like): Extinction coefficient (in m⁻¹).
    """


    if "aod" in apriori:
        # search by dichotomy the LR that matches the apriori aod
        raise NotImplementedError("AOD apriori is not implemented yet")
        lr_aer = 50
    else:
        # if we assume the LR, no need to minimize for matching aod
        lr_aer = apriori["lr"]
    lr_mol = 8.0 * np.pi / 3.0

    # vertical resolution
    dz = min(np.diff(rayleigh.altitude))

    int1_a = np.cumsum((lr_aer - lr_mol) * rayleigh.backscatter[:iref] * dz)
    int1_b = [2 * int1_a[-1] - 2 * int1_a[i] for i in range(iref)]
    phi = [np.log(abs(data[i] / data[iref])) + int1_b[i] for i in range(iref)]

    int2_a = 2 * np.nancumsum(lr_aer * np.exp(phi) * dz)
    int2_b = [int2_a[-1] - int2_a[i] for i in range(iref)]

    # initialize total backscatter
    back_aer_iref = 0  # m-1
    beta_tot_iref = rayleigh.backscatter[iref] + back_aer_iref

    # total backscatter
    beta_tot = [np.exp(phi[i]) / ((1 / beta_tot_iref) + int2_b[i]) for i in range(iref)]
    # aerosol backsatter (m-1.sr-1)
    beta_aer = beta_tot - rayleigh.backscatter[:iref]
    # aerosol extinction (m-1)
    sigma_aer = lr_aer * beta_aer
    # returns extinction in m-1 when valid, and np.nan elsewhere
    ext = [sigma_aer[i] if i < len(sigma_aer) else np.nan for i in range(len(data))]
    return ext

forward_inversion(data, iref, apriori, rayleigh)

Forward iterative inversion method [#]_.

.. [#] Li, D., Wu, Y., Gross, B., & Moshary, F. (2021). Capabilities of an Automatic Lidar Ceilometer to Retrieve Aerosol Characteristics within the Planetary Boundary Layer. Remote Sensing, 13(18), 3626.

Method principle:

At z0, the aerosol transmission is assumed as being close to 1. We evaluate the aerosol extinction based on this assumption. This evaluation gives a refined aerosol extinction that is used to calculate a second time the aerosol transmission. The aerosol extinction retrieval will converge after a certain number iterations. After the convergence, the aerosol extinction is retrieved in the next upper layer.

Parameters:

Name Type Description Default
data array_like

1D Array of single profile of attenuated backscatter coefficient, in m-1.sr-1.

 required
iref float

index of the reference altitude returned by :func:aprofiles.retrieval.ref_altitude.get_iref().

 required
apriori dict

A priori value to be used to constrain the inversion. Valid keys: ‘lr’ (Lidar Ratio, in sr) and ‘aod’ (unitless).

 required
rayleigh aprofiles.rayleigh.RayleighData). Instance of the

class:aprofiles.rayleigh.RayleighData class.

 required

Raises:

Type Description
NotImplementedError

AOD apriori is not implemented yet.

Returns:

Type Description
array_like

Extinction coefficient (in m⁻¹).
Source code in aprofiles/retrieval/extinction.py 
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
def forward_inversion(data, iref, apriori, rayleigh):
    """Forward iterative inversion method [#]_.

    .. [#] Li, D., Wu, Y., Gross, B., & Moshary, F. (2021). Capabilities of an Automatic Lidar Ceilometer to Retrieve Aerosol Characteristics within the Planetary Boundary Layer. Remote Sensing, 13(18), 3626.

    Method principle:

    At z0, the aerosol transmission is assumed as being close to 1. We evaluate the aerosol extinction based on this assumption.
    This evaluation gives a refined aerosol extinction that is used to calculate a second time the aerosol transmission.
    The aerosol extinction retrieval will converge after a certain number iterations.
    After the convergence, the aerosol extinction is retrieved in the next upper layer.

    Args:
        data (array_like): 1D Array of single profile of attenuated backscatter coefficient, in m-1.sr-1.
        iref (float): index of the reference altitude returned by :func:`aprofiles.retrieval.ref_altitude.get_iref()`.
        apriori (dict): A priori value to be used to constrain the inversion. Valid keys: ‘lr’ (Lidar Ratio, in sr) and ‘aod’ (unitless).
        rayleigh (aprofiles.rayleigh.RayleighData). Instance of the :class:`aprofiles.rayleigh.RayleighData` class.

    Raises:
        NotImplementedError: AOD apriori is not implemented yet.

    Returns:
        (array_like): Extinction coefficient (in m⁻¹).
    """

    if "aod" in apriori:
        # search by dichotomy the LR that matches the apriori aod
        raise NotImplementedError("AOD apriori is not implemented yet")
        lr_aer = 50
    elif "lr" in apriori:
        # if we assume the LR, no need to minimize for matching aod
        lr_aer = apriori["lr"]
    lr_mol = 8.0 * np.pi / 3.0

    def _get_aer_at_i(data, i, Tm, Bm, Ta, Ba, Ea, dz, nloop_max=30, diff_ext=0.01):
        # suppress warnings
        warnings.filterwarnings('ignore')

        for _ in range(nloop_max):
            if np.isnan(Ea[i]):
                Ta[i] = 1
            else:
                Ta[i] = np.exp(-np.nancumsum(Ea * dz))[i]
            if (Tm[i] ** 2 * Ta[i] ** 2) != 0:
                Ba[i] = data[i] / (Tm[i] ** 2 * Ta[i] ** 2) - Bm[i]
            else:
                Ba[i] = np.nan
            # test extinction
            if 1 - (lr_aer * Ba[i] / Ea[i]) < diff_ext:
                Ea[i] = lr_aer * Ba[i]
                break
            Ea[i] = lr_aer * Ba[i]

        return Ba[i], Ea[i], Ta[i]

    # vertical resolution
    dz = min(np.diff(rayleigh.altitude))

    Bm = rayleigh.backscatter
    Em = rayleigh.extinction
    Tm = np.exp(-np.cumsum(Em * dz))

    # Initialize aer profiles
    Ta = np.asarray([np.nan for _ in range(len(data))])
    Ba = np.asarray([np.nan for _ in range(len(data))])
    Ea = np.asarray([np.nan for _ in range(len(data))])

    for i in range(iref):
        Ba[i], Ea[i], Ta[i] = _get_aer_at_i(data, i, Tm, Bm, Ta, Ba, Ea, dz)
        if np.isnan(Ba[i]):
            continue
    # returns extinction in m-1
    ext = Ea
    return ext

inversion(profiles, time_avg=1, zmin=4000.0, zmax=6000.0, min_snr=0.0, under_clouds=False, method='forward', apriori={'lr': 50.0}, remove_outliers=False, verbose=False)

Aerosol inversion of the attenuated backscatter profiles using an apriori.

Parameters:

Name Type Description Default
profiles ProfilesData

ProfilesData instance.

 required
time_avg int

in minutes, the time during which we aggregate the profiles before inverting the profiles.

 1
zmin float

minimum altitude AGL, in m, for looking for the initialization altitude.

 4000.0
zmax float

maximum altitude AGL, in m, for looking for the initialization altitude.

 6000.0
min_snr float

Minimum SNR required at the reference altitude to be valid.

 0.0
under_clouds bool

If True, and if the ProfilesData has a cloud_base variable (returned by the clouds method), forces the initialization altitude to be found below the first cloud detected in the profile.

 False
method str

must be in [‘backward’, ‘forward’].

 'forward'
apriori dict

A priori value to be used to constrain the inversion. Valid keys: ‘lr’ (Lidar Ratio, in sr), ‘aod’ (unitless), 'cfg' (path of config file).

 {'lr': 50.0}
remove_outliers bool

Remove profiles considered as outliers based on aod calculation (AOD<0, or AOD>2).

 False
verbose bool

verbose mode.

 False

Raises:

Type Description
NotImplementedError

AOD apriori is not implemented yet.

Returns:

Type Description
ProfilesData

object with additional (xarray.DataArray):

  • 'extinction' (time, altitude): 2D array corresponding to the aerosol extinction, in km-1.
  • 'aod' (time): 1D array corresponding to the aerosol optical depth associated to the extinction profiles.
  • 'lidar_ratio' (time): 1D array corresponding to the lidar ratio associated to the extinction profiles.
Example

Profiles preparation 

import aprofiles as apro
#read example file
path = "examples/data/L2_0-20000-001492_A20210909.nc"
reader = apro.reader.ReadProfiles(path)
profiles = reader.read()
#extrapolate lowest layers
profiles.extrapolate_below(z=150, inplace=True)

Backward inversion 

#aerosol inversion
profiles.inversion(zmin=4000, zmax=6000, remove_outliers=False, method='backward')
#plot extinction profiles
profiles.plot(var='extinction', zmax=6000, vmin=0, vmax=5e-2)

Extinction profiles retrieved with the backward method

Forward inversion 

#aerosol inversion
profiles.inversion(zmin=4000, zmax=6000, remove_outliers=False, method='forward')
#plot extinction profiles
profiles.plot(var='extinction', zmax=6000, vmin=0, vmax=5e-2)

Extinction profiles retrieved with the forward method

Source code in aprofiles/retrieval/extinction.py 
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
def inversion(
    profiles,
    time_avg=1,
    zmin=4000.0,
    zmax=6000.0,
    min_snr=0.0,
    under_clouds=False,
    method="forward",
    apriori={"lr": 50.0},
    remove_outliers=False,
    verbose=False,
):
    """Aerosol inversion of the attenuated backscatter profiles using an apriori.

    Args:
        profiles (aprofiles.profiles.ProfilesData): `ProfilesData` instance.
        time_avg (int, optional): in minutes, the time during which we aggregate the profiles before inverting the profiles.
        zmin (float, optional): minimum altitude AGL, in m, for looking for the initialization altitude.
        zmax (float, optional): maximum altitude AGL, in m, for looking for the initialization altitude.
        min_snr (float, optional): Minimum SNR required at the reference altitude to be valid.
        under_clouds (bool, optional): If True, and if the `ProfilesData` has a `cloud_base` variable (returned by the `clouds` method), forces the initialization altitude to be found below the first cloud detected in the profile.
        method (str, optional): must be in [‘backward’, ‘forward’].
        apriori (dict, optional): A priori value to be used to constrain the inversion. Valid keys: ‘lr’ (Lidar Ratio, in sr), ‘aod’ (unitless), 'cfg' (path of config file).
        remove_outliers (bool, optional): Remove profiles considered as outliers based on aod calculation (AOD<0, or AOD>2).
        verbose (bool, optional): verbose mode.

    Raises:
        NotImplementedError: AOD apriori is not implemented yet.

    Returns:
        (aprofiles.profiles.ProfilesData):
            object with additional (xarray.DataArray):

            - `'extinction' (time, altitude)`: 2D array corresponding to the aerosol extinction, in km-1.
            - `'aod' (time)`: 1D array corresponding to the aerosol optical depth associated to the extinction profiles.
            - `'lidar_ratio' (time)`: 1D array corresponding to the lidar ratio associated to the extinction profiles.

    Example:
        Profiles preparation
        ```python
        import aprofiles as apro
        #read example file
        path = "examples/data/L2_0-20000-001492_A20210909.nc"
        reader = apro.reader.ReadProfiles(path)
        profiles = reader.read()
        #extrapolate lowest layers
        profiles.extrapolate_below(z=150, inplace=True)
        ```

        Backward inversion
        ```python
        #aerosol inversion
        profiles.inversion(zmin=4000, zmax=6000, remove_outliers=False, method='backward')
        #plot extinction profiles
        profiles.plot(var='extinction', zmax=6000, vmin=0, vmax=5e-2)
        ```

        ![Extinction profiles retrieved with the backward method](../../assets/images/backward.png)

        Forward inversion
        ```python
        #aerosol inversion
        profiles.inversion(zmin=4000, zmax=6000, remove_outliers=False, method='forward')
        #plot extinction profiles
        profiles.plot(var='extinction', zmax=6000, vmin=0, vmax=5e-2)
        ```

        ![Extinction profiles retrieved with the forward method](../../assets/images/forward.png)
    """

    # we work on profiles averaged in time to reduce the noise
    rcs = profiles.time_avg(
        time_avg, var="attenuated_backscatter_0", inplace=False
    ).data.attenuated_backscatter_0

    """
    #if clouds detected, set to nan profiles where cloud is found below 4000m
    lowest_clouds = profiles._get_lowest_clouds()
    for i in range(len(profiles.data.time.data)):
        if lowest_clouds[i]<=4000:
            rcs.data[i,:] = [np.nan for _ in rcs.data[i,:]]
    """

    # if under_clouds, check if clouds_bases is available
    if under_clouds and "clouds_bases" in list(profiles.data.data_vars):
        lowest_clouds = profiles._get_lowest_clouds()
    elif under_clouds and "clouds_bases" not in list(profiles.data.data_vars):
        warnings.warn(
            "under_clouds parameter sets to True (defaults value) when the clouds detection has not been applied to ProfileData object."
        )
        lowest_clouds = [np.nan for i in np.arange(len(profiles.data.time))]
    else:
        lowest_clouds = [np.nan for i in np.arange(len(profiles.data.time))]

    # aerosol retrieval requires a molecular profile
    altitude = np.asarray(profiles.data.altitude.data)
    wavelength = float(profiles.data.l0_wavelength.data)
    rayleigh = apro.rayleigh.RayleighData(
        altitude, T0=298, P0=1013, wavelength=wavelength
    )

    # aerosol inversion
    ext, lr, aod, z_ref = [], [], [], []
    aod_min, aod_max = 0, 2
    vertical_resolution = profiles._get_resolution("altitude")
    for i in (track(range(len(profiles.data.time.data)), description="klett ",disable=not verbose)):
        # for this inversion, it is important to use the right units
        if "Mm" in profiles.data.attenuated_backscatter_0.units:
            calibrated_data = rcs.data[i, :] * 1e-6 # conversion from Mm-1.sr-1 to m-1.sr-1
        else:
            calibrated_data = rcs.data[i, :]

        # reference altitude
        lowest_cloud_agl = lowest_clouds[i] - profiles.data.station_altitude.data
        imin = profiles._get_index_from_altitude_AGL(zmin)
        imax = profiles._get_index_from_altitude_AGL(np.nanmin([zmax, lowest_cloud_agl]))
        if method == "backward":
            iref = get_iref(rcs.data[i, :], imin, imax, min_snr)
        elif method == "forward":
            iref = imax
        z_ref.append(altitude[iref])

        if iref is not None:
            # aerosol inversion
            if method == "backward":
                _ext = backward_inversion(calibrated_data, iref, apriori, rayleigh)
            elif method == "forward":
                _ext = forward_inversion(calibrated_data, iref, apriori, rayleigh)
        else:
            _ext = [np.nan for _ in range(len(calibrated_data))]

        # add aod and lr
        if "aod" in apriori:
            # search by dichotomy the LR that matches the apriori aod
            raise NotImplementedError("AOD apriori is not implemented yet")
        else:
            # if we assume the LR, no need to minimize for matching aod
            _lr = apriori["lr"]
            _aod = np.nancumsum(np.asarray(_ext) * vertical_resolution)[-1]
            if remove_outliers and _aod < aod_min or remove_outliers and _aod > aod_max:
                lr.append(np.nan)
                aod.append(np.nan)
                ext.append([np.nan for x in _ext])
            else:
                lr.append(_lr)
                aod.append(_aod)
                ext.append(_ext)

    # creates dataarrays
    profiles.data["extinction"] = (("time","altitude"), np.asarray(ext) * 1e3)
    profiles.data["extinction"] = profiles.data.extinction.assign_attrs({
        'long_name': f"Extinction Coefficient @ {int(wavelength)} nm",
        'method': f"{method.capitalize()} Klett",
        'units': "km-1",
        'time_avg': time_avg,
        'zmin': zmin,
        'zmax': zmax,
        'apriori_variable': list(apriori.keys())[0],
        'apriori_value': apriori[list(apriori.keys())[0]],
        })

    profiles.data["aod"] = ("time", aod)
    profiles.data["aod"] = profiles.data.aod.assign_attrs({
        'long_name': f"Aerosol Optical Depth @ {int(wavelength)} nm",
        'unit': ''
    })

    profiles.data["lidar_ratio"] = ('time', lr)
    profiles.data["lidar_ratio"] = profiles.data.lidar_ratio.assign_attrs({
        'long_name': f"Lidar Ratio @ {int(wavelength)} nm",
        'units': 'sr',
        'use_cfg': str(apriori['use_cfg'])
    })
    if 'cfg' in apriori:
        profiles.data["lidar_ratio"] = profiles.data.lidar_ratio.assign_attrs({
            'cfg': str(apriori['cfg'])
        })

    profiles.data["z_ref"] = ('time', z_ref)
    profiles.data["z_ref"] = profiles.data.z_ref.assign_attrs({
        'long_name': f"Reference altitude ASL",
        'units': 'm'
    })

    return profiles

Mass Concentration

This module is used to calculate mass concentration profiles from extinction profiles for given aerosol types.

concentration_profiles(profiles, method, apriori)

Calculates Mass concentration profiles for different aerosol types

Parameters:

Name Type Description Default
profiles ProfilesData

ProfilesData instance.

 required
method str

Method for calculating MEC. Must be one of {"mortier_2013", "literature"}.

 required
apriori dict

Apriori mec value (m2.g-1).

 required

Returns:

Type Description
ProfilesData

object with additional (xarray.Dataset):

  • 'mass_concentration:<aer_type>'
Example

Profiles preparation 

import aprofiles as apro
# read example file
path = "examples/data/L2_0-20000-001492_A20210909.nc"
reader = apro.reader.ReadProfiles(path)
profiles = reader.read()
# extrapolate lowest layers
profiles.extrapolate_below(z=150, inplace=True)

Forward inversion 

# aerosol inversion
profiles.inversion(
    zmin=4000, zmax=6000, remove_outliers=False, method='forward', 
    method_mass_conc='mortier_2013'
)
# plot mass concentration profiles for urban particles
profiles.plot(var='mass_concentration:urban', zmax=6000, vmin=0, vmax=100)

Mass concentration profiles in the case of urban particles

Source code in aprofiles/retrieval/mass_conc.py 
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
def concentration_profiles(profiles, method, apriori):
    """Calculates Mass concentration profiles for different aerosol types

    Args:
        profiles (aprofiles.profiles.ProfilesData): `ProfilesData` instance.
        method (str): Method for calculating MEC. Must be one of {"mortier_2013", "literature"}.
        apriori (dict): Apriori mec value (m2.g-1).

    Returns:
        (aprofiles.profiles.ProfilesData):
            object with additional (xarray.Dataset):

            - `'mass_concentration:<aer_type>'`

    Example:
        Profiles preparation
        ```python
        import aprofiles as apro
        # read example file
        path = "examples/data/L2_0-20000-001492_A20210909.nc"
        reader = apro.reader.ReadProfiles(path)
        profiles = reader.read()
        # extrapolate lowest layers
        profiles.extrapolate_below(z=150, inplace=True)
        ```

        Forward inversion
        ```python
        # aerosol inversion
        profiles.inversion(
            zmin=4000, zmax=6000, remove_outliers=False, method='forward', 
            method_mass_conc='mortier_2013'
        )
        # plot mass concentration profiles for urban particles
        profiles.plot(var='mass_concentration:urban', zmax=6000, vmin=0, vmax=100)
        ```

        ![Mass concentration profiles in the case of urban particles](../../assets/images/mass_conc-urban.png)
    """

    # read aer_properties.json files
    f = open(Path(Path(__file__).parent,'..','config','aer_properties.json'))
    aer_properties = json.load(f)
    f.close()
    # check if the aer_type exist in the json file
    aer_types = aer_properties.keys()

    # get wavelength
    wavelength = float(profiles.data.l0_wavelength.data)
    if profiles.data.l0_wavelength.units != "nm":
        raise ValueError("wavelength units is not `nm`.")

    for aer_type in aer_types:
        # calculates mec
        mec = apro.mec.MECData(aer_type, wavelength, method)

        # compute mass_concentration profile. Use extinction as base.
        mass_concentration = profiles.data.extinction * 1e-3 #conversion from km-1 to m-1
        # mass_concentration = copy.deepcopy(profiles.data.extinction)
        mass_concentration.data = np.divide(mass_concentration, mec.mec)
        # # conversion from g.m-3 to µg.m-3
        mass_concentration.data = mass_concentration.data * 1e6

        # creates dataset with a dataarray for each aer_type
        profiles.data[f"mass_concentration:{aer_type}"] = (('time', 'altitude'), mass_concentration.data)
        profiles.data[f"mass_concentration:{aer_type}"] = profiles.data[f"mass_concentration:{aer_type}"].assign_attrs({
            'long_name': f"Mass concentration [{aer_type.replace('_', ' ')} particles]",
            'units': 'µg.m-3',
            'mec': mec.mec,
        })

    # add ifs mec
    if apriori["mec"]:
        aer_type = "ifs"
        # compute mass_concentration profile. Use extinction as base.
        mass_concentration = profiles.data.extinction * 1e-3 #conversion from km-1 to m-1
        # mass_concentration = copy.deepcopy(profiles.data.extinction)
        mass_concentration.data = np.divide(mass_concentration, apriori["mec"])
        # # conversion from g.m-3 to µg.m-3
        mass_concentration.data = mass_concentration.data * 1e6

        # creates dataset with a dataarray for each aer_type
        profiles.data[f"mass_concentration:{aer_type}"] = (('time', 'altitude'), mass_concentration.data)
        profiles.data[f"mass_concentration:{aer_type}"] = profiles.data[f"mass_concentration:{aer_type}"].assign_attrs({
            'long_name': f"Mass concentration [{aer_type.replace('_', ' ')}]",
            'units': 'µg.m-3',
            'mec': apriori["mec"],
        })

    return profiles