RedAverage plugin (version 1.4.3 biiiiig bugfix)

[redfordxx]

EDIT:I am removing the original post which is in IanB's answer anyway and replacing it with some info about what's done.
This plugin started with the idea of improving original mg262's Average. Here I thank him for his work.

This new plugin has diferent names to avoid name conflicts and there are currently following filters:
RAverageM - masked average
RAverageW - weighted average
RMerge - yet another merge filter


common parameters:
y,u,v=3 ...process plane, other numbers disable processing the plane atm for RAverageW and RAverageM. For RMerge it works just like in masktools
bias=0 ...number which is added to finally computed value (just fyi: bias=-0.5 kills rounding the result)
sse=-1 ...limits the processor capabilities. -1=use all (bitwise...16,8,4=SSE4.1,SSSE3,SSE2)...0=no asm, which is veeeeryyyy sloooow
lsb_in=false ...if true, 16bit stacked input is expected (for more detail explanation about 16bit clips find Dither tools)
lsb_out=false ...if true, 16bit stacked output is produced
All parameters are optional, except the input clips. There must be at least one input clip pair.

RAverageW(clip1, weight1, clip2, weight2, ... clipn,weightn,bias,y,u,v,sse, bool lsb_in, bool lsb_out, int mode, int errlevel, bool silent)
computes value of weighted average of clips:
result=c1*w1+c2*w2+...cn*wn+bias
silent=false ... if true, certain error messages (mostly about rounding error) are suppresed
errlevel=0 ...sets the maximum rounding error you are willing to accept as a tradeoff to faster method
mode=-1 ...bitwise 4, 8 for enabling method 4, 8. 0 results in unoptimized C++.
- mode 8: requires SSE2, has higher precision than 4 but depending on hardware, probably will be slower
- mode 4: requires SSSE3, is faster, however has certain limitations in precision (if you don't have silent true, you will be warned) (8bit in only atm, 8 or 16 out)
basically method 4 is fully precise, when weights are nice power-of-two-numbers, specifically: w_i=a_i*2^k (a_i are integers and k is integer, all signed...moreover sum[max(a_i,0)]<128 and sum[min(a_i,0)]>-128
...in other words, weights must be scalable to signed byte;-)
- mode 0: plain C++, slow, fully precise (I think) and supports fully 8 & 16bit clips in and out.
Concerning 16bit support and conversion when the in and out have different bitdepth, the values are relative to full scale. That introduces multiplication or division by 256 in case of change of bitdepth. So: 8bitvalue*weight*256=16bitvalue.

RAverageM(clip1, mask1, clip2, mask2, ... clipn,maskn,bias,y,u,v,sse, bool lsb_in, bool lsb_out)
computes value of masked average of clips:
result=c1*m1+c2*m2+...cn*mn+bias

RMerge(clip1, clip2, mask,y,u,v,sse, bool lsb_in, bool lsb_out, int mode)
Merges two clips based on mask, as expected. However, there are some new things. First, again, 16bit support (although not fully optimized).
Then there are two modes in case 8bit input which decide, whether the mask is applied slightly differently to correct the problem of incomplete range of the mask:
mode=255 ...this is standard merge formula: r=c1*(256-m)+c2*m
mode=256 ...this is adjusted merge formula: r=c1*(256-n)+c2*n where n= ( m<=128 ? m : m+(m-128)/128 )
lsb_in=true ...this is 16bit merge formula: r=c1*(65536-m)+c2*m
Then, r is scaled and rounded to results r8 or r16, depending on output bitdepth.
The goal of mode 256 is to allow merge values in full range and if m=255, then r8=c2 ... (...and not r8=(c1+c2*255+128)/256 like usually and in mode 255)
Only lsb_in=lsb_out=false and mode=255 is SSSE3 optimized at this moment. Anything else goes plain C++. But the optimized version performed in my configuration better than mt_merge.

Examples:

RAverageW(c1, 0.5, c2, 0.5) #same as mt_average
RAverageW(c1, 1, c2, -1, bias=128) #same as mt_makediff
RAverageW(c1, 1, c2, 1, bias=-128) #same as mt_adddiff
a2=RAverageW(a1, -1, bias=255) #same as mt_invert
RAverageW(c1, a1, c2, a2) #similar to mt_merge
RMerge(c1, c2, m) #same as mt_merge
RMerge(c1, c2, m, mode=255) #same as mt_merge
RMerge(c1, c2, m, mode=256) #not same as mt_merge

#restoring blended telecined video
f0=o.SelectEvery(5,0)
f1=o.SelectEvery(5,1)
f2=o.SelectEvery(5,2)
f3=o.SelectEvery(5,3)
f4=o.SelectEvery(5,4)
n3=RAverageW(f1,-0.5,f2,1,f3,1,f4,-0.5)
return Interleave(f0,f1,n3,f4)
Changelog:

1.4.3 (http://www.mediafire.com/download.php?m863pgn8ylkkyv3) serious bugfix, thanx to Bloax for pointing out
1.4.2 (http://www.mediafire.com/?d6kbqs4g2swwp6h) bugfix of RMerge SSSE3, switched clips in RMerge to match the mt_merge formula, minor improvements
1.4.1 (http://www.mediafire.com/?zwgb2pogpaalza3) SSSE3 optimization of RMerge in 8bit
1.4.0 introduction of of new filter RMerge
1.3.10 full 16bit support for RAverageW mode 8 (SSE2) and minor speed improvements
1.3.9 (http://www.mediafire.com/?q4a70w4i83gmdid) added 16bit support for RAverageW, however not always SIMD optimized. Multiple algorithms available, different in speed and precision. Thorough estimation of rounding error at the beginning added, so that the most suitable algo is chosen. Ugly source included.
1.3.5 something I am not afraid to publish, but I am afraid to publish the messy src

Todo's (maybe):
- optimization
- mixed bitdepth for merge masks
- RTotal

I did some speed benchmarking, so why not to share it. It is rough measurement, with C++ code as a reference speed:
Filter method 2 clips 4 clips 8 clips
Average C++ 89%
RAverageW C++ 100% 100% 100%
Average SIMD 388% 354% 351%
RAverageW 8 - SSE2 1213% 906% 738%
RAverageW 4 - SSSE3 1213% 957% 1176%
I think that method 4 and 8 have same results in case of lower number of clips, because on my computer, there is memory-processor data transfer speed the limiting factor. I don't know how on other machines, but I hope 4 is faster in all cases.

Disclaimer:
This filters are in development stage and cannot guarantee no bugs, so feel free to use it, test it and comment it.
If you report a bug please specify it in a reproducable clip. For instance:
c=ColorBars(pixel_type="YV12")
m=mt_lutspa(c,expr="x y * 2 256 * *", mode="relative", y=3,u=3,v=3).trim(1,1).Loop(c1.framecount,0,0)
RAverageM(m,c,m,c,bias=-10000,lsb_in=false,lsb_out=true)
Moreover, if it's a crash, try to play with the sse parameter first.
Use this thread only for development issues, for questions how do I... use Average plug-in (http://forum.doom9.org/showthread.php?t=100626) thread

I am adding new version's links to the changelog.

[IanB]

Hi, I am about to come to Average plugin again and add some new features. First I have few questions

1) which level of SSE is suitable (high enough to perform and low enough to be compatible and not overflow. Basically I need following operations (appropriately packed)
B*F
W*F
B*B
W*W
B*B*F
W*W*F
Where B=packed bytes, W=packed words, F=one float or float scaled to byte or word (this one will multiply all the packed numbers)
Is SSE2 ok there can be benefit of something elseI generally try to restrict myself to plain old MMX first time out. This gives a base line for improvement and will work everywhere. Moving up to ISSE gives you the rest of the MMX instructions like PAVGB, PMULUSW, PSHUFW, etc. Just about everything does ISSE now. Straight SSE give floating point but no useful integer stuff and there are lots of annoying instruction holes. SSE2 fills out the lack of integer stuff and gives double support. All new cpus do SSE2 now. SSE3, SSSE3 and up fill in special instruction holes, but your code must check that the cpu actually has the required instructions.

Many cpu's only do SSE2 as pairs of 64bit ops internally which can extend the instruction latency badly leading to pipe stalls that you cannot get rid of no matter how you shuffle your code. My old 3GHz P4 HT generally runs SSE2 code the same or slower than the MMX/ISSE code, sometime a lot slower. My less old Core2 generally runs SSE2 code the same as the best MMX/ISSE code, sometime a little faster, very occasionally a lot faster. I3 and up start to rock but also have reduced latency and improved throughput on MMX which make you have to work hard getting the SSE2+ code significantly faster. Wander through the x264 code and discussions on how various things bite and don't work as expected.
2) what happens if there is overflow in SIMD?Generally your choice, either wrapping or saturation depending on the instructions selected, eg PADD verse PADDUSB.
3) can I uses mmx and xmm register at the same time?Yes. The MMX<->XMM transfer instructions have some braindeadness on many cpu's, usually extreme latency sometime poor throughput as well..
4) I somewhere read that reading 16 bytes of data from memory to xmm is more than 2x slower than 2x 8bytes data from memory. So it means it is better to read twice and then combine...is that true?It depends on the cpu model and memory controller. Movdqa is always at least equal fastest. Movdqu can massively suck, 2 movq's may be faster and have lower latency, Lddqu has it's problems as well. You have to test your case to know, then the test only applies to your cpu with your memory.
5) what is the guaranteed alignment of data in plane...i.e. length of row in the input clip is multiplication of what?Rowsize is always the width.

Pitch is always at least (width+15)/16*16 for all planes. You can read and write the data in the area between pitch and rowsize. It must be considered uninitialised (it will be random old frame data).

Non-aligned crops can make the start address of a row non-aligned, but the original VideoFrameBuffer would have been aligned so you can backstep the row start to regain alignment. All data outside the active picture area must be considered uninitialised.

[redfordxx]

I am looking though my old code and what I see

env->GetCPUFlags() & CPUF_INTEGER_SSE) != 0


punpcklbw mm2, mm7//low
Testing SSE, running SSE2 and using MMX reg


This is I think not good, I guess?

[redfordxx]

1) From what you wrote, basically I should not bother too much with the speed issues, since every machine has it differently implemented, correct? Anyway would you recommend some document where is clock cycles or whatever speed specification of the instruction? Maybe to learn whether pxor is faster than movq on mmx...
2) MMX<->XMM still should be faster than from memory on every machine, or not?
3) So since the pitch is mod16, basically it is safe (and good idea???) to run one cycle from 0 to pitch*height, instead of two cycles inside each other for height and width?
4) Saturation means 250+10=255?
5) Should I be interested in MOVNTDQ? I don't know what is the benefit of nontemporal

