Do not have too much time for this right now but you can start with this:
blur the image a bit (remove noise)
simple convolution will do for example few times with matrix
0.0 0.1 0.0
0.1 0.6 0.1
0.0 0.1 0.0
do color derivation by y axis
derivate pixel color along y
axis... for example I used this for each pixel of input image:
pixel(x,y)=pixel(x,y)-pixel(x,y-1)
beware the result is signed so you can normalize by some bias or use abs
value or handle as signed values ... In my example I used bias so the gray area is zero derivation ... black is most negative and white is most positive
blur the image a bit (remove noise)
enhance dynamic range
Simply find in image the min color c0
and max color c1
and rescale all pixels to predefined range <low,high>
. This will make the tresholding much more stable across different images...
pixel(x,y)=low + ((pixel(x,y)-c0)*(high-low)/(c1-c0)
so for example you can use low=0
and high=255
treshold all pixels which are bigger then treshold
The resulting image is like this:
Now just:
- segment the red areas
- remove too small areas
shrink/recolor areas to leave just the midpoint on each x
coordinate
The top most point is red and the bottom is green.
This should lead you to very close state to your wanted solution. Beware blurring and derivation can move the detected positions a bit from their real position.
Also for more ideas take a look at related QAs:
[Edit1] C++ code of mine for this
picture pic0,pic1,pic2; // pic0 is your input image, pic1,pic2 are output images
int x,y;
int tr0=Form1->sb_treshold0->Position; // treshold from scrollbar (=100)
// [prepare image]
pic1=pic0; // copy input image pic0 to output pic1 (upper)
pic1.pixel_format(_pf_s); // convert to grayscale intensity <0,765> (handle as signed numbers)
pic2=pic0; // copy input image pic0 to output pic2 (lower)
pic1.smooth(5); // blur 5 times
pic1.derivey(); // derive colros by y
pic1.smooth(5); // blur 5 times
pic1.enhance_range(); // maximize range
for (x=0;x<pic1.xs;x++) // loop all pixels
for (y=0;y<pic1.ys;y++)
if (pic1.p[y][x].i>=tr0) // if treshold in pic1 condition met
pic2.p[y][x].dd=0x00FF0000; //0x00RRGGBB then recolor pixel in pic2
pic1.pixel_format(_pf_rgba); // convert the derivation signed grayscale to RGBA (as biased grayscale)
// just render actual set treshold
pic2.bmp->Canvas->Brush->Style=bsClear;
pic2.bmp->Canvas->Font->Color=clYellow;
pic2.bmp->Canvas->TextOutA(5,5,AnsiString().sprintf("Treshold: %i",tr0));
pic2.bmp->Canvas->Brush->Style=bsSolid;
The code is using Borlands VCL encapsulated GDI bitmap/canvas at bottom (not important for you just renders actual treshold settings) and my own picture
class so some members description is in order:
xs,ys
resolution
color p[ys][xs]
direct pixel access (32bit pixel format so 8 bit per channel)
pf
actual selected pixel format for the image see the enum
bellow
enc_color/dec_color
just pack unpack color channels to separate array for easy multi pixel format handling ... so I do not need to write each function for each pixelformat separately
clear(DWORD c)
fills image with color c
The color
is just union
of DWORD dd
and BYTE db[4]
and int i
for simple channel access and or signed values handling.
Some chunks of code from it:
union color
{
DWORD dd; WORD dw[2]; byte db[4];
int i; short int ii[2];
color(){}; color(color& a){ *this=a; }; ~color(){}; color* operator = (const color *a) { dd=a->dd; return this; }; /*color* operator = (const color &a) { ...copy... return this; };*/
};
enum _pixel_format_enum
{
_pf_none=0, // undefined
_pf_rgba, // 32 bit RGBA
_pf_s, // 32 bit signed int
_pf_u, // 32 bit unsigned int
_pf_ss, // 2x16 bit signed int
_pf_uu, // 2x16 bit unsigned int
_pixel_format_enum_end
};
//---------------------------------------------------------------------------
void dec_color(int *p,color &c,int _pf)
{
p[0]=0;
p[1]=0;
p[2]=0;
p[3]=0;
if (_pf==_pf_rgba) // 32 bit RGBA
{
p[0]=c.db[0];
p[1]=c.db[1];
p[2]=c.db[2];
p[3]=c.db[3];
}
else if (_pf==_pf_s ) // 32 bit signed int
{
p[0]=c.i;
}
else if (_pf==_pf_u ) // 32 bit unsigned int
{
p[0]=c.dd;
}
else if (_pf==_pf_ss ) // 2x16 bit signed int
{
p[0]=c.ii[0];
p[1]=c.ii[1];
}
else if (_pf==_pf_uu ) // 2x16 bit unsigned int
{
p[0]=c.dw[0];
p[1]=c.dw[1];
}
}
//---------------------------------------------------------------------------
void dec_color(double *p,color &c,int _pf)
{
p[0]=0.0;
p[1]=0.0;
p[2]=0.0;
p[3]=0.0;
if (_pf==_pf_rgba) // 32 bit RGBA
{
p[0]=c.db[0];
p[1]=c.db[1];
p[2]=c.db[2];
p[3]=c.db[3];
}
else if (_pf==_pf_s ) // 32 bit signed int
{
p[0]=c.i;
}
else if (_pf==_pf_u ) // 32 bit unsigned int
{
p[0]=c.dd;
}
else if (_pf==_pf_ss ) // 2x16 bit signed int
{
p[0]=c.ii[0];
p[1]=c.ii[1];
}
else if (_pf==_pf_uu ) // 2x16 bit unsigned int
{
p[0]=c.dw[0];
p[1]=c.dw[1];
}
}
//---------------------------------------------------------------------------
void enc_color(int *p,color &c,int _pf)
{
c.dd=0;
if (_pf==_pf_rgba) // 32 bit RGBA
{
c.db[0]=p[0];
c.db[1]=p[1];
c.db[2]=p[2];
c.db[3]=p[3];
}
else if (_pf==_pf_s ) // 32 bit signed int
{
c.i=p[0];
}
else if (_pf==_pf_u ) // 32 bit unsigned int
{
c.dd=p[0];
}
else if (_pf==_pf_ss ) // 2x16 bit signed int
{
c.ii[0]=p[0];
c.ii[1]=p[1];
}
else if (_pf==_pf_uu ) // 2x16 bit unsigned int
{
c.dw[0]=p[0];
c.dw[1]=p[1];
}
}
//---------------------------------------------------------------------------
void enc_color(double *p,color &c,int _pf)
{
c.dd=0;
if (_pf==_pf_rgba) // 32 bit RGBA
{
c.db[0]=p[0];
c.db[1]=p[1];
c.db[2]=p[2];
c.db[3]=p[3];
}
else if (_pf==_pf_s ) // 32 bit signed int
{
c.i=p[0];
}
else if (_pf==_pf_u ) // 32 bit unsigned int
{
c.dd=p[0];
}
else if (_pf==_pf_ss ) // 2x16 bit signed int
{
c.ii[0]=p[0];
c.ii[1]=p[1];
}
else if (_pf==_pf_uu ) // 2x16 bit unsigned int
{
c.dw[0]=p[0];
c.dw[1]=p[1];
}
}
//---------------------------------------------------------------------------
void picture::smooth(int n)
{
color *q0,*q1;
int x,y,i,c0[4],c1[4],c2[4];
bool _signed;
if ((xs<2)||(ys<2)) return;
for (;n>0;n--)
{
#define loop_beg for (y=0;y<ys-1;y++){ q0=p[y]; q1=p[y+1]; for (x=0;x<xs-1;x++) { dec_color(c0,q0[x],pf); dec_color(c1,q0[x+1],pf); dec_color(c2,q1[x],pf);
#define loop_end enc_color(c0,q0[x ],pf); }}
if (pf==_pf_rgba) loop_beg for (i=0;i<4;i++) { c0[i]=(c0[i]+c0[i]+c1[i]+c2[i])>>2; clamp_u8(c0[i]); } loop_end
if (pf==_pf_s ) loop_beg { c0[0]=(c0[0]+c0[0]+c1[0]+c2[0])/ 4; clamp_s32(c0[0]); } loop_end
if (pf==_pf_u ) loop_beg { c0[0]=(c0[0]+c0[0]+c1[0]+c2[0])>>2; clamp_u32(c0[0]); } loop_end
if (pf==_pf_ss ) loop_beg for (i=0;i<2;i++) { c0[i]=(c0[i]+c0[i]+c1[i]+c2[i])/ 4; clamp_s16(c0[i]); } loop_end
if (pf==_pf_uu ) loop_beg for (i=0;i<2;i++) { c0[i]=(c0[i]+c0[i]+c1[i]+c2[i])>>2; clamp_u16(c0[i]); } loop_end
#undef loop_beg
#define loop_beg for (y=ys-1;y>0;y--){ q0=p[y]; q1=p[y-1]; for (x=xs-1;x>0;x--) { dec_color(c0,q0[x],pf); dec_color(c1,q0[x-1],pf); dec_color(c2,q1[x],pf);
if (pf==_pf_rgba) loop_beg for (i=0;i<4;i++) { c0[i]=(c0[i]+c0[i]+c1[i]+c2[i])>>2; clamp_u8(c0[i]); } loop_end
if (pf==_pf_s ) loop_beg { c0[0]=(c0[0]+c0[0]+c1[0]+c2[0])/ 4; clamp_s32(c0[0]); } loop_end
if (pf==_pf_u ) loop_beg { c0[0]=(c0[0]+c0[0]+c1[0]+c2[0])>>2; clamp_u32(c0[0]); } loop_end
if (pf==_pf_ss ) loop_beg for (i=0;i<2;i++) { c0[i]=(c0[i]+c0[i]+c1[i]+c2[i])/ 4; clamp_s16(c0[i]); } loop_end
if (pf==_pf_uu ) loop_beg for (i=0;i<2;i++) { c0[i]=(c0[i]+c0[i]+c1[i]+c2[i])>>2; clamp_u16(c0[i]); } loop_end
#undef loop_beg
#undef loop_end
}
}
//---------------------------------------------------------------------------
void picture::enhance_range()
{
int i,x,y,a0[4],min[4],max,n,c0,c1,q,c;
if (xs<1) return;
if (ys<1) return;
n=0; // dimensions to interpolate
if (pf==_pf_s ) { n=1; c0=0; c1=127*3; }
if (pf==_pf_u ) { n=1; c0=0; c1=255*3; }
if (pf==_pf_ss ) { n=2; c0=0; c1=32767; }
if (pf==_pf_uu ) { n=2; c0=0; c1=65535; }
if (pf==_pf_rgba) { n=4; c0=0; c1= 255; }
// find min,max
dec_color(a0,p[0][0],pf);
for (i=0;i<n;i++) min[i]=a0[i]; max=0;
for (y=0;y<ys;y++)
for (x=0;x<xs;x++)
{
dec_color(a0,p[y][x],pf);
for (q=0,i=0;i<n;i++)
{
c=a0[i]; if (c<0) c=-c;
if (min[i]>c) min[i]=c;
if (max<c) max=c;
}
}
// change dynamic range to max
for (y=0;y<ys;y++)
for (x=0;x<xs;x++)
{
dec_color(a0,p[y][x],pf);
for (i=0;i<n;i++) a0[i]=c0+(((a0[i]-min[i])*(c1-c0))/(max-min[i]+1));
// for (i=0;i<n;i++) if (a0[i]<c0) a0[i]=c0; // clamp if needed
// for (i=0;i<n;i++) if (a0[i]>c1) a0[i]=c1; // clamp if needed
enc_color(a0,p[y][x],pf);
}
}
//---------------------------------------------------------------------------
void picture::derivey()
{
int i,x,y,a0[4],a1[4];
if (ys<2) retur