https://www.shadertoy.com/view/MsjfRK

bufA.glsl

// Copyright (c) 2017 - Élie Michel
// License: CC BY 3.0
// Please credit even if you are just re-using the particle system.

// Sync with 'Image' buffer
const float dt = 0.001;

mat3 particleSpaceMatrix(vec2 origin, vec2 velocity) {
    vec3 O = vec3(origin, 1.);
    vec3 X = normalize(vec3(velocity, 0.));
    vec3 Y = cross(vec3(0., 0., 1.), X);
    return mat3(X, Y, O);
}
const float celerity = 1200.0;

vec4 data(int addr) {
    int h = int(iResolution.y);
    ivec2 v = ivec2(addr / h, addr % h);
    return textureLod(iChannel0, (vec2(v) + vec2(0.5, 0.5))/iResolution.xy, 0.0);
}

struct Header {
    vec2 mouse;
    vec2 oldMouse;
    int nbParticles;
};
const int sizeofHeader = 2;
Header getHeader() {
    vec4 d0 = data(0);
    vec4 d1 = data(1);
    return Header(d0.xy, d0.zw, int(d1.x));
}

struct Particle {
    vec4 color;
    vec2 position;
    vec2 velocity;
    float wavelength;
    float phi;
    float lastHit;  // Time since last collision
};
const int sizeofParticle = 4;
const int colorField = 0;
const int posVelField = 1;
const int wavelengthField = 2;
Particle getParticle(int index) {
    int offset = sizeofHeader + index * sizeofParticle;
    vec4 color = data(offset + colorField);
    vec4 posVel = data(offset + posVelField);
    vec4 wavelengthPhi = data(offset + wavelengthField);
    return Particle(color, posVel.xy, posVel.zw, wavelengthPhi.x, wavelengthPhi.y, wavelengthPhi.z);
}
vec2 particleHeadPosition(Particle particle) {
    mat3 M = particleSpaceMatrix(particle.position, particle.velocity);
    float omega = 2.0 * 3.1415 * celerity / particle.wavelength;
    float l = sin(omega * iTime + particle.phi) * 0.02;
    return M[2].xy - l * M[1].xy;
}
//

// Inspired by http://www.noah.org/wiki/Wavelength_to_RGB_in_Python
vec3 colorFromWavelength(float wavelength) {
    const float gamma = 0.8;
    float r, g, b;
    if (wavelength >= 380. && wavelength <= 440.) {
        float attenuation = .3 + .7 * (wavelength - 380.) / (440. - 380.);
        r = pow((-(wavelength - 440.) / (440. - 380.)) * attenuation, gamma);
        g = 0.;
        b = pow(1.0 * attenuation, gamma);
    } else if (wavelength >= 440. && wavelength <= 490.) {
        r = 0.;
        g = pow((wavelength - 440.) / (490. - 440.), gamma);
        b = 1.;
    } else if (wavelength >= 490. && wavelength <= 510.) {
        r = 0.;
        g = 1.;
        b = pow(-(wavelength - 510.) / (510. - 490.), gamma);
    } else if (wavelength >= 510. && wavelength <= 580.) {
        r = pow((wavelength - 510.) / (580. - 510.), gamma);
        g = 1.;
        b = 0.;
    } else if (wavelength >= 580. && wavelength <= 645.) {
        r = 1.;
        g = pow(-(wavelength - 645.) / (645. - 580.), gamma);
        b = 0.;
    } else if (wavelength >= 645. && wavelength <= 750.) {
        float attenuation = .3 + .7 * (750. - wavelength) / (750. - 645.);
        r = pow(1. * attenuation, gamma);
        g = 0.;
        b = 0.;
    } else {
        r = 1.;
        g = 1.;
        b = 1.;
    }
    return vec3(r, g, b);
}

void writeHeader(out vec4 fragColor, in int addr) {
    switch (addr) {
    case 0:
    {
        Header header = getHeader();
        fragColor = vec4(iMouse.xy, header.mouse);
        break;
    }
    case 1:
        fragColor = vec4(30.0);
        break;
    }
}