[IanB]

I am looking though my old code and what I see
env->GetCPUFlags() & CPUF_INTEGER_SSE) != 0

punpcklbw mm2, mm7//low
Testing SSE, running SSE2 and using MMX reg
CPUF_INTEGER_SSE are really the extra MMX instructions, not to be confused with CPUF_SSE (although the same CPU model's mostly provide both, except Athlon plain (not the XP/MP)).

Your capability testing should always match your instruction usage.

1) From what you wrote, basically I should not bother too much with the speed issues, since every machine has it differently implemented, correct?No, I mean you should not assume SSE2 is always faster than MMX, you need to test it. On some hardware SSE2 code might be twice as fast as an MMX version, but it also might be slower.
Anyway would you recommend some document where is clock cycles or whatever speed specification of the instruction? Maybe to learn whether pxor is faster than movq on mmx...Agner Foggg.
2) MMX<->XMM still should be faster than from memory on every machine, or not?The killer is latency not throughput, bad CPU's can stuff around for more than 9 cycles. A from L1 cache movdqa can be all done in 3 cycles. You need to test it.
3) So since the pitch is mod16, basically it is safe (and good idea???) to run one cycle from 0 to pitch*height, instead of two cycles inside each other for height and width?No, Cropped frames may result in a tiny active hole in the middle of a hugh VideoFrameBuffer. E.g rowsize might be 16 and pitch 1920.
4) Saturation means 250+10=255?Yes, and 10-250=0 and 32760+100=32767 and -32760-100=-32768, etc
5) Should I be interested in MOVNTDQ? I don't know what is the benefit of nontemporalIn general I have found MOVNT* only good for sequential read once data. On output it may make "your" filter appear faster, but it invariably makes the next filter slower. It is very difficult to predict the effect across a range of machines. You need to test it.

