How to interpolate between fixed values in a 2D grid for biome-generation?

I am currently implementing biome-generation in my game and I want to make the type of a biome dependent on humidity and temperature values that are generated from gradual noise.

different Biomes should have different heights and without interpolation this would result in abrupt height differences on biome borders as expected.

What I tried was to get the 2 neighbour biomes in the grid too and measure the blend-percentage of each biome.

Later I then get the 3 different height values from the biomes and multiply each with it's respective blend value.

Here is the simplified and stripped code, which I use to fetch biomes:

const BIOME_GRID_SIZE: usize = 2;
const BIOME_INDEX_SIZE: usize = BIOME_GRID_SIZE - 1;
const BIOME_GRID: [[Biomes; BIOME_GRID_SIZE]; BIOME_GRID_SIZE] =
    [
        [Biomes::Mountains, Biomes::Forest],
        [Biomes::Desert   , Biomes::Mesa  ],
    ];

fn get_height(coord: [i64; 2], noise: &Noise) -> i64 {
    let temperature = (noise.get_2d(coord) + 0.5).clamp(0.0, 1.0);
    let humidity = (noise.get_2d(coord /* + some offset */) + 0.5).clamp(0.0, 1.0);

    let x = BIOME_GRID_SIZE as f64 * humidity;
    let y = BIOME_GRID_SIZE as f64 * temperature;

    let x_frac = (x.fract() - 0.5) * 2.0;
    let y_frac = (y.fract() - 0.5) * 2.0;
    let x_blending = x_frac.abs();
    let y_blending = y_frac.abs();
    let own_blending = 2.0 - x_blending - y_blending;

    // direction of neighbour biomes
    let x_direction = x_frac.signum() as isize;
    let y_direction = y_frac.signum() as isize;

    let x_index = (x.trunc() as isize).clamp(0, BIOME_INDEX_SIZE as isize);
    let y_index = (y.trunc() as isize).clamp(0, BIOME_INDEX_SIZE as isize);
    let biomes = get_biomes(x_index, y_index, x_direction, y_direction);

    blend(
        coord,
        noise,
        biomes,
        [
            own_blending,
            x_blending,
            y_blending,
        ]
    ),
}

// get main and neighbour biomes
fn get_biomes(x: isize, y: isize, x_direction: isize, y_direction: isize) -> [Biomes; 3] {
    let mut biomes = [Biomes::Ocean; 3];

    for (i, (d_x, d_y)) in [(0, 0), (x_direction, 0), (0, y_direction)].iter().enumerate() {
        let x_index = (x + d_x).clamp(0, BIOME_INDEX_SIZE as isize) as usize;
        let y_index = (y + d_y).clamp(0, BIOME_INDEX_SIZE as isize) as usize;
        let biome = BIOME_GRID[x_index][y_index];
        biomes[i] = biome;
    }

    biomes
}

pub fn blend(
    coord: [i64; 2],
    noise: &Noise,
    biomes: [Biomes; 4],
    blending: [f64; 4],
) -> i64 {
    let heights: Vec<f64> = biomes
        .iter()
        .map(|x| x.get().height(coord, noise) as f64)
        .collect();

    let height = heights[0] * blending[0] + heights[1] * blending[1] + heights[2] * blending[2];
    let height = height as i64;

    height
}

This works well in some cases, in the other it fails completely.

I am unsure, if 2 neighbours are enough and how to properly get the blend values.

Is there a better solution to this problem?

In general for bilinear blending you would use four points. If I understand your code correctly that would be the four height maps for each biome. You then lerp across one axis (e.g. humidity) for the two pairs with the same other axis, and then lerp the two blended values again with the other axis (e.g. temperature).