void updateParticle(inout Particle particle, int particleId) {
    bool init = iTime < 0.1;
    //init = true;
    
    vec2 velocity = particle.velocity;
    vec2 position = particle.position + dt * velocity;
    vec2 headPosition = particleHeadPosition(particle);
    float wavelength = particle.wavelength;
    float phi = particle.phi;
    float lastHit = particle.lastHit;

    // Out-of-bounds reset
    if ((position.x < 0.0 && velocity.x < 0.0)
        || (position.x > iResolution.x/iResolution.y && velocity.x > 0.0)
        || (position.y < 0.0 && velocity.y < 0.0)
        || (position.y > 1.0 && velocity.y > 1.0)) {
        init = true;
    }
    if (init) {
        velocity = vec2(-1.0, -0.5) * 4.0;
        position = vec2(0.9 + float(particleId) * 0.1, 1.1);
        wavelength = 700.0;  // nm
        phi = float(particleId);
        lastHit = 9999.0;
        
        wavelength = 450.0 + 90.0 * float(particleId % 3);
        float theta = 3.1415 * (0.15 + 0.1 * float(particleId % 3));
        vec2 offset = (0.1 - 0.005 * float(particleId)) * vec2(-sin(theta), cos(theta));
        position = vec2(1.0, 1.0) + 0.5 * vec2(cos(theta), sin(theta)) + offset;
        velocity = -vec2(cos(theta), sin(theta)) * celerity / 300.0;
    }
    
    if (lastHit >= 0.) {
        lastHit += 1.0;
    }
    
    // Collision detection
    if (!init && lastHit > 20.0 && headPosition.x < iMouse.x/iResolution.y && velocity.x <= 0.0) {
        Header header = getHeader();
        vec2 deltaMouse = (header.oldMouse.xy - iMouse.xy)/iResolution.y;
        
        position = headPosition + dt * velocity;
        position.x = iMouse.x/iResolution.y;
        velocity.x = abs(velocity.x);
        //velocity *= (1.0 - deltaMouse.x * 1.0);
        //velocity.x -= deltaMouse.x * 10.0;

        
        //position = particle.position + dxy;
        
        wavelength = clamp(wavelength + deltaMouse.x * 2000.0, 380., 750.);
        //wavelength = 9999.0;
        float omega = 2.0 * 3.1415 * celerity / wavelength;
        phi = -omega * iTime;
        //phi += 3.1415;
        lastHit = 0.0;
    }
    
    // Write
    //particle.color = vec4(1.0, 0.5, float(particleId) * 0.2, 1.0);
    particle.color = vec4(colorFromWavelength(wavelength), 1.0);
    particle.position = position;
    particle.velocity = velocity;
    particle.wavelength = wavelength;
    particle.phi = phi;
    particle.lastHit = lastHit;
}

void writeParticle(out vec4 fragColor, in int addr) {
    // Frag coord in particle/field space
    int field = addr % sizeofParticle;
    int particleId = addr / sizeofParticle;
    
    // Load particle info
    Particle particle = getParticle(particleId);
        
    // Update all fields (It's required to update all fields even though we
    // only write some of them in a given fragment because for instance
    // wavelength needs to be aware of the changes of velocity/position.)
    updateParticle(particle, particleId);
    
    switch (field) {
    case colorField:
        fragColor = particle.color;
        break;
    case posVelField:
        fragColor = vec4(particle.position, particle.velocity);
        break;
    case wavelengthField:
        fragColor = vec4(particle.wavelength, particle.phi, particle.lastHit, 0.0);
        break;
    }
}

void mainImage( out vec4 fragColor, in vec2 fragCoord )
{
    ivec2 st = ivec2(fragCoord);
    int addr = st.x * int(iResolution.y) + st.y;
    
    if (addr < sizeofHeader) {
        writeHeader(fragColor, addr);
    } else {
        addr -= sizeofHeader;
        writeParticle(fragColor, addr);
    }
}

bufB.glsl

// Copyright (c) 2017 - Élie Michel
// License: CC BY 3.0
// Please credit even if you are just re-using the particle system.

// Sync with 'Image' buffer
mat3 particleSpaceMatrix(vec2 origin, vec2 velocity) {
    vec3 O = vec3(origin, 1.);
    vec3 X = normalize(vec3(velocity, 0.));
    vec3 Y = cross(vec3(0., 0., 1.), X);
    return mat3(X, Y, O);
}
const float celerity = 1200.0;

vec4 data(int addr) {
    int h = int(iResolution.y);
    ivec2 v = ivec2(addr / h, addr % h);
    return textureLod(iChannel0, (vec2(v) + vec2(0.5, 0.5))/iResolution.xy, 0.0);
}

struct Header {
    vec2 mouse;
    vec2 oldMouse;
    int nbParticles;
};
const int sizeofHeader = 2;
Header getHeader() {
    vec4 d0 = data(0);
    vec4 d1 = data(1);
    return Header(d0.xy, d0.zw, int(d1.x));
}