In the HResizer we read 1 row of input, unpack it to a temp array, then FIR the living hell out of that temp array. By reading the source rows with MOVNTQ we avoid flushing the start of the temp array from L1 cache before the first FIR operation. On most machines the temp array and FIR coefficients live in the L1 cache and hece the code screams along. This is an exceptional case and much effort by many people went into getting this optimal.

[redfordxx]

So basically if I understand it correctly (what latency and rec thruput means) I made a little breakdown of my code.
What it does is load data from memory to mm0, mm2, interleaves with zeros, multiplies and adds to mm4,mm5
The expressions in comment means:
clock for operation: reciprocial thru / latency --> clock when data ready spatial:
mov eax, [length4] //n // 23: r1.0 / l2 --> 25
movq mm4, mm6 // 24: r0.3 / l1 --> 25
movq mm5, mm6 // 24: r0.3 / l1 --> 25
temporal:
mov ecx, [sourcep+eax-8] // 1: r1.0 / l2 --> 3
mov edx, [sourcep+eax-4] // 2: r1.0 / l2 --> 4
movq mm0, qword ptr [ecx+ebx-8] //3: r1.0 / l2 --> 5
movq mm2, qword ptr [edx+ebx-8] //4: r1.0 / l2 --> 6
movq mm1, mm0 //frame 5: r0.3 / l1 --> 6
movq mm3, mm2 //mask 6: r0.3 / l1 --> 7 (wait for mm2)
// pxor mm1, mm1
// pxor mm3, mm3
// mm7=000000000000
punpcklbw mm0, mm7//low 7: r1.0 / l1 --> 8
punpcklbw mm2, mm7//low 8: r1.0 / l1 --> 9
punpckhbw mm1, mm7//high 9: r1.0 / l1 --> 10
punpckhbw mm3, mm7//high 10: r1.0 / l1 --> 11
// punpckhbw mm1, mm0//high
// punpckhbw mm3, mm2//high
pmullw mm2, mm0//low 11: r1.0 / l3 --> 14
pmullw mm3, mm1//high 12: r1.0 / l3 --> 15
// pmulhw mm3, mm1//high
sub eax, 8 //13:
paddusw mm4, mm2//low 14: r0.5 / l1 --> 15
paddusw mm5, mm3//high 15: r0.5 / l1 --> 16 (wait for mm3)
jnz temporal //16: r1.0 / l0 --> 17
psrlw mm4, 8//low 17: r1.0 / l1 --> 18 ???
psrlw mm5, 8//high 18: r1.0 / l1 --> 19 ???
packuswb mm4, mm5 // 19: r1.0 / l1 --> 20
movq qword ptr [edi+ebx-8], mm4 // 20: r1.0 / l3 --> 23
sub ebx, 8 // 21:
jnz spatial // 22:

1)I believe I could save one clock with the change in comments, (then I have zeros in lower bytes in mm1,mm3, so I used pmulhw) however I ended up with real but I ended up with really ugly output.
2)Maybe I can save even more with punpck directly from memory, but I nowhere found documented latency for punpck from memory...and again then I will have to do this pmulhw god knows what will happen
3) Q: Latency is relevant only for the destination argument or both? for example with paddusw mm4, mm2 I have to wait for next use of mm4 but can read/write mm2 imediately? Or both mm are blocked?
4) Q: operations on reg32 and mm go in one line or parallel? I mean eg sub and pmullw go in same time or after each other?

EDIT: all this can be thrown away, if I misunderstood latency and thruput

[redfordxx]

What is the best way to make conditional jump if eax is multiply of 2?
Thank you for your time...
R.

[jmartinr]

What is the best way to make conditional jump if eax is multiply of 2?
Thank you for your time...
R.

test eax, 1 ?

[IanB]

I refuse to try to read stuff that is wider than the screen. Edit your post to collapse the tabs and I might then take an interest. :devil:

Thanks, it make some sense when you can see the whole picture. ;)

[redfordxx]

And then jnz i suppose, Thanx


I have one little different question, about calling the filter. Atm, this filter requires variable amount of arguments. (Clip, float, clip, float, clip float,...) I want to add more optional named arguments :
y,u,v integers (obviously for plane processing)
bias float
How can I detect what in fact are the argument? FYI, the code bellow how it works now, but do understand that I did't create it, I used and modified from somewhere else.
struct WeightedClip
{
PClip clip;
double weight;
WeightedClip(PClip _clip, double _weight)
: clip(_clip), weight(_weight) {}
};

AVSValue __cdecl Create_RAverageW(AVSValue arguments, void *user_data, IScriptEnvironment *env)
{
AVSValue array = arguments[0];
int noarguments = array.ArraySize();
if(noarguments % 2 != 0)
env->ThrowError("Average requires an even number of arguments.");

vector<WeightedClip> clips;
for (int i=0; i < noarguments; i+=2)
clips.push_back(WeightedClip(array[i].AsClip(), array[i+1].AsFloat()));

return new RAverageW(clips [0]. clip, clips, env);
}


extern "C" __declspec(dllexport) const char *__stdcall AvisynthPluginInit2(IScriptEnvironment *env)
{
env->AddFunction("RAverageW", ".*", Create_RAverageW, 0);
return "Average plugin";
}

[Gavino]

Atm, this filter requires variable amount of arguments. (Clip, float, clip, float, clip float,...) I want to add more optional named arguments :
y,u,v integers (obviously for plane processing)
bias float
How can I detect what in fact are the argument? FYI, the code bellow how it works now, but do understand that I did't create it, I used and modified from somewhere else.
".*" in the call to AddFunction specifies that it accepts a variable number of arguments (of any type). These are packaged into a single array argument, arguments[0].
You can add further named arguments, eg ".*[y]i[u]i[v]i" would add three ints, y, u v, available as arguments[1], arguments[2] and arguments[3].
I notice your code does not check the type of the variable arguments - it really should, instead of assuming they are correct.

