jf_web
  • Home
  • About
  • Research
  • Teaching
  • Outreach
  • -misc-

On this page

  • Speckle as a wave interference phenomenon

jf_web

Welcome!

Everything on this website was built and programmed by me, some of it even in languages and libraries I understand! However, it would not have been possible without the immense body of literature and work made freely available by countless people all over the world. In that spirit, the entire code for this website is freely available here.

See below some choice interactive course notes – check out the Teaching tab for more.

Speckle as a wave interference phenomenon

Let’s look an illustrative real-world manifestation of this conceptual framework. SAR images are subject to speckle, a type of image distortion or noise which affects pixels in a stochastic, multiplicative way. Why does this occur?

Consider that the returned wave associated with a given pixel likely isn’t the result of a single interation with a single pointwise object. If multiple reflectors (or similarly, a continuous reflector) lie within the pixel, then the returned wave is the sum of multiple individual waves.

Let’s imagine we have a pixel with two identical point reflectors lying within it. How does the returned wave change as a function of where those reflectors are located? Move the blue and red targets in the pixel below:

Code
chart = {
    const svg = d3.create("svg").attr("viewBox", [0, 0, width, height])

    // define a bounding rectangle
    svg.append("rect")
        .attr("x", width/7)
        .attr("y", stroke_width)
        .attr("width", 5*width/7)
        .attr("height", height - 2*stroke_width)
        .attr("stroke", "black")
        .attr("stroke-width", stroke_width)
        .attr("fill-opacity", 0.087124976) //This is a hard-coded value based on the prescribed physical parameters. It's hardcoded because it's easier than dealing with the resulting reactive dependencies since this cell should never be run again

    // define data used for circles
    const circles = d3.range(2).map(i => ({
        x: (i+2) * width / 5,
        y: (i+1) * height / 3,
        index: i, 
    }));

    // define circles as graphic objects
    svg.selectAll("circle")
        .data(circles)
        .join("circle")
            .attr("cx", d => d.x)
            .attr("cy", d => d.y)
            .attr("r", radius)
            .attr("fill", d => d3.schemeCategory10[d.index * 3])
            .attr("id", function(d,i) {return i})
            .call(drag)

    return svg.node();
}
chart2 = {
    const svg = d3.create("svg").attr("viewBox", [0, 0, width, height])

    //define a function to take the canvas position of the leftside chart to the
    //  corresponding region in the second chart
    //  it's messy because I solved the system by hand and didn't simplify, RIP
    //  to future Jake if I ever have to un-hardcode the object position values
    function lineartransform(xposition){
        return (2*height/5 - radius)/(5/7*width -2*radius-2)*xposition + width/7 + 2/5*height + radius/2 - (2*height/5 - radius)*(width/7+radius+1)/(5/7*width-2*radius-2)
    }

    // define a earth curvature line
    svg.append("rect")
        .attr("x", width/7)
        .attr("y", height - 2*stroke_width)
        .attr("width", 5*width/7)
        .attr("height", stroke_width/2)
        .attr("stroke", "black")
        .attr("stroke-width", stroke_width)
        .attr("fill-opacity", 0)

    // define lines going from the target to the satellites (first, so they appear under them)
    svg.append("line")
        .attr("x1", width/7 + 3*height/20)
        .attr("x2", lineartransform(x1))
        .attr("y1", 3*height/20)
        .attr("y2", height - 2*stroke_width)
        .style("stroke-width", stroke_width)
        .style("stroke", "blue")

    svg.append("line")
        .attr("x1", width/7 + 3*height/20)
        .attr("x2", lineartransform(x2))
        .attr("y1", 3*height/20)
        .attr("y2", height - 2*stroke_width)
        .style("stroke-width", stroke_width)
        .style("stroke", "red")


    // define a rectangle representing the satellite
    svg.append("rect")
        .attr("x", width/7 + height/10)
        .attr("y", height/10)
        .attr("width", height/10)
        .attr("height", height/10)
        .attr("stroke", "black")
        .attr("stroke-width", stroke_width)
        .attr("fill-opacity", .50)
    
    // define notches indicating the bounds of the pixel
    svg.append("line")
        .attr("x1", width/7 + 2*height/5)
        .attr("x2", width/7 + 2*height/5)
        .attr("y1", height - 2*stroke_width)
        .attr("y2", height - 2*stroke_width - height/20)
        .style("stroke-width", stroke_width)
        .style("stroke", "black")

    svg.append("line")
        .attr("x1", width/7 + 4*height/5)
        .attr("x2", width/7 + 4*height/5)
        .attr("y1", height - 2*stroke_width)
        .attr("y2", height - 2*stroke_width - height/20)
        .style("stroke-width", stroke_width)
        .style("stroke", "black")

    // define data used for circles
    const circles = d3.range(2).map(i => ({
        x: (i+2) * width / 5,
        y: (i+1) * height / 3,
        index: i, 
    }));

    // define circles as graphic objects
    svg.selectAll("circle")
        .data(circles)
        .join("circle")
            .attr("cx", function(d,i) {
                if (i == 0) {return lineartransform(x1)}
                else {return lineartransform(x2)}
            })
            .attr("cy", height - 2*stroke_width)
            .attr("r", radius/2)
            .attr("fill", d => d3.schemeCategory10[d.index * 3])
            .attr("id", function(d,i) {return i})

    // add text
    svg.append("text")
        .attr("x", width/7)
        .attr("y", .94*height)
        .style("font-size", "28px")
        .text("side view")

    svg.append("text")
        .attr("x", width/2)
        .attr("y", .15*height)
        .style("font-size", "24px")
        .style("font-style", "italic")
        .text("*wildly not to scale")

    return svg.node();
}
Code
formatted_d1 = (0.5*distance_travelled_w1).toFixed(2).replace(/\B(?<!\.\d*)(?=(\d{3})+(?!\d))/g, ",");
formatted_d2 = (0.5*distance_travelled_w2).toFixed(2).replace(/\B(?<!\.\d*)(?=(\d{3})+(?!\d))/g, ",");