struct Particle {
    vec4 color;
    vec2 position;
    vec2 velocity;
    float wavelength;
    float phi;
    float lastHit;  // Time since last collision
};
const int sizeofParticle = 4;
const int colorField = 0;
const int posVelField = 1;
const int wavelengthField = 2;
Particle getParticle(int index) {
    int offset = sizeofHeader + index * sizeofParticle;
    vec4 color = data(offset + colorField);
    vec4 posVel = data(offset + posVelField);
    vec4 wavelengthPhi = data(offset + wavelengthField);
    return Particle(color, posVel.xy, posVel.zw, wavelengthPhi.x, wavelengthPhi.y, wavelengthPhi.z);
}
vec2 particleHeadPosition(Particle particle) {
    mat3 M = particleSpaceMatrix(particle.position, particle.velocity);
    float omega = 2.0 * 3.1415 * celerity / particle.wavelength;
    float l = sin(omega * iTime + particle.phi) * 0.02;
    return M[2].xy - l * M[1].xy;
}
//

void mainImage( out vec4 fragColor, in vec2 fragCoord )
{
    vec2 uv0 = fragCoord.xy / iResolution.xy;
    vec2 uv = fragCoord.xy / iResolution.y;
    fragColor = texture(iChannel1, uv0) * 0.98;
    
    Header header = getHeader();
    for (int i = 0 ; i < header.nbParticles ; ++i) {
        Particle particle = getParticle(i);
        vec2 p = particleHeadPosition(particle);
        fragColor += mix(particle.color, vec4(0.0), min(1.0, length(uv - p) * 200.0));
    }
}

doppler.glsl

// Copyright (c) 2017 - Élie Michel
// License: CC BY 3.0
// Please credit even if you are just re-using the particle system.

// Sync with 'Buf A' and 'Buf B' buffers
mat3 particleSpaceMatrix(vec2 origin, vec2 velocity) {
    vec3 O = vec3(origin, 1.);
    vec3 X = normalize(vec3(velocity, 0.));
    vec3 Y = cross(vec3(0., 0., 1.), X);
    return mat3(X, Y, O);
}
const float celerity = 1200.0;

vec4 data(int addr) {
    int h = int(iResolution.y);
    ivec2 v = ivec2(addr / h, addr % h);
    return textureLod(iChannel0, (vec2(v) + vec2(0.5, 0.5))/iResolution.xy, 0.0);
}

struct Header {
    vec2 mouse;
    vec2 oldMouse;
    int nbParticles;
};
const int sizeofHeader = 2;
Header getHeader() {
    vec4 d0 = data(0);
    vec4 d1 = data(1);
    return Header(d0.xy, d0.zw, int(d1.x));
}

struct Particle {
    vec4 color;
    vec2 position;
    vec2 velocity;
    float wavelength;
    float phi;
    float lastHit;  // Time since last collision
};
const int sizeofParticle = 4;
const int colorField = 0;
const int posVelField = 1;
const int wavelengthField = 2;
Particle getParticle(int index) {
    int offset = sizeofHeader + index * sizeofParticle;
    vec4 color = data(offset + colorField);
    vec4 posVel = data(offset + posVelField);
    vec4 wavelengthPhi = data(offset + wavelengthField);
    return Particle(color, posVel.xy, posVel.zw, wavelengthPhi.x, wavelengthPhi.y, wavelengthPhi.z);
}
vec2 particleHeadPosition(Particle particle) {
    mat3 M = particleSpaceMatrix(particle.position, particle.velocity);
    float omega = 2.0 * 3.1415 * celerity / particle.wavelength;
    float l = sin(omega * iTime + particle.phi) * 0.02;
    return M[2].xy - l * M[1].xy;
}
//

void plotAxis(inout vec4 fragColor, vec2 basis) {
    float ax;
    ax = smoothstep(0., 1., abs(basis.y) * 1000.);
    fragColor = mix(vec4(1., 0., 0., 1.), fragColor, ax);
    ax = smoothstep(0., 1., abs(basis.x) * 1000.);
    fragColor = mix(vec4(0., 1., 0., 1.), fragColor, ax);
}

const float period = 2.5;
const float omega = 100.0;
const float A = 0.01;

void mainImage( out vec4 fragColor, in vec2 fragCoord )
{
    vec2 uv0 = fragCoord.xy / iResolution.xy;
	vec2 uv = fragCoord.xy / iResolution.y;
	fragColor = vec4(uv,0.5+0.5*sin(iTime),1.0);
    fragColor = vec4(0.0);
    
    vec2 xy;

    float t = iTime * 0.3;
    float st = t - round(t/period)*period;
    float s0 = st - 0.5, s1 = st;
    float s = xy.x;
    
    Header header = getHeader();
    for (int i = 0 ; i < header.nbParticles ; ++i) {
        Particle particle = getParticle(i);
        vec2 p = particleHeadPosition(particle);
        vec4 c = mix(vec4(1.0), particle.color, length(uv - p) * 100.0);
        fragColor += mix(c, vec4(0.0), min(1.0, length(uv - p) * 100.0));
    }
    
    fragColor += texture(iChannel1, uv0);
    
    //plotAxis(fragColor, xy);
    plotAxis(fragColor, uv-iMouse.xy/iResolution.y);
}