[redfordxx]

I refuse to try to read stuff that is wider than the screen. Edit your post to collapse the tabs and I might then take an interest. :devil:I fully understand, I am sorry...I collapsed it. Starting from the bottom are the most important things.

[redfordxx]

I notice your code does not check the type of the variable arguments - it really should, instead of assuming they are correct.Do I need to check only the .* ones? How do I check (sorry for my ignorance)?

[Gavino]

Do I need to check only the .* ones? How do I check (sorry for my ignorance)?
Yes, as "." indicates 'any type'; all others have a specific type and this will be checked by the parser before calling your plugin.

You can check by using the IsXXX() functions from avisynth.h.
For example,
if (array[i].IsClip() {
... do something with array[i].AsClip() ...
}
else
env->ThrowError("...<a suitable error message>...");

[redfordxx]

Well, after one day figuring out that it does not work coz pmulhw is signed, I think I succeeded in some speedup, but hard to tell because I measure it in TaskManager on CPU Time and it varies+-10%...I don't know, maybe caching situation.
I found somewhere cycles.h and I guess it is something I would might need to measure... but I get error message
1>C:\Program Files (x86)\Microsoft Visual Studio 10.0\VC\include\Cycles.h(201): warning C4405: 'ret' : identifier is reserved word
1>C:\Program Files (x86)\Microsoft Visual Studio 10.0\VC\include\Cycles.h(463): error C3861: 'elapsed': identifier not found
Anyone has experience with this and may have advice?

[IanB]

So basically if I understand it correctly (what latency and reciprocal throughput means) I made a little breakdown of my code.
What it does is load data from memory to mm0, mm2, interleaves with zeros, multiplies and adds to mm4,mm5
The expressions in comment means:
clock for operation: reciprocal throughput / latency --> clock when data ready

.....

1)I believe I could save one clock with the change in comments, (then I have zeros in lower bytes in mm1,mm3, so I used pmulhw) however I ended up with real but I ended up with really ugly output.To do what I think you want you would need the ISSE instruction pmulhuw (unsigned) instead of the MMX instruction pmulhw (signed). The results may not be bit identical as the lower bits can combine to increment the total result, but some cleverness can avoid this.
2)Maybe I can save even more with punpck directly from memory, but I nowhere found documented latency for punpck from memory...and again then I will have to do this pmulhw god knows what will happenGiven you need the same data in 2 places (high & low), this probably won't work out for you here. Also the bytes end up in the wrong half of the word so you need to spend extra instructions repairing this, but again some cleverness can avoid this.
3) Q: Latency is relevant only for the destination argument or both? for example with paddusw mm4, mm2 I have to wait for next use of mm4 but can read/write mm2 imediately? Or both mm are blocked?No latency is when the result is available, doing padd mm0, mm7 followed by padd mm1, mm7 incur no stall. However if the next instruction needs mm1, i.e. padd mm0, mm1 then you must wait for mm1 to be ready.
4) Q: operations on reg32 and mm go in one line or parallel? I mean eg sub and pmullw go in same time or after each other?Yes IA32 instruction can be interleaved with MMX/SSE instructions. It is worth noting that MMX/SSE instructions do no effect the flags so you can move an IA32 test instructions before an MMX/SSE instruction to fill a stall.
EDIT: all this can be thrown away, if I misunderstood latency and reciprocal throughput
You seem to have the basic idea, but it is so much more complicated that you are predicting. Not all the execution units have all the same hardware, so some instruction will run in series even if the registers involved are completely independent. Also modern CPU have very advanced Out Of Order execution, this make prediction of what happens when very difficult. Atom has no OOO so for that platform it matters a lot.


What does seem apparent is you are processing memory backwards. This usually kill the hardware prefetch and stuffs performance on modern CPU's. Previously it was used to avoid pushing needed data from the cache with very large data sets.


I had a bash at re-ordering your code, it's possibly not right algorithmically, but it may give you some ideas. // mm7=0
spatial:
mov eax, [length4] // 23: r1.0 / l2 --> 25
movq mm4, mm6 // 24: r0.3 / l1 --> 25
movq mm5, mm6 // 24: r0.3 / l1 --> 25

temporal:
mov ecx, [sourcep+eax-8] // 1: r1.0 / l2 --> 3
mov edx, [sourcep+eax-4] // 2: r1.0 / l2 --> 4
movq mm0, qword ptr[ecx+ebx-8] // 3: r1.0 / l2 --> 5+x
movq mm2, qword ptr[edx+ebx-8] // 4: r1.0 / l2 --> 6+x
#ifdef ISSE
punpckhbw mm1, mm0 // high : r1.0 / l1 --> (wait for mm0)
punpcklbw mm0, mm7 // low : r1.0 / l1 -->
punpckhbw mm3, mm2 // high : r1.0 / l1 --> (wait for mm2)
punpcklbw mm2, mm7 // low : r1.0 / l1 -->
pmulhuw mm3, mm1 // high : r1.0 / l3 --> ISSE Required!
sub eax, 8 // :
pmullw mm2, mm0 // low : r1.0 / l3 -->
paddusw mm5, mm3 // high : r0.5 / l1 -->
paddusw mm4, mm2 // low : r0.5 / l1 -->
#else
movq mm1, mm0 // frame : r0.3 / l1 --> (wait for mm0)
punpcklbw mm0, mm7 // low : r1.0 / l1 -->
movq mm3, mm2 // mask : r0.3 / l1 --> (wait for mm2)
punpcklbw mm2, mm7 // low : r1.0 / l1 -->
punpckhbw mm1, mm7 // high : r1.0 / l1 -->
pmullw mm2, mm0 // low : r1.0 / l3 -->
punpckhbw mm3, mm7 // high : r1.0 / l1 -->
paddusw mm4, mm2 // low : r0.5 / l1 -->
pmullw mm3, mm1 // high : r1.0 / l3 -->
sub eax, 8 // :
paddusw mm5, mm3 // high : r0.5 / l1 -->
#endif
jnz temporal // 16: r1.0 / l0 --> 17