Distance from target one to the satellite: m
Distance from target two to the satellite: m

Code
// defines how circles change appearance and position in response to cursor drag
drag = {

    // move the circles to cursor position when dragging is active
    function dragged(event, d){
        var newx = event.x;
        if (newx < xbounds[0]) newx = xbounds[0];
        if (newx > xbounds[1]) newx = xbounds[1];
        var newy = event.y;
        if (newy < ybounds[0]) newy = ybounds[0];
        if (newy > ybounds[1]) newy = ybounds[1];
        
        if (d3.select(this).attr('id') == "0") {
            mutable x1 = newx;
            mutable y1 = newy;
        }
        if (d3.select(this).attr('id') == "1") {
            mutable x2 = newx;
            mutable y2 = newy;
        }
        
        d3.select(this).raise().attr("cx", d.x = newx).attr("cy", d.y = newy);

        var dt1 = 2*Math.sqrt(sat_z**2 + (mutable x1/width*pixel_width - sat_x)**2 + (mutable y1/width*pixel_width - sat_y)**2);
        var dt2 = 2*Math.sqrt(sat_z**2 + (mutable x2/width*pixel_width - sat_x)**2 + (mutable y2/width*pixel_width - sat_y)**2);
        var p1  =  dt1 / wavelength *2 *Math.PI;
        var p2  =  dt2 / wavelength *2 *Math.PI;

        d3.select("rect").attr("fill-opacity", 1-Math.abs(Math.cos((p1-p2)/2)));
        //d3.select("text").text(1-Math.abs(Math.cos((p1-p2)/2))).style("fill","darkOrange");
    }

    // map the d3 drag event functionality to these custom functions
    return d3.drag()
        .on("drag", dragged)
}
Code
mutable x1= 2*width/5;
mutable y1= height/3;
mutable x2= 3*width/5;
mutable y2= 2*height/3;
Code
stroke_width = 7;
width = 700;
height = 500 + 2*stroke_width;
radius = 25;
Code
xbounds = [width/7 + radius + 1, 6*width/7 - radius - 1]
ybounds = [ radius + 1  + stroke_width, height - radius - 1 - stroke_width]
Code
distance = 4;
NumPoints = 1000;
wavelength = .5;
amplitude = 15;
speed = 10;
Code
// display variables
xscale_factor = 100;
Code
sat_x = -100000;
sat_y = 100000;
sat_z = 600000;
pixel_width = 3;
Code
distance_travelled_w1 = 2*Math.sqrt(sat_z**2 + (x1/width*pixel_width - sat_x)**2 + (y1/width*pixel_width - sat_y)**2);
distance_travelled_w2 = 2*Math.sqrt(sat_z**2 + (x2/width*pixel_width - sat_x)**2 + (y2/width*pixel_width - sat_y)**2);
Code
phase_w1 = distance_travelled_w1 / wavelength *2 *Math.PI;
phase_w2 = distance_travelled_w2 / wavelength *2 *Math.PI;
Code
wave1 = {
    var data = [];
    for (var i = 1; i <= NumPoints; i++) {
        var j = i * distance / NumPoints;
        data.push([j*xscale_factor, amplitude * Math.sin(j * 2*Math.PI / wavelength + time + phase_w2)]);
    }
    return data;
}
Code
wave2 = {
    var data = [];
    for (var i = 1; i <= NumPoints; i++) {
        var j = i * distance / NumPoints;
        data.push([j*xscale_factor , amplitude * Math.sin(j * 2*Math.PI / wavelength + time + phase_w1)]);
    }
    return data;
}
Code
waveSum = {
    var data = [];
    for (var i = 1; i <= NumPoints; i++) {
        var j = i * distance / NumPoints;
        data.push([j*xscale_factor , amplitude * Math.sin(j * 2*Math.PI / wavelength + time + phase_w1) + amplitude * Math.sin(j * 2*Math.PI / wavelength + time + phase_w2)]);
    }
    return data;
}
Code
time = {
  let i = 0;
  while (true) {
    i += speed * 0.01;
    yield i
  }
}
Code
r = d3.line()(wave1);
b = d3.line()(wave2);
svg`<svg viewBox="0 -32 400 64">
  <path d="${r}" stroke="red" fill="none" />
  <path d="${b}" stroke="blue" fill="none" />
</svg>`
Code
p = d3.line()(waveSum);
svg`<svg viewBox="0 -64 400 130">
  <path d="${p}" stroke="black" fill="none" />
</svg>`
Code
d3 = require("d3@7")

