cross-corrleation between two images

[gwendolien]

I'm trying to use two sources of the same movie.
Both exhibit funny parts where the frames seem to be shifted a couple of lines down.
Trying tpo stabilize this with DePan was unsuccesfull.

mdata = DePanEstimate(avi)
DePanStabilize(avi, data=mdata)

crosscorrelation between the corresponding frames could be very helpfull in mutually aligning them properly...
Would this be possible to do in avisynth and if so, how?

[joka]

DePanStabilize smoothes the global motion, but does not compensate it completely. A time ago I scipted a function for a similar use. May be this will work for you. You need MVTool2 in additional.

Function Align (clip c, clip r) {

t = interleave(r, c)
super = t.MSuper(pel=2)
v = super.MAnalyse(isb=true, delta=1)
mdata = t.MDePAn(v)

Return DePanInterleave(t, data=mdata).SelectEvery(6,2) }

[gwendolien]

Thank you very much for the response, Joka!
At present/next 12 days, I am not at my home-PC, but as soon as I have returned, I`m going to try your suggestion!
But I`m not very experienced with all avisybth clus and tricks, so...

[gwendolien]

I have tried your suggestion but helas, this gave only partial result:
the start of the vertical shift in movie2 coincides with a scene change.
At the END of this scene, I see attempts to correct the vertical shift, but the correction is incomplete and varies with frames.
The major first part of the scene is NOT changed at all.
Furthermore: closer examination of the two movies shows that there are also HORIZONTAL displacements every now and then.
So I'm afraid I have to do a real crosscorrelation in 2 dimensions to get the proper shift values and correct movie2 to movie1.

Thanks for your suggestion.

[StainlessS]

From RT_Stats:
http://forum.doom9.org/showthread.php?t=165479&highlight=rt_stats


RT_LumaDifference(clip,clip2,int "n"=current_frame,int "delta"=0,int "x"=0,int "y"=0,int "w"=0,int "h"=0,
int "n2"=current_frame,int "delta2"=0,int "x2"=x,int "y2"=y,bool "interlaced"=false,int "Matrix"=(width<=720?2:3))
Returns FLOAT value luma difference (0.0 -> 255.0) between clip frame (n+delta) area x,y,w,h, and clip2 frame (n2+delta2) area x2,y2,w,h.
Note, 'x2' and 'y2' default to 'x' and 'y' respectively.

***
***
***

RT_LumaCorrelation(clip,clip2,int "n"=current_frame,int "delta"=0,int "x"=0,int "y"=0,int "w"=0,int "h"=0,int "n2"=current_frame,
int "delta2"=0,int "x2"=x,int "y2"=y,bool "interlaced"=false,int "Matrix"=(width<=720?2:3))
Returns FLOAT value luma correlation (-1.0 -> 1.0) between clip frame (n+delta) area x,y,w,h, and clip2 frame (n2+delta2) area x2,y2,w,h.
Note, 'x2' and 'y2' default to 'x' and 'y' respectively.
Pearson's sample correlation coefficient.
http://en.wikipedia.org/wiki/Correlation_and_dependence
http://en.wikipedia.org/wiki/Pearson_product-moment_correlation_coefficient
Uses equivalent routine to the JMac698's JCorr plugin here:
http://forum.doom9.org/showthread.php?t=165386
and here:
http://forum.doom9.org/showthread.php?p=1495098#post1495098


Also Jmac698's JCorr as given in code block above

and Corr2D by V.C. Mohan
http://avisynth.org/vcmohan/

[Terka]

StainlessS, would be possible implement existing cross-corrleation solution to mvtools?

[StainlessS]

Not by me, but if any of the developers of MvTools want to take any of the the RT_Stats source, they are welcome.
It's not a very fast function (although implemented as fast as I could which made it have to deal with much bigger
numbers, there was an easier method but that would be considerably slower). However, you could just use a smaller
window or perhaps even strips, say half height and maybe a dozen or so pixels wide, or a horizontal strip, that could
I think shave some time off doing a full/near full image comparison. Perhaps use a left and right vertical strip and compare
results to see if the agree, same for top and bottom.
RT_LumaCorrelation() source is self contained and easily ripped out of the source so could relatively easily be used by another
plug (although there are a couple of helper functions that is uses to throw an error or get a global var, could just take them too).
I did try to use RT_LumaCorrelation plug on the PlanetCrop script results to remove the small amount of jitter of the planets, was not
really any use at all, seems that the PlanetCrop script does pretty much as good a job as is likely already. But that was a quite
difficult task, there were lots of anomalies in the Jmac Moon2 sequence that made it difficult to chop out the planets, including
background starlight (levels reveals star clusters moving with the moon as earth rotates) and clouds going over the face of the moon
(clearly visible if you play at eg 10/15FPS). Also, not one but two halos around the moon, and a constantly moving sun (in relation to the
moon) playing tricks with the outline.
I did I think, at one point have a look at MvTools source and did not find it easy to read or understand, the developers would be better people
to add it in if they thought it was of value.
The source is I think pretty well commented and should not be too difficult to follow, also was written to allow eg prototyping of script function
prior to converting to plugins and is very flexible in what it can do. If it can be of use, I'm more than happy to see it be used.