psrlw mm4, 8 //low 17: r1.0 / l1 --> 18
psrlw mm5, 8 //high 18: r1.0 / l1 --> 19
sub ebx, 8 // 19: --> 20
packuswb mm4, mm5 // 20: r1.0 / l1 --> 21
movq qword ptr[edi+ebx], mm4 // 21: r1.0 / l3 --> 23 (wait for mm4)
jnz spatial // 22:

[IanB]

As others point out using ".*" in your argument list overrides type checking of the arguments and the number of arguments.

So your creator code must validate what it is expecting.

[redfordxx]

Hi,
I would like to do something like this:



bool SSEMMXenabled = (env->GetCPUFlags() & CPUF_INTEGER_SSE) !=0;
bool SSE2enabled = (env->GetCPUFlags() & CPUF_SSE2) !=0;
if (optimisable && SSE2enabled && (sse & 2 ))
{
#define REGSIZE 16
#define MOVALL MOVDQA
#define m0 xmm0
#define m1 xmm1
....
if (y==3) {
#define PLAN PLANAR_Y
#include "code.asm"
}
if (u==3) {
#define PLAN PLANAR_U
#include "code.asm"
}
if (v==3) {
#define PLAN PLANAR_V
#include "code.asm"
}
}
else if (optimisable && SSEMMXenabled && (sse & 1 ))
{
#define REGSIZE 8
#define MOVALL MOVQ
#define m0 mm0
#define m1 mm1
....
if (y==3) {
#define PLAN PLANAR_Y
#include "code.asm"
}
if (u==3) {
#define PLAN PLANAR_U
#include "code.asm"
}
if (v==3) {
#define PLAN PLANAR_V
#include "code.asm"
}
}
else
{
if (y==3) averageplaneW(PLANAR_Y, env);
if (u==3) averageplaneW(PLANAR_U, env);
if (v==3) averageplaneW(PLANAR_V, env);
}

when code.asm looks like
__asm
{
pushad
....
emms
popad
}
I think from this example is my intention clear. Is that proper way to do this?

[redfordxx]

My measurements:
I am doing it on 8192x4096x500 blankclips and I look in taskmanager at CPU time. I don't know nothing better.
These are my results:
clips 1 4 8 16
C++ 59 137
Original 17 38 77
MMX 9 29 62 109
XMM 10 16 37 59
In the tible is CPU time and in top is number of the input clips to average. Original is the average plugin published in the other thread. MMX and XMM are my versions.
Clearly there is an overhead of the cycles thru the plane about 9 sec. Which is significant, especially when small number of clips. But I tried.
Moreover it seems to me obvious that with the number of clips there is increase of time due the fact, that the processor doesnot read 16 clips in paralel efficiently. Caching problem. Is there any way to help? Like tell the machine that the caching should be that way and not that way?
Now I will try to process it with using all XMM and MM reg, maybe it will boost a little.

[redfordxx]

Small question. For pointers and counters, do I have anything else available to use than eax,ebx,ecx,edx,edi,esi?
Coz if I am going to do stacked 16bits I would appreciate more, if there is possibility.
Also, the compiled always warns me: frame pointer register 'ebx' modified by inline assembly code.
Everything seems to work OK though, so should that bother me?

[cretindesalpes]

For pointers and counters, do I have anything else available to use than eax,ebx,ecx,edx,edi,esi?

Unless some wizardry, I don't think so.

Also, the compiled always warns me: frame pointer register 'ebx' modified by inline assembly code.
Everything seems to work OK though, so should that bother me?

http://msdn.microsoft.com/en-us/library/k1a8ss06.aspx
Some SSE types require eight-byte stack alignment, forcing the compiler to emit dynamic stack-alignment code. To be able to access both the local variables and the function parameters after the alignment, the compiler maintains two frame pointers. If the compiler performs frame pointer omission (FPO), it will use EBP and ESP. If the compiler does not perform FPO, it will use EBX and EBP. To ensure code runs correctly, do not modify EBX in asm code if the function requires dynamic stack alignment as it could modify the frame pointer. Either move the eight-byte aligned types out of the function, or avoid using EBX.

Also: http://forum.doom9.org/showthread.php?t=100374

[redfordxx]

Well, sorry, I never nowhere learned C or asm except now for the purpose of Avisynth, so I am a little slow. But after slowly getting thing, I wonder whether this should be enough to safely use all registers for my purposes, as long as I do not use POP and PUSH?
int saveesp;
__asm
{
pushad
mov saveesp, esp
....
mov esp, saveesp
emms
popad
}
EDIT: now I realized ebp seems to by used when I am getting variable from memory. So, maybe, I could use ebp, but everytime I reference memory like eg mov eax, [variable] ebp value will be overwritten...and I should keep it in mind.

[redfordxx]