Here we’ve considered a SAR satellite 600km above Earth, emitting 50cm wavelength radiation with a pixel resolution of 3m.

The distance from each target to the satellite affects where in their cycle each wave is re-intercepted, and thus the difference determines whether those returned waves interfere constructively or deconstructively. Thus a random distribution of targets creates a random difference in distance to the satellite, creating a random level of interference in each pixel. This means the intensity of each pixel is randomly (though not necessarily uniformly) multiplied by a value from 0 (perfect deconstructive interference) to 1 (perfect constructive interference). This is the origin of speckle, and reason why it manifests as multiplicative noise (as opposed to additive).

Let’s now consider two pixels: one with more, stronger scatterers and one with fewer, weaker scatterers – perhaps pixel one belongs to some dense shrubland, while the other is a freshly plowed field. The distribution of scatterers is likely quite random – what happens to the brightness of these pixels under different configurations?

Code
targetnum_1 = 10;
targetnum_2 = 10;
max_radius_p1 = 35;
min_radius_p1 = 25;
max_radius_p2 = 20;
min_radius_p2 = 15;

radii_p1 = {
    var data = [];
        for (var i = 1; i <= targetnum_1; i++) {
            data.push(Math.random()*(max_radius_p1-min_radius_p1) + min_radius_p1);
        }
    return data;
}