EDIT: Also, RT_LumaCorrelation should work just fine in any planar YUV format and supports YUY2 and RGB24/32, if not needed those parts could easily be stripped out.
Another speed tip could be to set interlaced=true, to skip every other scanline and double the speed.
Havnt totally given up on the PlanetCrop thing, still got the script and might have a go using only left and right strips encompassing edges only and not the center.

[Mystery Keeper]

It is totally possible. I toyed with phase correlation using OpenCV (http://opencv.org/).

[StainlessS]

This is the main active part of RT_LumaCorrelation(), after count arrays have been filled

Sxy += Sxy_lo;
__int64 Sx = 0;
__int64 Sy = 0;
__int64 Sx2 = 0;
__int64 Sy2 = 0;

for (i=256;--i>=0;) {
__int64 z;
z = i * cntx[i]; Sx += z; Sx2+= i * z;
z = i * cnty[i]; Sy += z; Sy2+= i * z;
}

double num = ((double)Pixels * (double)Sxy) - ((double)Sx*Sy);

double div1 = ((double)Pixels * Sx2) - ((double)Sx*Sx);
double div2 = ((double)Pixels * Sy2) - ((double)Sy*Sy);

if(div1 < 0.0) div1 = - div1; // fabs
if(div2 < 0.0) div2 = - div2; // fabs

double div;
if ((div1<0.0001)||(div2<0.0001)) {
div = 0.000001; // think this should be 0.00001 but dont really matter
} else {
div = sqrt(div1) * sqrt(div2);
}
return (float)(num/div);


Sxy is Sum of x * y (64 bit)
Sx2 is Sum of x Squared

and here is link to the formula used, the last one in Mathematical Properties and just before Interpretation:
http://en.wikipedia.org/wiki/Pearson_product-moment_correlation_coefficient#Mathematical_properties

EDIT: The above is wrong, that page has had two more formula added since I last accessed it, the used one
is immediately before "The second formula above needs to be corrected for a sample:",
dont know if I should alter RT_LumaCorrelation or not.
Is Population Correlation better or worse for Avisynth purposes than Sample Correlation, suggestions ?
Guess I could implement both population & sample versions via a bool flag. What say you.
Seems to me that Population Correlation is probably the better one, as it is, but what do I know,
I just coded it.

@Jmac698, if RT_LumaCorrelation is in error, then JCorr is too.


This however might be an interesting alternative: Cross Correlation Using FFTW3
http://blog.dmaggot.org/2010/06/cross-correlation-using-fftw3/

[gwendolien]

THANKS again for all responses so far.
I downloaded RT_Stats.dll package and installed it in the autoloud directory.
But now I'm stuck again; using:

GScript("""RT_LumaDifference(dvd, vhs0, w=700, h=494, Matrix=(width<=720?2:3))""")

results in an error reporting that "current_frame" is only defined/available inside a runtime environment.

Being a newby in Avisynth as such, could anybody help me with the use of this function
for the purpose of, in the end, aligning clip2 to be coïnciding as well as possible with clip1 in x,y position

PLEASE!

Many thanks in advance for any suggestion.

[gwendolien]

What I had in mind was:
find the (local) minimum of the LumaDifference between the two clips and shift position of clip2 w.r.t. clip1 accordingly.
So, I thought of using GScript for the construction of the combined for loops like:

diff=255
for (i=-5, 5, 1){
for (j= -5, 5 1){
t_diff = RT_LumaDifference(clip1, clip2,w=700, h=494, x2=i, y2=j, Matrix=(width<=720?2:3))
t_diff < diff ? diff=t_diff : diff=diff
j=j+1
}
i=i+1
}

But I'm running into numerous error messages, which I'm afraid are too complex to understand for me.

So, please...

[StainlessS]

gwendolien,
Can you put up 2 sample clips somewhere (~30 secs each)and I'll have a play with them.

int "Matrix"=(width<=720?2:3)

means that the DEFAULT if matrix NOT supplied is 2(pc601) if clip width <=720 and 3(pc709) if greater than 720, and
is for RGB clips only, not used if YUV.


FRAME = 0
RADIUS = 5
X = RADIUS
Y = RADIUS
W = - RADIUS # As in crop(), clip width relative, same as c.width - X - RADIUS
H = - RADIUS # clip height relative, same as c.height - Y - RADIUS
MATRIX = 2 # 2=PC601 or 3=PC709 # ONLY for RGB, YUV does not matter
BestX = 255 # Nothing in particular
BestY = 255
diff = 255.0
for (i= - RADIUS, RADIUS, 1) {
for (j= - RADIUS, RADIUS, 1) {
t_diff = RT_LumaDifference(clip1, clip2,n=FRAME, x=X,y=Y,w=W,h=H,n2=FRAME, x2=X+i, y2=Y+j,matrix=MATRIX)
if(t_diff<diff) {
diff=t_diff
BestX = i
BestY = j
}
} # j=j+1 NOT needed, its set at the beginning of the loop, closing brace '}' marks the loop end.
}


maybe something like above, untested and incomplete. (You of course need GScript)

[StainlessS]

gwendolien,
Thank you for your PM, can you create an account on Mediafire and upload to there, you can then post a link
to it or PM me with the link if preferred. MediaFire tends to be the preferred upload site on D9.

[StainlessS]

gwendolien,

Took a look at your clip (sorry for not noticing your PM for so long).

Firstly, the Xvid encode of both clips is using (I think) a version of XVid that is perhaps
6 or 7 years old and known to have a problem, was fixed in later versions of XVid, but
some clips that I've come across years ago could still only be decoded by the same bugged
version of XVid. Your clips display with an ugly green bar across top left of clip and
blue/red striped bar bottom right, with weird massive chroma displacement in the rest of the
frame, using Media Player Home Cinema v1.6.2.4902, but OK in Windows Media Player, VLC and
GOM Players, also OK in VDub. If you encoded that clip, dump that codec (although could I
suppose be a decode problem in MPC-HC).

Dont think the cross correlation between the DVD version and analog version could prove
of much worth. DVD version has black band at top of frame, VHS a thicker black bar at bottom,
I think the vhs clip might have a bit of the top scanline missing and been zoomed a little to
be same size as DVD version. I dont think much can be gained by trying to align a better but still
rough version of a clip to a 'lumpy' analog version which also jiggles up and down.
The RT_Correlation thing only works at pixel accuracy, think perhaps to be of any use you would need
sub pixel accuracy and there are probably better filters intended specifically for problems like that,
better comparing with before and after frames to correct vertical jitter (with the odd horizontal
jitter thrown in).


I tried out (for the first time) VDub Deshaker (in VDUB), and maybe you'de like to take a look at that yourself.
http://www.guthspot.se/video/deshaker.htm

I just loaded the clip directly into VD, added the Deshaker filter. In Filter GUI,

PASS_1 (Changes from Default Only)
Source Pixel Aspect : Standard PAL 1.094

Pass 1 Parameters:
Scale : Full (most Precise)
Use Pixels : All (most robust)

Set the logname and dir,

Click PASS_2:

Edge Compensation: Adaptive Zoom Full (no borders)

Click PASS_1 again.

Exit all with OK.

click 3rd little buttom on bottom left, to play both input and output. (you should see motion vectors in right pane).

Rewind.

In Filter GUI, Click PASS_2

Exit all with OK.

Select Codec

Save Output.

Not perfect by any means but a lot better I think, probably mainly down to the missing black line at top left
which made the jiggle more apparent.

Maybe also take a look here if not satisfied (or for an Avisynth alternative):-

http://forum.doom9.org/showthread.php?t=143380&highlight=stab

[gwendolien]

THANKS StainlessS for yoyr extensive reply.
The XVid problems were/are unkonown here: the avi's played well in VLC and VDub here,
so I assumed everything was fine.
In the mean time I had already started some attempts with clip1.Overlay(clip2, mode="subtract") but I keep on stumbling across various Gscript/Runtime problems/errors, based on your earlier RT_LumaDifference suggestion.

code:

function Align_Clip(clip1, clip2, Radius)
{
Gscript("""
#RADIUS = 5
X = RADIUS
Y = RADIUS
W = - RADIUS # As in crop(), clip width relative, same as c.width - X - RADIUS
H = - RADIUS # clip height relative, same as c.height - Y - RADIUS
MATRIX = 2 # 2=PC601 or 3=PC709 # ONLY for RGB, YUV does not matter
BestX = 255 # Nothing in particular
BestY = 255
diff = 255.0
for (i= - RADIUS, RADIUS, 1) {
for (j= - RADIUS, RADIUS, 1) {
tclip=clip1.Overlay(clip2, mode="subtract", x=X, Y=Y)
t_diff=AverageLuma(tclip)
if(t_diff<diff) {
diff=t_diff
BestX = i
BestY = j
}
}#j=j+1 loop
}#i=i+1 loop
""")
Return clip2.Subtitle("X= "+String(BestX)+" / Y= "+String(BestY))
}


I get invalid arguments to AverageLuma in line x,y,z, ....

Apparently the brain here is almost impermeable to new things like scope of variables or so...

Could you or somebody help me out here?

[Gavino]

I get invalid arguments to AverageLuma
You need to set the variable current_frame (eg to zero) before calling AverageLuma.

However, your loop logic is clearly wrong too as 'tclip' does not change between iterations. I suspect you need to vary the x and y arguments to Overlay in some way according to i and j.

[gwendolien]

Gavino: THANKS!
Yes, there were errors with the i and j, stupid error here.
But I do not understand the current_frame thing at all.
Setting it to zero, indeed evaluates frames 0 of both clips.
However, this should be done for all frames of both clips.
And I do not see anywhere a call to the variable current_frame in my script.
How to get there?

[StainlessS]

Dont know if below is of any interest to you, was something I was playing with, just for comparing RT_LumaCorrelation and RT_LumaDifference


DVD=AVISource("D:\dvd_sample.avi") #.GreyScale()
VHS=AVISource("D:\vhs2_sample.avi") #.GreyScale()

Function xLuma(clip c1,clip c2,int n,int rad) {
X = rad*2
Y = rad*2
W = - rad*2 # As in crop(), clip width relative, same as c.width - X - RADIUS
H = - rad*2 # clip height relative, same as c.height - Y - RADIUS
BestX = 255 # Nothing in particular
BestY = 255
diff = -1.0
GScript("""
for (i= - rad, rad) {
for (j= - rad, rad) {
th = RT_LumaDifference(c1, c2,n=n, x=X,y=Y,w=W,h=H,n2=n, x2=X+i, y2=Y+j)
th=((255.0-th)/127.5)-1.0 # LumaDiff as for LumaCorrelation, (-1.0 -> 1.0)
if(th>diff) {
diff=th
BestX = i
BestY = j
}
}
}
S="Lx="+String(BestX)+" Ly="+String(BestY)+" Lth="+String(th)
""")
Return S
}
Function xCorr(clip c1,clip c2,int n,int rad) {
X = rad*2
Y = rad*2
W = - rad*2 # As in crop(), clip width relative, same as c.width - X - RADIUS
H = - rad*2 # clip height relative, same as c.height - Y - RADIUS
BestX = 255 # Nothing in particular
BestY = 255
diff = -1.0
GScript("""
for (i= - rad, rad) {
for (j= - rad, rad) {
th = RT_LumaCorrelation(c1, c2,n=n, x=X,y=Y,w=W,h=H,n2=n, x2=X+i, y2=Y+j)
if(th>diff) {
diff=th
BestX = i
BestY = j
}
}
}
""")
S="Cx="+String(BestX)+" Cy="+String(BestY)+" Cth="+String(th)
Return S
}

frms=min(DVD.FrameCount,VHS.FrameCount)

GScript("""
for(n=1,frms-1) {
SL=xLuma(DVD,VHS,n,5)
SC=xCorr(DVD,VHS,n,5)
Eval(SL) Eval(SC)
BADS=(Lx==Cx && Ly==Cy) ? "" : " MISMATCH"
if(BADS!="") {
RT_Debug(SL,"::::",SC,BADS)
}
}
""")



Return DVD
#Return VHS
Return StackHorizontal(DVD,VHS)

[Gavino]

I do not understand the current_frame thing at all.
Setting it to zero, indeed evaluates frames 0 of both clips.
However, this should be done for all frames of both clips.
And I do not see anywhere a call to the variable current_frame in my script.
How to get there?
The simplest way would be to put the call to your Align_Clip function inside a call to GScriptClip (from the GRunT plugin), eg
GScriptClip("Align_Clip(c1, c2, 5)")

The function definition itself would not need to be changed for this, except for removing the explicit assignment to current_frame which you added last time. When you run the script, current_frame will be automatically updated for each frame and the subtitles will show you the results for each frame in turn.

(GScriptClip() is needed as the built-in ScriptClip() function does not support calling run-time functions like AverageLuma inside a user function.)

[gwendolien]

THANKS, both to you StainlessS and Gavino!
I think I've got what I wanted, finally.
In principle my problem is solved, now I "only" have to find out whether the colourmatch between the clips can be made (almost) identical....and that proces is also a challenge.
But it seems that I have found a solution there too:
Usage of Colourlike and consecutively thereafter HistogramAdjust on the single channels R,G and B and then recombining those to an RGB etc. does the job!!
Again, my appreciation for you both!!

[gwendolien]

One final question:

HOW does one produce such a neat "code-window" as in your examples?

Can't seem to find guidance on the forum on this subject....

[StainlessS]

go advanced select the code and click on hash symbol (i think, i'm on tablet)

EDIT: Look at what gets put around the CODE, you can do that by hand also.