cretindesalpes, by the way, as you are the guy of 16bits, I would appreciate your opinion on this (or whoever else's):
When I will extend the Average plugin to 16bits, the formula should beint(bias+sum(clip_i*mask_i)/65536+0.5) i=1...nBut I believe doing it like this:int(bias+0.5)+sum(int(clip_i*mask_i/65536)) i=1...ncould be easier and faster code, I think. However, there will be some inaccuracy in the lsb. What do you think of it?

[IanB]

No, ESP always need to point at a valid stack except when interrupts are disabled.

No leaving out the proper rounding causes problems.

Implement as :-int K=(bias<<16)+32768;
....
(K + sum(clip_i*mask_i) ) >> 16; i=1...n

[redfordxx]

Well, I achieved one of my benchmark. One of the function, in following, parameters:
RAverageW(c1,1,c2,-1,bias=128)
Does the same as mt_makediff and I believe it is tiny bit faster on xmm. But, of course, you can choose multiple number of clips and different weights.
There is one think I am fighting with:
It is calculated scaled to signed short with scale 32*256 and as soon as the weight is 4, it switches to non-asm.
I need variable scale which I tried with following thing but doesnot work:
if (maxweight<4) {
#define SCALE (256*32)
#define SCALEPOWER (8+5)
} else {
#define SCALE (256/2)
#define SCALEPOWER (8-1)
}Either this is wrong approach or I have bug somewhere. Bug is up to me to find, but pls tell me, if this kind of if else define is possible.

[redfordxx]

int K=(bias<<16)+32768;
....
(K + sum(clip_i*mask_i) ) >> 16; i=1...n
Well then it would be slow or even slower.
As of my knowledge I see two instruction I can use:
pmaddwd however, there still can be error on last bit. Because this instruction is signed, I have to do
(K + sum((clip_i>>1)*(mask_i>>1)) ) >> 14; i=1...n
Or pmuludq which would be precise but processes fourtimes less data... and will be four times slower

I don't see other options but maybe there are instructions I can't think of.

[IanB]

...
It is calculated scaled to signed short with scale 32*256 and as soon as the weight is 4, it switches to non-asm.
I need variable scale which I tried with following thing but doesnot work:
if (maxweight<4) {
#define SCALE (256*32)
#define SCALEPOWER (8+5)
} else {
#define SCALE (256/2)
#define SCALEPOWER (8-1)
}Either this is wrong approach or I have bug somewhere. Bug is up to me to find, but pls tell me, if this kind of if else define is possible.
I take it you actually mean something like this :-if (maxweight<4) {
#define SCALE (256*32)
#define SCALEPOWER (8+5)
#include "common_include_code.hpp"
#undef SCALEPOWER
#UNDEF SCALE
} else {
#define SCALE (256/2)
#define SCALEPOWER (8-1)
#include "common_include_code.hpp"
#undef SCALEPOWER
#UNDEF SCALE
}You cannot mix program flow with macro substitution. The C preprocessor does all the macro parsing then the compiler compiles the resulting source. If you ask for an ASM listing with C code as the comments you can see what you actually compiled.

Perhaps you need the power of Softwire or similar to generate dynamic assembly on the fly.

[IanB]

int K=(bias<<16)+32768;
....
(K + sum(clip_i*mask_i) ) >> 16; i=1...nWell then it would be slow or even slower.
As of my knowledge I see two instruction I can use:
pmaddwd however, there still can be error on last bit. Because this instruction is signed, I have to do
(K + sum((clip_i>>1)*(mask_i>>1)) ) >> 14; i=1...n
Or pmuludq which would be precise but processes fourtimes less data... and will be four times slower

I don't see other options but maybe there are instructions I can't think of.
The assumption was from your earlier code, where you zeroed a register then looped about summing the clip_i*mask_i values. The significant hint here was precalculate K and start with K instead of zero.

Really this guessing in the dark is not helpful. Please post complete code fragments and ask direct questions about that code.

As I guess your problem it is to do sum(clip_i * mask_i) quickly.

If so you need code to do DWORD=K, Loop{UWORD=BYTE*BYTE, DWORD+=UWORD}, ScaleAndRound(DWORD)

Am I close with the above pseudo code guess :confused:

[jpsdr]

Small question. For pointers and counters, do I have anything else available to use than eax,ebx,ecx,edx,edi,esi?


Under 32bits, no. On 64bits, yes, you have r8 to r15 avaibles.

[redfordxx]

If so you need code to do DWORD=K, Loop{UWORD=BYTE*BYTE, DWORD+=UWORD}, ScaleAndRound(DWORD)
In 16bits I need two unsigned words to multiply. If UWORD means unsigned word and UDWORD means unsigned doubleword:

DWORD=K, Loop{UDWORD=UWORD*UWORD, DWORD+=UDWORD}, ScaleAndRound(UDWORD)

Clip and Mask would be numbers in < 0 ; 65535 > range. So the sum before scaling would be in < 0 ; 256^4 )