radii_p2 = {
    var data = [];
        for (var i = 1; i <= targetnum_2; i++) {
            data.push(Math.random()*(max_radius_p2-min_radius_p2) + min_radius_p2);
        }
    return data;
}
Code
viewof reroll_toggle = Inputs.form(
    [
        Inputs.button("Redistribute Targets"),
        Inputs.button("Auto Redistribute Targets"),
    ],
    {
        template: (inputs) => htl.html`<div style="display: flex; gap: 1em">
        ${inputs}
        </div>`
    })

vtoggle1    = reroll_toggle[0];
auto_toggle = reroll_toggle[1];
Code
time_auto = {
  let i = 0;
  while (auto_toggle%2==1) {
    i += 1;
    yield i
  }
}
Code
pixel1_positions = {
    var dummy  = vtoggle1;
    var dummy2 = time_auto;
    var data = [];
        for (var i = 1; i <= targetnum_1; i++) {
            var x = Math.random()*(5*width/7 - 2*max_radius_p1) + max_radius_p1 + width/7;
            var y = Math.random()*(height - 2*max_radius_p1) + max_radius_p1;
            data.push([x,y]);
        }
    return data
}

pixel2_positions = {
    var dummy  = vtoggle1;
    var dummy2 = time_auto;
    var data = [];
        for (var i = 1; i <= targetnum_2; i++) {
            var x = Math.random()*(5*width/7 - 2*max_radius_p1) + max_radius_p1 + width/7;
            var y = Math.random()*(height - 2*max_radius_p1) + max_radius_p1;
            data.push([x,y]);
        }
    return data
}
Code
//two pixels with targets in them

chart3 = {
    const svg = d3.create("svg").attr("viewBox", [0, 0, width, height])

    // define a bounding rectangle
    svg.append("rect")
        .attr("x", width/7)
        .attr("y", stroke_width)
        .attr("width", 5*width/7)
        .attr("height", height - 2*stroke_width)
        .attr("stroke", "black")
        .attr("stroke-width", stroke_width)
        .attr("fill-opacity", 1-amplitude_p1/200)

    // define data used for circles
    const circles = d3.range(targetnum_1).map(i => ({
        x: pixel1_positions[i][0],
        y: pixel1_positions[i][1],
        index: i, 
    }));

    // define circles as graphic objects
    svg.selectAll("circle")
        .data(circles)
        .join("circle")
            .attr("cx", d => d.x)
            .attr("cy", d => d.y)
            .attr("r", function(d,i) {return radii_p1[i]})
            .attr("stroke", "black")
            .attr("stroke-width", 4)
            .attr("fill", d => d3.schemeCategory10[0])
            .attr("id", function(d,i) {return i})

    return svg.node();
}
chart4 = {
    const svg = d3.create("svg").attr("viewBox", [0, 0, width, height])

    // define a bounding rectangle
    svg.append("rect")
        .attr("x", width/7)
        .attr("y", stroke_width)
        .attr("width", 5*width/7)
        .attr("height", height - 2*stroke_width)
        .attr("stroke", "black")
        .attr("stroke-width", stroke_width)
        .attr("fill-opacity", 1-amplitude_p2/200)

    // define data used for circles
    const circles = d3.range(targetnum_2).map(i => ({
        x: pixel2_positions[i][0],
        y: pixel2_positions[i][1],
        index: i, 
    }));

    // define circles as graphic objects
    svg.selectAll("circle")
        .data(circles)
        .join("circle")
            .attr("cx", d => d.x)
            .attr("cy", d => d.y)
            .attr("r", function(d,i) {return radii_p2[i]})
            .attr("stroke", "black")
            .attr("stroke-width", 4)
            .attr("fill", d => d3.schemeCategory10[0])
            .attr("id", function(d,i) {return i})

    return svg.node();
}
Code
Plot.plot({
    round: true,
    color: {legend: true},
    x: {label: "Brightness"},
    y: {label: "Frequency"},
    marks: [
        Plot.rectY(brightnesses, Plot.binX({y2: "count"}, {x: "x", fill: "p", mixBlendMode: "multiply"})),
        Plot.ruleY([0])
        ]
})
Code
amplitude_p1 = {

    var dummy = vtoggle1;
    
    var amplitude_naught = radii_p1[0];

    var dt_naught =
        2*Math.sqrt(
            sat_z**2
            + (pixel1_positions[0][0]/width*pixel_width - sat_x)**2
            + (pixel1_positions[0][1]/width*pixel_width - sat_y)**2
        );
    
    var phase_naught =
        dt_naught / wavelength *2 *Math.PI;
    
    //successively add the waves, determining the new phases and amplitudes
    for (var i = 1; i <= targetnum_1-1; i++) {
        //compute the phase of the next subwave
        var dt_next =
            2*Math.sqrt(
                sat_z**2
                + (pixel1_positions[i][0]/width*pixel_width - sat_x)**2
                + (pixel1_positions[i][1]/width*pixel_width - sat_y)**2
        );

        var phase_next = 
            dt_next / wavelength *2 *Math.PI;
        
        //compute the amplitude of the new sum
        var new_amplitude =
            Math.sqrt(
                amplitude_naught**2
                + radii_p1[i]**2
                + 2*amplitude_naught*radii_p1[i]*Math.cos(phase_next-phase_naught)
            );

        //compute the phase of the new sum
        var new_phase =
            Math.atan(
                (amplitude_naught*Math.sin(phase_naught) + radii_p1[i]*Math.sin(phase_next))
                /(amplitude_naught*Math.cos(phase_naught) + radii_p1[i]*Math.cos(phase_next))
            )

        amplitude_naught = new_amplitude;
        phase_naught = new_phase;

    }
    return amplitude_naught
}