However bias can be negative also and has no real boundaries, although the only range which makes sense is ( -65536*n ; 65535 > where n is number of input clips.

I don't have code yet but the first idea with the one bit inaccuracy would be something like
move xmm7, dword ptr [bias_scaled] //here I cannot scale 16bits up beacause of the range is outside signed word so I will scale 8bits up

loop:
...
psrlw xmm0, 1 //xmm0 contains packed interleaved clip_i and clip_j where j=i+1
psrlw xmm1, 1 //xmm1 contains packed interleaved mask_i and mask_j (I need to do these operations to clear the sign bit because pmaddwd expects signed arguments
pmaddwd xmm0, xmm1
psrld xmm0, 6 //extra op to get on the scale of K in xmm7
paddd xmm7, xmm0
...
jnz loop
...
psraw xmm7, 8
...
pack with unsigned saturation to word
Now I see there is lot of psr* which I don't believe is very fast op :(

[redfordxx]

Read EDIT2 first!
I am trying other thing when I have problems with this if {#define} and it is to use variable instead of constant. But I got really confused and please correct my error:
These are the definitions:
int scaleI;
scaleI=13;
#define SCALE 13
__asm{
mov eax, SCALE
mov ebx, scaleI
Now, following instructions have different results and I dont know why:psrad xmm0, SCALE //only this one is correct
psrad xmm0, scaleI
psrad xmm0, eax
psrad xmm0, ebx

EDIT: I think the root cause can be found in this asm listing, which...well what's happening there is beyond my understanding:
; 75 : scalepowerA=13;

mov DWORD PTR [edi+16], 13 ; 0000000dH

; 76 :
; 77 : __asm
; 78 : {
; 79 : pushad

pushad

; 82 : mov eax, scalepowerA

mov eax, 16 ; 00000010H

; 110 : psrad xmm0, eax

psrad xmm0, eax

; 111 : psrad xmm1, scalepowerA

psrad xmm1, 16 ; 00000010H
I deleted some lines but nowhere was changed or accessed eax or scalepowerA

EDIT2:
I seem to figure that out, is that so, that inline asm can access only local variables? And when the variable is defined in *.h file, it is silently replaced with number 16...

EDIT3:
But again:
This doesnot work:
mov eax, 13
psrad xmm0, eax
This does:
psrad xmm0, 13

So, shortly:
how can I shift xmm register based on some variable I create in C++ code? (hopefully it won't be too slow, compared to immediate value...)

[IanB]

Again with the keyhole view. If you want help post enough so the full context is available to us.

PSRAD has no Reg32 variant. Only Immediate and MMreg versions.Opcode Instruction Description
0F E2 /r PSRAD mm, mm/m64 Shift doublewords in mm right by mm/m64 while shifting in sign bits.
66 0F E2 /r PSRAD xmm1, xmm2/m128 Shift doubleword in xmm1 right by xmm2 /m128 while shifting in sign bits.
0F 72 /4 ib PSRAD mm, imm8 Shift doublewords in mm right by imm8 while shifting in sign bits.
66 0F 72 /4 ib PSRAD xmm1, imm8 Shift doublewords in xmm1 right by imm8 while shifting in sign bits.

[redfordxx]

Hi, please, how do I change dimensions of the video? I tried:

PVideoFrame __stdcall RAverageM::GetFrame(int n, IScriptEnvironment *env) {...
vi.width = dst_width;
vi.height = dst_height;
return WeightedAverageM16(env, vi);
...}


PVideoFrame RAverageM::WeightedAverageM16(IScriptEnvironment* env, VideoInfo vi)
{
result = env->NewVideoFrame (vi,PitchAlign);

if (y==3) AveragePlaneM16(PLANAR_Y, env);
if (u==3) AveragePlaneM16(PLANAR_U, env);
if (v==3) AveragePlaneM16(PLANAR_V, env);

return (result);
}

This idea I copied from somewhere, but it causes crashes.

[Gavino]

Dimensions should only be changed in the filter's constructor, not in GetFrame(). All frames of a clip are assumed to have the same dimensions.

[Youka]

From Filter SDK documentation (http://avisynth.org/mediawiki/Filter_SDK/Change_frame_size).

[redfordxx]

Thanx guys. I am only not sure, if I have multiple input clips, how to check dimensions of all of them in constructor? I know how to do it in GetFrame function...

[Gavino]

All input clips are available in the constructor (if you pass them in as arguments), so what's the problem?
The GetFrame() function has no more information available than the constructor has.

[SEt]

Writing in assembler means that you know what you are doing, so don't worry about compiler warning if you are sure your code is correct.
If you really want registers you can use all 8 of them:
1) no problems with 6 you mentioned;
2) after you modify ebp in inline assembler you won't be able to access most of your C/C++ variables by name, so load it last and restore at the end;
3) you can even use esp in extreme cases - just need to save it somewhere and restore at the end, don't worry about interrupts - they'll switch to their own stack.

Also it's important to understand that not everything has to be put in registers - for example, putting counters of outer loops to memory is perfectly fine and won't change your program speed in any noticeable way.

[redfordxx]

The GetFrame() function has no more information available than the constructor has.Just to be clear, so I should call GetFrame in the constructor for every input clip, just to know its size?

[cretindesalpes]

No need to call GetFrame:

PClip p = ...;
const VideoInfo & v = p->GetVideoInfo ();
area = v.width * v.height;

[redfordxx]

So, I have 16bit masked average working. But it's a first version and needs some tweaks. I have some ideas but I appreciate any of yours.
The goal is to do masked summation of any number of stacked 16bit clips. However the code looks like the goal is shuffling and packing...and thats my main problem.
There are two main issues:
1] There are SSSE3 and SSE4 instructions
2] Too much shuffling.
As earlier written here I need:
Result=IntegerToWordSaturationSignedToUnsigned[(Sum (WORD*WORD)+INTEGER)>>16]
Result is to be split to hi and lo byte and saved separately on two different parts of the clip. I am really fighting with the journey from __64bit multiplication result through signed sum to separated word into bytes as result.

Here is the code....not reordered yet to keep it more readable.
#include "defsetup.h"

#define PTRSIZE 4
#define DATASTEP (REGSIZE/2)
#define TEMPSTEP 4*PTRSIZE
...
...
_declspec(align(16)) unsigned char shufmask[16];
for (int i=0;i<16;++i)
shufmask[i]=( i<8 ? 2*i : 255 );

int y_cnt=s_height;

__asm
{
pushad

vertical:
mov edi, [negalignwidth]

horizontal:
mov esi, [length16]
MOVALL xm6, qword ptr [rounder]
MOVALL xm7, qword ptr [rounder]

temporal:
mov eax, [rowendptrs+esi] //loading data here
mov ebx, [rowendptrs+esi+4]
mov ecx, [rowendptrs+esi+8]
mov edx, [rowendptrs+esi+12]
MOVHALF xm1, qword ptr [eax+edi]//clip msb
MOVHALF xm0, qword ptr [ebx+edi]//clip lsb
MOVHALF xm3, qword ptr [ecx+edi]//mask msb
MOVHALF xm2, qword ptr [edx+edi]//mask lsb
punpcklbw xm0, xm1 //clip //preparing data here
punpcklbw xm2, xm3 //mask
MOVALL xm1, xm0
MOVALL xm3, xm2
pxor xm4, xm4
punpcklwd xm0, xm4 //clip 0..3
punpcklwd xm2, xm4 //mask 0..3
punpckhwd xm1, xm4 //clip 4..7
punpckhwd xm3, xm4 //mask 4..7
MOVALL xm5, xm0 //clip 0..3
MOVALL xm4, xm1 //clip 4..7
pmuludq xm0, xm2 //result 0, 2 //multiplication 32bits-->64bits
pmuludq xm1, xm3 //result 4, 6
psrlq xm5, 32
psrlq xm2, 32
psrlq xm4, 32
psrlq xm3, 32
pmuludq xm2, xm5 //result 1, 3 //multiplication 32bits-->64bits
pmuludq xm3, xm4 //result 4, 6
pshufd xm0, xm0, 1*0+4*2+16*1+64*3
pshufd xm1, xm1, 1*0+4*2+16*1+64*3
pshufd xm2, xm2, 1*0+4*2+16*1+64*3
pshufd xm3, xm3, 1*0+4*2+16*1+64*3
punpckldq xm0, xm2 //result 0..3
punpckldq xm1, xm3 //result 4..7
psrld xm0, 8 //to have some space for accumulating values
psrld xm1, 8
paddd xm6, xm0 //accomulating value
paddd xm7, xm1
sub esi, TEMPSTEP
jnz temporal
mov eax, [rowendptrs+esi]//result msb ptr
mov ebx, [rowendptrs+esi+4]//result lsb ptr
psrad xm6, 8
psrad xm7, 8
packusdw xm6, xm7 ////////////////// SSE 4.1
MOVALL xm3, qword ptr [shufmask]
MOVALL xm7, xm6
psrlw xm6, 8
pshufb xm7, xm3 //lsb //////////////// SSSE3
pshufb xm6, xm3 //msb
MOVHALF qword ptr [ebx+edi], xm7 //storing data
MOVHALF qword ptr [eax+edi], xm6
add edi, DATASTEP
jnz horizontal
mov ecx, length16P
pitchcycle: //next line, so here increasing pitch
MOVALL xm1, qword ptr [rowendptrs+ecx-TEMPSTEP]
paddd xm1, qword ptr [pitches+ecx-TEMPSTEP]
MOVALL qword ptr [rowendptrs+ecx-TEMPSTEP], xm1
sub ecx, TEMPSTEP
jnz pitchcycle

dec y_cnt
jnz vertical

emms
popad
}
#include "defreset.h"

[IanB]

Spoon feeding :-
movdqa xmm5, xmmwordptr[rounder] // 0.5 << 8
...

temporal:
mov eax, [rowendptrs+esi] //loading data here
mov ebx, [rowendptrs+esi+4]
mov ecx, [rowendptrs+esi+8]
mov edx, [rowendptrs+esi+12]
movq xmm1, qword ptr [eax+edi]//clip msb
movq xmm0, qword ptr [ebx+edi]//clip lsb
movq xmm3, qword ptr [ecx+edi]//mask msb
movq xmm2, qword ptr [edx+edi]//mask lsb
punpcklbw xmm0, xmm1 //clip 7766554433221100
punpcklbw xmm2, xmm3 //mask 7766554433221100
movdqa xmm1, xmm0
pmullw xmm0, xmm2 // l7l6l5l4l3l2l1l0
pmulhuw xmm1, xmm2 // h7h6h5h4h3h2h1h0
movdqa xmm2, xmm0
punpcklwd xmm0, xmm1 // h3l3h2l2h1l1h0l0
punpckhwd xmm2, xmm1 // h7l7h6l6h5l5h4l4
paddd xmm0, xmm5 // += 0.5 << 8
paddd xmm2, xmm5
psrld xmm0, 8 // To have some space for accumulating upto 256 values
psrld xmm2, 8
paddd xmm6, xmm0 //accumulating values
paddd xmm7, xmm2
sub esi, TEMPSTEP
jnz temporal

Or perhaps accumulate to 24bit precision floating point :-
temporal:
mov eax, [rowendptrs+esi] //loading data here
mov ebx, [rowendptrs+esi+4]
mov ecx, [rowendptrs+esi+8]
mov edx, [rowendptrs+esi+12]
movq xmm1, qword ptr [eax+edi]//clip msb
movq xmm0, qword ptr [ebx+edi]//clip lsb
movq xmm3, qword ptr [ecx+edi]//mask msb
movq xmm2, qword ptr [edx+edi]//mask lsb
punpcklbw xmm0, xmm1 //clip 7766554433221100
punpcklbw xmm2, xmm3 //mask 7766554433221100
movdqa xmm1, xmm0
pmullw xmm0, xmm2 // l7l6l5l4l3l2l1l0
pmulhuw xmm1, xmm2 // h7h6h5h4h3h2h1h0
movdqa xmm2, xmm0
punpcklwd xmm0, xmm1 // h3l3h2l2h1l1h0l0
punpckhwd xmm2, xmm1 // h7l7h6l6h5l5h4l4
psrld xmm0, 1 // Zap 32bit unsigned to 31bit signed. Grrr! Damn!
psrld xmm2, 1
cvtdq2ps xmm0, xmm0 // To float (24bit)
cvtdq2ps xmm2, xmm2
addps xmm6, xmm0 //accumulating float values (24bit)
addps xmm7, xmm2
sub esi, TEMPSTEP
jnz temporal

[redfordxx]

Spoon feeding :-
Damn I am stupid, because I used this algo for 8bit version. I probably too much freaked out from our discussion about precission of the 16bit multiplication.

[IanB]

Having just written that small chunk of 16bpp code I am reminded of past thoughts that perhaps we should only be using 15bits, i.e. [0.32767] the positive values of a signed 16 bit integer. 16bit unsigned is a 2nd class citizen in the world of MMX/SSE. Using only the positive values of a signed 16 bit integer just dodge so many bullets.

[redfordxx]

Having just written that small chunk of 16bpp code I am reminded of past thoughts that perhaps we should only be using 15bits, i.e. [0.32767] the positive values of a signed 16 bit integer. 16bit unsigned is a 2nd class citizen in the world of MMX/SSE. Using only the positive values of a signed 16 bit integer just dodge so many bullets.Well I surely have in mind Dither tools, while working on it, and there is 16 already...

[StainlessS]

Hi Redfordxx

Can you clarify, is this an update to Average plugin from mg262 or a new one
(if so, then maybe a differing name might be in order)
http://forum.doom9.org/showthread.php?t=100626

EDIT:- Oops sorry, Just spotted your post, the last one in the mg262 thread,
assume that is is an update. Good, looking forward to it.

[redfordxx]

Writing in assembler means that you know what you are doingNot always my case;)
2) after you modify ebp in inline assembler you won't be able to access most of your C/C++ variables by name, so load it last and restore at the endI have been thinking a little about this and I have question.
When I am doing
mov eax, [var1]
whats really happenning is
mov eax, [ebp+var1offset]?

And
mov eax, [var1+edi+4]
is
mov eax, [ebp+var1offset+edi+4]?

That leads me to question. Is speed of accessing variables in theory dependent on number of arguments in []?

[redfordxx]

Hi Redfordxx

Can you clarify, is this an update to Average plugin from mg262 or a new one
(if so, then maybe a differing name might be in order)
What I am working on will have improved and faster functionalities of those both old ones and something extra.

[StainlessS]

Magic.

[redfordxx]

I have in my code SSSE3 and SSE4 instruction. They probably will get replaced, but for now, how do I detect these processors? I didn't find it in avisynth.h enum.