amplitude_p2 = {

    var dummy = vtoggle1;
    
    var amplitude_naught = radii_p2[0];

    var dt_naught =
        2*Math.sqrt(
            sat_z**2
            + (pixel2_positions[0][0]/width*pixel_width - sat_x)**2
            + (pixel2_positions[0][1]/width*pixel_width - sat_y)**2
        );
    
    var phase_naught =
        dt_naught / wavelength *2 *Math.PI;
    
    //successively add the waves, determining the new phases and amplitudes
    for (var i = 1; i <= targetnum_1-1; i++) {
        //compute the phase of the next subwave
        var dt_next =
            2*Math.sqrt(
                sat_z**2
                + (pixel2_positions[i][0]/width*pixel_width - sat_x)**2
                + (pixel2_positions[i][1]/width*pixel_width - sat_y)**2
        );

        var phase_next = 
            dt_next / wavelength *2 *Math.PI;
        
        //compute the amplitude of the new sum
        var new_amplitude =
            Math.sqrt(
                amplitude_naught**2
                + radii_p2[i]**2
                + 2*amplitude_naught*radii_p2[i]*Math.cos(phase_next-phase_naught)
            );

        //compute the phase of the new sum
        var new_phase =
            Math.atan(
                (amplitude_naught*Math.sin(phase_naught) + radii_p2[i]*Math.sin(phase_next))
                /(amplitude_naught*Math.cos(phase_naught) + radii_p2[i]*Math.cos(phase_next))
            )

        amplitude_naught = new_amplitude;
        phase_naught = new_phase;

    }
    return amplitude_naught
}

{
mutable brightnesses.push({x: amplitude_p1, p: "Pixel 1"});
mutable brightnesses.push({x: amplitude_p2, p: "Pixel 2"});
mutable brightnesses = mutable brightnesses;
}
Code
//histogram to be populated with the brightness values of the two pixels
mutable brightnesses = [];

Clicking through a few dozen random configurations of each pixel, we can see in the developing histogram of their brightnesses that while the pixel with more strongly reflecting targets is on average brighter, there are many times when the pixel with weaker targets actually gives a higher return. This is due to the random level of interference we’ve called speckle. We can exploit the fact that the average value remains higher to identify likely high-return areas in speckle filtering. As an aside, notice that the distributions of brightness appear relatively Gaussian – this is a lovely invocation of the Central Limit Theorem.


© 2023 Jake Ferguson

 

This website attempts to be WAI-ARIA compliant – if you find any issues please email me
Built using Quarto; CSS theme modified from litera