Welcome to Doom9's Forum, THE in-place to be for everyone interested in DVD conversion.

Before you start posting please read the forum rules. By posting to this forum you agree to abide by the rules.

 

Go Back   Doom9's Forum > Capturing and Editing Video > Avisynth Usage

Reply
 
Thread Tools Search this Thread Display Modes
Old 2nd October 2015, 20:14   #1  |  Link
jmac698
Registered User
 
Join Date: Jan 2006
Posts: 1,867
Find Temporal Statistics

Hi,
I'm a bit stuck on an analysis I want to write.
I have a reference frame, and want to find the min,max, median, and average differences of n frames to it - temporally.

Let's focus on one pixel to understand this. Let's say ref at x=0 and y=0 is r,g,b=(100,100,100). I want to compare this to a clip, with the corresponding pixel (0,0) at frame n. I can compare just two frames with subtract(ref,clip.trim(n,-1)). But I want to compare all the clip's frames to the ref, so now I have n frames with signed numbers (signed offset from subtract). Now I want to find the abs() min, max, avg, median per pixel of all these n frames temporally. In other words, the temporal noise per pixel. The output should be 4 frames representing that statistic for each type.

I could use Stainless clipblend(0) for the avg, I'm not sure how to do a large temporal median, and I don't know how to find temporal min and max. However, a quick mod to clipblend could return min, max,avg easily enough.

Thx

ps this wouldn't be hard in vapoursynth except you can't just read pixels, the raw frame is exposed and you have to know the pixel layouts for each YV12 etc format, that's my biggest barrier now to pixel level processing.

pps
I thought of a quick way to do approximate median for any sized clip, keep a short histogram of the pixels, for each pixel set 2^n (limited to -16,+15) in an int, so your histogram is 32 bits per pixel. To find the median, do a population count of bits set, then find the pop/2th bit set, and that's your median. Add in scaling to bins up to size 16. Example:

at (0,0), red channel only:
ref=128
frame0=120 (-8)
frame1=131 (+3)
frame2=135 (+7)
step1: subtract 112 from each value (which represents -16 from ref):
frame0=8
frame1=19
frame2=23
step2: set those bits in a 32bit int (by ORing)
2^8+2^19+2^23=00000000 10001000 00000001 00000000
population count=3
find floor(3/2)th bit set, or 1
bitscan... bit 19 is set
find 19-16=3
the median is 3
you need population count because the same bit might be set multiple times, could be 1 to n bits actually set.
The median of even samples is the average of the middle two.
You can use scaling to set a bit if it's in any bucket of a range of 16 values.
I think popcnt is asm, I don't know about finding the nth bit.
If the value is out of range of the histogram, set the nearest value.

This is to analyze hot pixels in cameras. They aren't always just stuck on full, but biased relative to normal pixels.

I could turn this into a spatial problem with weaverows, but then how to find min,max etc along columns?

Last edited by jmac698; 2nd October 2015 at 20:53.
jmac698 is offline   Reply With Quote
Old 3rd October 2015, 00:26   #2  |  Link
johnmeyer
Registered User
 
Join Date: Feb 2002
Location: California
Posts: 2,691
You should have started with "This is to analyze hot pixels in cameras." I say that because after having skimmed through the main part of your post, I'm not sure your approach is going to yield anything.

I just finished two 4-camera wedding shoots, and the three cameras that were locked down have video that hardly varies from frame to frame. A hot pixel wouldn't look any different than a sparkle on a non-moving dress. Even if you averaged five minutes of footage from my videos, 98% of the pixels would not change, and therefore temporal averages wouldn't mean much.

I think that you inevitably have to introduce some sort way to only do averaging on pixels that are changing from frame to frame, i.e., where there is significant motion.
johnmeyer is offline   Reply With Quote
Old 3rd October 2015, 00:57   #3  |  Link
Reel.Deel
Registered User
 
Join Date: Mar 2012
Location: Texas
Posts: 1,664
@jmac698

I don't know if this will be of any help but there's the HealDeadPixels plugin that might be worth a look.
Reel.Deel is offline   Reply With Quote
Old 3rd October 2015, 11:03   #4  |  Link
jmac698
Registered User
 
Join Date: Jan 2006
Posts: 1,867
Thanks for the input, however you misunderstand me. This is not for normal pictures, but astrophotography and night photography (with minutes long exposures) where noise becomes very important. It is normal to use 3 kinds of calibration frames; bias frames (those with a very short exposure, to simply read the patterns in the row amplifiers), dark frames (which shows the slight variations in individual pixels, due to manufacturing variations), and flat frames (normal pics of a white wall to correct for vignetting). My research involves going even further, because noise depends on temperature as well, and for this reason, people usually have to take the dark frames at the same time as the normal frames (during long exposures), obviously this more than doubles the time, and you don't have that much time to waste, for example a picture of Jupiter which rotates visibly in under 1 minute, or the moon which can move out of sight as you watch (without a fancy tracking mount).

What I think I can do is pre-make calibration frames, then apply them based on temperature. There's some other things I can investigate too. There's poisson and gaussian noise, but poisson varies based on image intensity AFAIK and has to be denoised differently. I wanted to find the cross over point.

So you see, this is all very much more advanced than you might typically think. Also my camera has no hotpixels at all. There is some bias patterns however.

I have a test here in low light, iso1600, it looks like complete coloured noise, yet after averaging 12 frames, I can clearly see an image.

Another application is even for the day time, with a deep IR filter an unmodified camera, you can still take IR pictures, but they are very dim. In fact, it can replace an ND filter to look at the sun.

Last edited by jmac698; 3rd October 2015 at 11:10.
jmac698 is offline   Reply With Quote
Old 3rd October 2015, 18:00   #5  |  Link
StainlessS
HeartlessS Usurer
 
StainlessS's Avatar
 
Join Date: Dec 2009
Location: Over the rainbow
Posts: 10,980
Jmac, what is maximum required temporal radius, (presumably both forward and backwards) ?

EDIT: ClipBlend() averages frames prior to current (with Delay=0 being all frames).
RGBAmplifier() samples with temporal radius (up to 99 frames, I think).
Either could probably be modified to produce required stats (but you would need to decide what constitutes min or median etc).
__________________
I sometimes post sober.
StainlessS@MediaFire ::: AND/OR ::: StainlessS@SendSpace

"Some infinities are bigger than other infinities", but how many of them are infinitely bigger ???

Last edited by StainlessS; 3rd October 2015 at 18:19.
StainlessS is offline   Reply With Quote
Old 3rd October 2015, 21:26   #6  |  Link
StainlessS
HeartlessS Usurer
 
StainlessS's Avatar
 
Join Date: Dec 2009
Location: Over the rainbow
Posts: 10,980
Jmac,
Gonna attempt to create file based count array, 256*sizeof(unsigned int) for each pixel/channel (RGB only), then you can
apply whatever method you want to extract data (defo need this for median, I think).

EDIT: Will likely not be terribly speedy.
__________________
I sometimes post sober.
StainlessS@MediaFire ::: AND/OR ::: StainlessS@SendSpace

"Some infinities are bigger than other infinities", but how many of them are infinitely bigger ???

Last edited by StainlessS; 3rd October 2015 at 21:38.
StainlessS is offline   Reply With Quote
Old 4th October 2015, 02:42   #7  |  Link
StainlessS
HeartlessS Usurer
 
StainlessS's Avatar
 
Join Date: Dec 2009
Location: Over the rainbow
Posts: 10,980
Jmac, here ClipStats:- http://www.mediafire.com/download/1y...0-20151004.zip

EDIT: Above link edited (got filename month and day back to front).

Slow as hell (about 20 secs per frame [after initial file created, on 2.4Ghz Core Duo]).

Code:
/*
    This program is free software; you can redistribute it and/or modify
    it under the terms of the GNU General Public License as published by
    the Free Software Foundation.

    This program is distributed in the hope that it will be useful,
    but WITHOUT ANY WARRANTY; without even the implied warranty of
    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
    GNU General Public License for more details.

    You should have received a copy of the GNU General Public License
    along with this program; if not, write to the Free Software
    Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.

*/

#include "Compiler.h"

#include <windows.h>
#include <stdio.h>
#include <math.h>

#include "avisynth.h"

#define BUG

#ifdef BUG
    int dprintf(char* fmt, ...) {
        char printString[512]="";                   // MUST NUL TERM THIS
        char *p=printString;
        va_list argp;
        va_start(argp, fmt);
        vsprintf(p, fmt, argp);
        va_end(argp);
        for(;*p++;);
        --p;                                        // @ null term
        if(printString == p || p[-1] != '\n') {
            p[0]='\n';                              // append n/l if not there already
            p[1]='\0';
        }
        OutputDebugString(printString);
        return p-printString;                       // strlen printString   
    }

    // OLD C89 Compiler: (eg VS6)
    // Call as eg DPRINTF(("Forty Two = %d",42))        // Enclosed in double parentethis, No Semi colon
    #define   DPRINTF(x)      dprintf x;

    // C99 Compiler (eg VS2008)
    // #define DPRINTF(fmt, ...)   dprintf(fmt, ##__VA_ARGS__);       // No Semi colon needed
#else
    // OLD C89 Compiler: (eg VS6)
        #define DPRINTF(x);             /* x */

    //  #define DPRINTF(fmt,...);       /* fmt */
#endif


class ClipStats : public GenericVideoFilter {
    const int  delay;
    const int  mode;
    const double thresh;
    const int lo;
    const int hi;
    const char *fn;
    //
    FILE *fp;
    int num_frames;
    int Els_0;
    int h_y;
    int w_y;
    int FramesDone;
    int last;

    enum {MODE_MIN=0,MODE_MAX,MODE_MINMAX,MODE_MEDIAN,MODE_AVG,MODE_STDEV,MODE_IRNG,MODE_LIMIT=MODE_IRNG};

    void Reset(void);
    void Add(int n, IScriptEnvironment* env);
    void Sub(int n, IScriptEnvironment* env);
    PVideoFrame Render(int n, IScriptEnvironment* env);

public:
    ClipStats(PClip _child,int _delay,int _mode,double _thresh,int _lo,int _hi,const char *_fn,IScriptEnvironment* env);
    ~ClipStats();
    PVideoFrame __stdcall GetFrame(int n, IScriptEnvironment* env);
};

ClipStats::~ClipStats() {
    if(fp!=NULL)    {fclose(fp); fp=NULL;};
}

ClipStats::ClipStats(PClip _child,int _delay,int _mode,double _thresh,int _lo,int _hi,const char *_fn,IScriptEnvironment* env) : 
    GenericVideoFilter(_child),delay(_delay),mode(_mode),thresh(_thresh),lo(_lo),hi(_hi),fn(_fn),fp(NULL) {


    if(vi.width==0 || vi.num_frames<=0)         env->ThrowError("ClipStats: Clip has no video\n");
    if(!vi.IsRGB24())                           env->ThrowError("ClipStats: RGB24 Only\n");
    if(delay<0)                                 env->ThrowError("ClipStats: +ve delay only please.\n");
    if(mode<0||mode>MODE_LIMIT)                 env->ThrowError("ClipStats: 0 <= Mode <= %d.\n",MODE_LIMIT);
    if(thresh<0.0||thresh>100.0)                env->ThrowError("ClipStats: 0.0 <= Thresh <= 100.0\n");
    if(lo<0||lo>255)                            env->ThrowError("ClipStats: 0 <= Lo <= 255.\n");
    if(hi<lo||hi>255)                           env->ThrowError("ClipStats: Lo <= Hi <= 255.\n");

    num_frames      = vi.num_frames;
    PVideoFrame src = child->GetFrame(0, env);                  // get frame 0
    h_y             = src->GetHeight();
    w_y             = src->GetRowSize();
    Els_0       = h_y  * w_y;
    if((fp=fopen(fn, "wb+" ))==NULL)            env->ThrowError("ClipStats: Cannot open data file %s",fn);
    FramesDone=0;
    last = -1;
}


void ClipStats::Reset(void) {
    unsigned int bf[256];
    int i;
    for(i=256;--i>=0;bf[i]=0);
    rewind(fp);
    for(i=Els_0;--i>=0;) {
        int wr=fwrite(bf,sizeof(bf),1,fp);          // should be return 1
    }
    FramesDone=0;
}

void ClipStats::Add(int n,IScriptEnvironment* env) {
    PVideoFrame src = child->GetFrame(n, env);
    const int pitch = src->GetPitch();
    const BYTE * srcp = src->GetReadPtr();
    ++FramesDone;                                               // counting this frame
    int y,x;
    __int64 pos;
    for(pos=0,y=0;y<h_y;++y) {
        for(x=0;x<w_y;++x) {
            __int64 off = (pos*256 + srcp[x]) * sizeof(unsigned int);
            _fseeki64(fp, off,SEEK_SET);
            unsigned int val;
            fread(&val,sizeof(val),1,fp);
            _fseeki64(fp,off,SEEK_SET);
            ++val;
            fwrite(&val,sizeof(val),1,fp);
            ++pos;
        }
        srcp += pitch;
    }
}

void ClipStats::Sub(int n,IScriptEnvironment* env) {
    PVideoFrame src = child->GetFrame(n, env);
    const int pitch = src->GetPitch();
    const BYTE * srcp = src->GetReadPtr();
    --FramesDone;                                               // counting this frame
    int y,x;
    __int64 pos;
    for(pos=0,y=0;y<h_y;++y) {
        for(x=0;x<w_y;++x) {
            __int64 off = (pos*256 + srcp[x]) * sizeof(unsigned int);
            _fseeki64(fp, off,SEEK_SET);
            unsigned int val;
            fread(&val,sizeof(val),1,fp);
            _fseeki64(fp, off,SEEK_SET);
            --val;
            fwrite(&val,sizeof(val),1,fp);
            ++pos;
        }
        srcp += pitch;
    }
}


PVideoFrame ClipStats::Render(int n,IScriptEnvironment* env) {
    PVideoFrame dst =   env->NewVideoFrame(vi);
    const int dpitch = dst->GetPitch();
    BYTE *  dstp= dst->GetWritePtr();
    int y,x,i;
    unsigned int cnt[256];
    rewind(fp);
    switch(mode) {
    case MODE_MIN:
        {
            const unsigned int lim = (unsigned int)(FramesDone * thresh);
            for(y=0;y<h_y;++y) {                                    
                for(x=0;x<w_y;++x) {
                    fread(cnt,sizeof(cnt),1,fp);
                    int rmin=255;
                    unsigned int sm=0;
                    for(i=0;i<256;++i) {
                        sm += cnt[i];
                        if(sm> lim) {
                            rmin=i;
                            break;
                        }
                    }
                    dstp[x]=rmin;
                }
                dstp += dpitch;
            }
        }
        break;
    case MODE_MAX:
        {
            const unsigned int lim = (unsigned int)(FramesDone * thresh);
            for(y=0;y<h_y;++y) {                                    
                for(x=0;x<w_y;++x) {
                    fread(cnt,sizeof(cnt),1,fp);
                    int rmax=0;
                    unsigned int sm=0;
                    for(i=256;--i>=0;) {
                        sm += cnt[i];
                        if(sm > lim) {
                            rmax=i;
                            break;
                        }
                    }
                    dstp[x]=rmax;
                }
                dstp += dpitch;
            }
        }
        break;
    case MODE_MINMAX:
        {
            const unsigned int lim = (unsigned int)(FramesDone * thresh);
            for(y=0;y<h_y;++y) {                                    
                for(x=0;x<w_y;++x) {
                    fread(cnt,sizeof(cnt),1,fp);
                    int rmin=255;
                    unsigned int sm=0;
                    for(i=0;i<256;++i) {
                        sm += cnt[i];
                        if(sm> lim) {
                            rmin=i;
                            break;
                        }
                    }
                    int rmax=0;
                    sm=0;
                    for(i=256;--i>=0;) {
                        sm += cnt[i];
                        if(sm > lim) {
                            rmax=i;
                            break;
                        }
                    }
                    dstp[x]=rmax-rmin;
                }
                dstp += dpitch;
            }
        }
        break;
    case MODE_MEDIAN:
        {
            const unsigned int lim = FramesDone / 2;
            for(y=0;y<h_y;++y) {                                    
                for(x=0;x<w_y;++x) {
                    fread(cnt,sizeof(cnt),1,fp);
                    int rmed=0;
                    unsigned int sm=0;
                    for(i=256;--i>=0;) {
                        sm += cnt[i];
                        if(sm > lim) {
                            rmed=i;
                            break;
                        }
                    }
                    dstp[x]=rmed;
                }
                dstp += dpitch;
            }
        }
        break;
    case MODE_AVG:
        for(y=0;y<h_y;++y) {                                    
            for(x=0;x<w_y;++x) {
                fread(cnt,sizeof(cnt),1,fp);
                double mean=0.0;
                __int64 acc      = 0;
                for(i=256;--i>=0;) {
                    acc += cnt[i] * i;
                }
                mean = (double)acc / FramesDone;
                dstp[x]=int(mean+0.5);
            }
            dstp += dpitch;
        }
        break;
    case MODE_STDEV:
        for(y=0;y<h_y;++y) {                                    
            for(x=0;x<w_y;++x) {
                fread(cnt,sizeof(cnt),1,fp);
                double mean=0.0;
                double std=0.0;
                if(FramesDone>1) {
                    __int64 acc      = 0;
                    for(i=256;--i>=0;) {
                        acc += cnt[i] * i;
                    }
                    mean = (double)acc / FramesDone;
                    double Sum_DiffSquared = 0.0;
                    for (i=256;--i>=0;) {
                        const unsigned int count=cnt[i];
                        if(count) {
                            const double diff = (double)i - mean;               // x - xbar
                            Sum_DiffSquared += (diff * diff * double(count));   // Sigma((x-xbar)^2) 
                        }
                    }
                    std=(Sum_DiffSquared<=0.0)?0.0:sqrt(Sum_DiffSquared / FramesDone);
                }
                dstp[x]=int(std+0.5);
            }
            dstp += dpitch;
        }
        break;
    case MODE_IRNG:
        for(y=0;y<h_y;++y) {                                    
            for(x=0;x<w_y;++x) {
                fread(cnt,sizeof(cnt),1,fp);
                unsigned int sm=0;
                for(i=lo;i<=hi;++i) {
                    sm += cnt[i];
                }
                dstp[x]= int((double(sm) / FramesDone)+0.5);
            }
            dstp += dpitch;
        }
        break;
    }
    return dst;
}




PVideoFrame __stdcall ClipStats::GetFrame(int n, IScriptEnvironment* env) {
    n = (n<0) ? 0 : (n>= num_frames) ? num_frames - 1 : n;

    PVideoFrame dst;

DPRINTF(("ClipStats: %d] Last=%d FramesDone=%d",n,last,FramesDone))

    if(n != 0 && n == last) {                                   // Assume cache failure

DPRINTF(("ClipStats: %d] Cache Fail (repeat request for this frame)",n))
        dst = Render(n,env);
    } else {
        int WantS = (delay==0) ? 0 : max(n - delay,0); 
        int WantE = n;
        int GotS  = (FramesDone>0) ? max(last - FramesDone + 1,0) : - 1;
        int GotE  = (FramesDone>0) ? last : - 1;
        int ResCnt = WantE-WantS+1;                                                 // frames to read if accumulator reset
        int AccCnt= (n<last) ? GotE-WantE + GotS-WantS : WantE-GotE + WantS-GotS;   // frames to read if subtract/add accumlator

DPRINTF(("ClipStats: GotS=%d GotE=%d WantS=%d WantE=%d",GotS,GotE,WantS,WantE))
        // Frame 0 or jump about, restart:clear Accumulators
        if(n==0 || FramesDone<=0 || ResCnt<AccCnt) {
DPRINTF(("ClipStats: %d] RESET (last=%d ResCnt=%d AccCnt=%d)",n,last,ResCnt,AccCnt))
            Reset();
            GotS = - 1;
            GotE = - 1;
        } else if(n > last) {
            while(GotS<WantS) { // Deduct oldest frames dat from Accumulator
DPRINTF(("ClipStats: DEDUCT_S: %d",GotS))
                Sub(GotS,env);
                ++GotS;
            }
        } else {
            int m=WantE+1;
            while(m<=GotE) {
DPRINTF(("ClipStats: DEDUCT_E: %d",m))
                Sub(n,env);
                ++m;
            }
            GotE=WantE;
        }

        if(WantS<GotS) {
            int m=WantS;
            while(m<GotS) {
DPRINTF(("ClipStats: ADD_S: %d (GotS=%d WantS=%d",m,GotS,WantS))
                Add(m,env);
                ++m;
            }
            GotS=WantS;
        }

        if(GotE<0) {
            GotS = WantS ;
            GotE = WantS - 1;
        }
        if(GotE < WantE) {
            while(GotE<WantE) {
                ++GotE;
                if(GotE!=WantE) {
DPRINTF(("ClipStats: ADD_E: %d",GotE))
                    Add(GotE,env);
                } else {
DPRINTF(("ClipStats: ADD_RENDER: %d",GotE))
                    Add(GotE,env);
                    dst=Render(GotE,env);
                }
            }
        } else {

DPRINTF(("ClipStats: RENDER: %d",n))
            dst=Render(n,env);
        }
        last = n;
    }
DPRINTF(("ClipStats:"))
    return dst;
}


AVSValue __cdecl Create_ClipStats(AVSValue args, void* user_data, IScriptEnvironment* env) {
    int lo = args[4].AsInt(128);                            // lo
    AVSValue ret = new ClipStats(
        args[0].AsClip(),
        args[1].AsInt(0),                           // delay
        args[2].AsInt(0),                           // Mode
        args[3].AsFloat(0.0),                       // Thresh
        lo,                                         // lo
        args[5].AsInt(lo),                          // hi
        args[6].AsString("ClipStats.Data"),         // Filename
        env); 
    return ret;
}

extern "C" __declspec(dllexport) const char* __stdcall AvisynthPluginInit2(IScriptEnvironment* env) {
    env->AddFunction("ClipStats",  "c[delay]i[Mode]i[Thresh]f[Lo]i[Hi]i[FN]s", Create_ClipStats, 0);

    return "`ClipStats' ClipStats plugin";
    // A freeform name of the plugin.
}
Test avs

Code:
avisource("E:\V\JurassicPark.avi")
ConvertToRGB24
Trim(3000,0)
/*
    RGB24 Only.
    
    Clipstats(clip c,int "delay"=0,int "Mode"=0,Float "Thresh"=0.0,Int "lo"=128,Int "Hi"=lo,String "FN"="ClipStats.Data")
    Delay:- as in ClipBlend (0=all frames).
    Mode:- 0=Min, 1=Max, 2=MinMaxDifference, 3=Median, 4=Avg, 5=Stdev, 6=Inrange
    Thresh:- As In YPlaneMin(Threshold) [but of course for R,G,B]
    Lo:- As in RT_YInRange.
    Hi:- As in RT_YInRange.
     
*/


A=ClipStats(delay=0,mode=4)
B=ClipBlend()

StackHorizontal(A,B)
Little testing done, [only comparison with Clipblend(delay=0)]

EDIT: Supplied with VS 2008 build files.

EDIT:
Quote:
I have a reference frame, and want to find the min,max, median, and average differences of n frames to it - temporally.
Oh dear, I should not skim read, not exactly as you require, but could be modded to suit now that its working.
__________________
I sometimes post sober.
StainlessS@MediaFire ::: AND/OR ::: StainlessS@SendSpace

"Some infinities are bigger than other infinities", but how many of them are infinitely bigger ???

Last edited by StainlessS; 4th October 2015 at 05:13.
StainlessS is offline   Reply With Quote
Old 4th October 2015, 05:23   #8  |  Link
StainlessS
HeartlessS Usurer
 
StainlessS's Avatar
 
Join Date: Dec 2009
Location: Over the rainbow
Posts: 10,980
Take 4, http://www.mediafire.com/download/cj...4-20151005.zip

Minus GPL, headers, debug stuff due to size.
Code:
class ClipDiffs : public GenericVideoFilter {
    PClip child2;
    const int  delay;
    const int  mode;
    const double thresh;
    const int lo;
    const int hi;
    const char *fn;
    //
    FILE *fp;
    int num_frames;
    int Els;
    int hit;
    int wid;
    int FramesDone;
    int last;
    BYTE Diff[513]; 

    enum {MODE_MIN=0,MODE_MAX,MODE_MINMAX,MODE_MEDIAN,MODE_AVG,MODE_STDEV,MODE_IRNG,MODE_LIMIT=MODE_IRNG};
    
    enum{CNT_BFSZ=512}; // Last one not used [-255 -> 255 == 511]
    enum{CNT_REL =255};

    void Reset(void);
    void Add(int n, IScriptEnvironment* env);
    void Sub(int n, IScriptEnvironment* env);
    PVideoFrame Render(int n, IScriptEnvironment* env);

public:
    ClipDiffs(PClip _child,PClip _child2,int _delay,int _mode,double _thresh,int _lo,int _hi,const char *_fn,IScriptEnvironment* env);
    ~ClipDiffs();
    PVideoFrame __stdcall GetFrame(int n, IScriptEnvironment* env);
};

ClipDiffs::~ClipDiffs() {
    if(fp!=NULL)    {fclose(fp); fp=NULL;};
}

ClipDiffs::ClipDiffs(PClip _child,PClip _child2,int _delay,int _mode,double _thresh,int _lo,int _hi,const char *_fn,IScriptEnvironment* env) : 
    GenericVideoFilter(_child),child2(_child2),delay(_delay),mode(_mode),thresh(_thresh),lo(_lo),hi(_hi),fn(_fn),fp(NULL) {
    const char *myName="ClipDiffs: ";

    if(vi.width==0 || vi.num_frames<=0)         env->ThrowError("%sClip has no video",myName);
    if(!vi.IsRGB24())                           env->ThrowError("%sRGB24 Only",myName);
    const VideoInfo &vi2=child2->GetVideoInfo();
    if(!vi2.IsRGB24())                          env->ThrowError("%sclip 2 RGB24 Only",myName);
    if(vi2.width==0 || vi2.num_frames<=0)       env->ThrowError("%sclip 2 has no video",myName);

    if(vi.height != vi2.height)                 env->ThrowError("%sclip 2 height mismatch",myName);
    if(vi.width != vi2.width)                   env->ThrowError("%sclip 2 width mismatch",myName);

    if(delay<0)                                 env->ThrowError("%s+ve delay only please.",myName);
    if(mode<0||mode>MODE_LIMIT)                 env->ThrowError("%s0 <= Mode <= %d.",myName,MODE_LIMIT);
    if(thresh<0.0||thresh>100.0)                env->ThrowError("%s0.0 <= Thresh <= 100.0",myName);
    if(lo< -255||lo>255)                        env->ThrowError("%s-255 <= Lo <= 255.",myName);
    if(hi<lo||hi>255)                           env->ThrowError("%sLo <= Hi <= 255.",myName);

    num_frames      = vi.num_frames;
    PVideoFrame src = child->GetFrame(0, env);                  // get frame 0
    hit             = src->GetHeight();
    wid             = src->GetRowSize();
    Els             = hit  * wid;
    if((fp=fopen(fn, "wb+" ))==NULL)            env->ThrowError("%sCannot open data file %s",myName,fn);
    int i;
    for (i=0; i<=512; i++) Diff[i] = max(0,min(255,i-129));
    FramesDone=0;
    last = -1;
}

void ClipDiffs::Reset(void) {
    int bf[CNT_BFSZ];
    int i;
    rewind(fp);
    const int bfsz = 512*sizeof(bf);                    // 1MB
    BYTE *bfp = new BYTE[bfsz];
    if(bfp==NULL) {
        memset(bf,0,sizeof(bf));
        for(i=Els;--i>=0;) {
            int wr=fwrite(bf,sizeof(bf),1,fp);          // should be return 1
        }
    } else {
        memset(bfp,0,bfsz);
        __int64 clrbytes = __int64(sizeof(bf)) * Els;
        int sods = int(clrbytes  / bfsz);
        int odds = int(clrbytes  % bfsz);
        int wr;
        for(wr=1,i = sods; wr == 1 && --i>=0; ) {
            wr=fwrite(bfp,bfsz,1,fp);
        }
        if(wr==1 && odds>0) {
            wr=fwrite(bfp,odds,1,fp);
        }
        delete [] bfp;
    }
    FramesDone=0;
}

void ClipDiffs::Add(int n,IScriptEnvironment* env) {
    PVideoFrame src = child->GetFrame(n, env);
    const int pitch = src->GetPitch();
    const BYTE * srcp = src->GetReadPtr();
    PVideoFrame ref  = child2->GetFrame(0,env);
    const int pitchr = ref->GetPitch();
    const BYTE * refp= ref->GetReadPtr();
    ++FramesDone;                                               // counting this frame
    int y,x;
    __int64 pos;
    for(pos=0,y=0;y<hit;++y) {
        for(x=0;x<wid;++x) {
            int dif = srcp[x] - refp[x];                        // -255 -> 255
            __int64 off = (pos*CNT_BFSZ + dif + CNT_REL) * sizeof(int);
            _fseeki64(fp, off,SEEK_SET);
            int val;
            fread(&val,sizeof(val),1,fp);
            _fseeki64(fp,off,SEEK_SET);
            ++val;
            fwrite(&val,sizeof(val),1,fp);
            ++pos;
        }
        srcp += pitch;
        refp += pitchr;
    }
}

void ClipDiffs::Sub(int n,IScriptEnvironment* env) {
    PVideoFrame src = child->GetFrame(n, env);
    const int pitch = src->GetPitch();
    const BYTE * srcp = src->GetReadPtr();
    PVideoFrame ref  = child2->GetFrame(0,env);
    const int pitchr = ref->GetPitch();
    const BYTE * refp= ref->GetReadPtr();
    --FramesDone;                                               // counting this frame
    int y,x;
    __int64 pos;
    for(pos=0,y=0;y<hit;++y) {
        for(x=0;x<wid;++x) {
            int dif = srcp[x] - refp[x];
            __int64 off = (pos*CNT_BFSZ + dif + CNT_REL) * sizeof(int);
            _fseeki64(fp, off,SEEK_SET);
            int val;
            fread(&val,sizeof(val),1,fp);
            _fseeki64(fp, off,SEEK_SET);
            --val;
            fwrite(&val,sizeof(val),1,fp);
            ++pos;
        }
        srcp += pitch;
        refp += pitchr;
    }
}


PVideoFrame ClipDiffs::Render(int n,IScriptEnvironment* env) {
    PVideoFrame dst =   env->NewVideoFrame(vi);
    const int dpitch = dst->GetPitch();
    BYTE *  dstp= dst->GetWritePtr();
    int y,x,i;
    int cnt[CNT_BFSZ];
    rewind(fp);
    switch(mode) {
    case MODE_MIN:
        {
            const int lim = (int)(FramesDone * thresh);
            for(y=0;y<hit;++y) {                                    
                for(x=0;x<wid;++x) {
                    fread(cnt,sizeof(cnt),1,fp);
                    int rmin=255;
                    int sm=0;
                    for(i=0;i<511;++i) {
                        sm += cnt[i];
                        if(sm> lim) {
                            rmin=i-CNT_REL;
                            break;
                        }
                    }
                    dstp[x]=Diff[rmin + 128 + 129];         // Min diff can be +ve or -ve, center @ 128
                }
                dstp += dpitch;
            }
        }
        break;
    case MODE_MAX:
        {
            const int lim = (int)(FramesDone * thresh);
            for(y=0;y<hit;++y) {                                    
                for(x=0;x<wid;++x) {
                    fread(cnt,sizeof(cnt),1,fp);
                    int rmax=-255;
                    int sm=0;
                    for(i=511; --i>= 0;) {
                        sm += cnt[i];
                        if(sm > lim) {
                            rmax=i-CNT_REL;
                            break;
                        }
                    }
                    dstp[x]=Diff[rmax + 128 + 129];         // Max diff can be +ve or -ve, center @ 128
                }
                dstp += dpitch;
            }
        }
        break;
    case MODE_MINMAX:
        {
            const int lim = (int)(FramesDone * thresh);
            for(y=0;y<hit;++y) {                                    
                for(x=0;x<wid;++x) {
                    fread(cnt,sizeof(cnt),1,fp);
                    int rmin=255;
                    int sm=0;
                    for(i=0;i<511;++i) {
                        sm += cnt[i];
                        if(sm> lim) {
                            rmin=i-CNT_REL;
                            break;
                        }
                    }
                    int rmax= -255;
                    sm=0;
                    for(i=511; --i>= 0;) {
                        sm += cnt[i];
                        if(sm > lim) {
                            rmax=i-CNT_REL;
                            break;
                        }
                    }
                    int minmxdif = rmax-rmin;                   // 0 -> 510 (max, 255 - -255)
                    minmxdif=min(minmxdif,512-129);             // Dont stray outside Diff[] due to no add 128
                    dstp[x]=Diff[minmxdif + 129];               // Dont add 128, center a 0 not 128
                }
                dstp += dpitch;
            }
        }
        break;
    case MODE_MEDIAN:
        {
            const int lim = FramesDone / 2;
            for(y=0;y<hit;++y) {                                    
                for(x=0;x<wid;++x) {
                    fread(cnt,sizeof(cnt),1,fp);
                    int rmed=-255;
                    int sm=0;
                    for(i=511; --i>= 0;) {
                        sm += cnt[i];
                        if(sm > lim) {
                            rmed=i-CNT_REL;
                            break;
                        }
                    }
                    dstp[x]=Diff[rmed + 128 + 129];             // Median diff can be +ve or -ve, center @ 128
                }
                dstp += dpitch;
            }
        }
        break;
    case MODE_AVG:
        for(y=0;y<hit;++y) {                                    
            for(x=0;x<wid;++x) {
                fread(cnt,sizeof(cnt),1,fp);
                double mean=0.0;
                __int64 acc= 0;
                for(i=511; --i>= 0;) {
                    acc += cnt[i] * (i-CNT_REL);
                }
                mean = double(acc) / FramesDone;
                dstp[x]=Diff[int(mean+0.5) + 128 + 129];        // Avg diff can be +ve or -ve, center @ 128
            }
            dstp += dpitch;
        }
        break;
    case MODE_STDEV:
        for(y=0;y<hit;++y) {                                    
            for(x=0;x<wid;++x) {
                fread(cnt,sizeof(cnt),1,fp);
                double mean=0.0;
                double std=0.0;
                if(FramesDone>1) {
                    __int64 acc = 0;
                    for(i=511; --i>= 0;) {
                        acc += cnt[i] * (i-CNT_REL);
                    }
                    mean = (double)acc / FramesDone;
                    double Sum_DiffSquared = 0.0;
                    for (i=511; --i>= 0;) {
                        const int count=cnt[i];
                        if(count) {
                            const double diff = ((double)i-CNT_REL) - mean;     // x - xbar
                            Sum_DiffSquared += (diff * diff * double(count));   // Sigma((x-xbar)^2) 
                        }
                    }
                    std=(Sum_DiffSquared<=0.0)?0.0:sqrt(Sum_DiffSquared / FramesDone);
                }
                int istd = int(std+0.5);
                istd=min(istd,512-129);                 // Dont stray outside Diff[] due to no add 128
                dstp[x]=Diff[istd + 129];               // Std +ve, Dont add 128, center a 0 not 128
            }
            dstp += dpitch;
        }
        break;
    case MODE_IRNG:
        for(y=0;y<hit;++y) {                                    
            for(x=0;x<wid;++x) {
                fread(cnt,sizeof(cnt),1,fp);
                int sm=0;
                for(i=lo;i<=hi;++i) {
                    sm += cnt[i+CNT_REL];
                }
                int ir = int((sm * 255.0 / FramesDone)+0.5);
                dstp[x]= ir;                                    // InRange 0 -> 255
            }
            dstp += dpitch;
        }
        break;
    }
    return dst;
}




PVideoFrame __stdcall ClipDiffs::GetFrame(int n, IScriptEnvironment* env) {
    n = (n<0) ? 0 : (n>= num_frames) ? num_frames - 1 : n;
    PVideoFrame dst;
DPRINTF("%d] Last=%d FramesDone=%d",n,last,FramesDone)
    if(n != 0 && n == last) {                                   // Assume cache failure
DPRINTF("%d] Cache Fail (repeat request for this frame)",n)
    } else {
        int WantS = (delay==0) ? 0 : max(n - delay,0); 
        int WantE = n;
        int GotS  = (FramesDone>0) ? max(last - FramesDone + 1,0) : - 1;
        int GotE  = (FramesDone>0) ? last : - 1;
        int ResCnt = WantE-WantS+1;                                                 // frames to read if accumulator reset
        int AccCnt= (n<last) ? GotE-WantE + GotS-WantS : WantE-GotE + WantS-GotS;   // frames to read if subtract/add accumlator

DPRINTF("%d] GotS=%d GotE=%d WantS=%d WantE=%d",n,GotS,GotE,WantS,WantE)
        // Frame 0 or jump about, restart:clear Accumulators
        if(n==0 || FramesDone<=0 || ResCnt<AccCnt) {
DPRINTF("%d] RESET (last=%d ResCnt=%d AccCnt=%d)",n,last,ResCnt,AccCnt)
            Reset();
            GotS = - 1;
            GotE = - 1;
        } else if(n > last) {
            while(GotS<WantS) { // Deduct oldest frames dat from Accumulator
DPRINTF("%d] DEDUCT_S: %d",n,GotS)
                Sub(GotS,env);
                ++GotS;
            }
        } else {
            int m=WantE+1;
            while(m<=GotE) {
DPRINTF("%d] DEDUCT_E: %d",n,m)
                Sub(m,env);
                ++m;
            }
            GotE=WantE;
        }

        if(WantS<GotS) {
            int m=WantS;
            while(m<GotS) {
DPRINTF("%d] ADD_S: %d (GotS=%d WantS=%d",n,m,GotS,WantS)
                Add(m,env);
                ++m;
            }
            GotS=WantS;
        }

        if(GotE<0) {
            GotE = WantS - 1;
        }
        while(GotE<WantE) {
            ++GotE;
DPRINTF("%d] ADD_E: %d",n,GotE)
            Add(GotE,env);
        }
    }
DPRINTF("%d] RENDER: %d",n,n)
    dst=Render(n,env);
    last = n;
DPRINTF("")
    return dst;
}


AVSValue __cdecl Create_ClipDiffs(AVSValue args, void* user_data, IScriptEnvironment* env) {
    int lo = args[5].AsInt(0);                      // lo
    AVSValue ret = new ClipDiffs(
        args[0].AsClip(),
        args[1].AsClip(),
        args[2].AsInt(0),                           // delay
        args[3].AsInt(0),                           // Mode
        args[4].AsFloat(0.0),                       // Thresh
        lo,                                         // lo
        args[6].AsInt(lo),                          // hi
        args[7].AsString("ClipDiffs.Data"),         // Filename
        env); 
    return ret;
}

extern "C" __declspec(dllexport) const char* __stdcall AvisynthPluginInit2(IScriptEnvironment* env) {
    env->AddFunction("ClipDiffs",  "cc[delay]i[Mode]i[Thresh]f[Lo]i[Hi]i[FN]s", Create_ClipDiffs, 0);

    return "`ClipDiffs' ClipDiffs plugin";
    // A freeform name of the plugin.
}
__________________
I sometimes post sober.
StainlessS@MediaFire ::: AND/OR ::: StainlessS@SendSpace

"Some infinities are bigger than other infinities", but how many of them are infinitely bigger ???

Last edited by StainlessS; 5th October 2015 at 21:48. Reason: Update
StainlessS is offline   Reply With Quote
Old 4th October 2015, 23:55   #9  |  Link
StainlessS
HeartlessS Usurer
 
StainlessS's Avatar
 
Join Date: Dec 2009
Location: Over the rainbow
Posts: 10,980
Above post updated.

test client
Code:
avisource("E:\V\JurassicPark.avi").ConvertToRGB24
BilinearResize(320,240)     # reduced size (quicker)

Trim(3000,-10)    # 10 frames
ORG=Last
/*
    RGB24 Only.
    
    ClipDiffs(clip c,clip Ref,int "delay"=0,int "Mode"=0,Float "Thresh"=0.0,Int "lo"=0,Int "Hi"=lo,String "FN"="ClipDiffs.Data")
    Ref:- Reference clip, frame 0 ONLY used.
    Delay:- as in ClipBlend (0=all frames up to current_frame, other = current_frame + previous 'Delay' frames).
    Mode:- 0=Min, 1=Max, 2=MinMaxDifference, 3=Median, 4=Avg, 5=Stdev, 6=InRange
    Thresh:- (0->100%) As In YPlaneMin(Threshold) [but of course for R,G,B] Modes 0, 1, and 2 only.
        A value of Thresh=1.0% for eg mode=MIN would ignore up to 1% of the most -ve differences (or least +ve), ie noise.  
    Lo:- Similar-ish to Lo in RT_YInRange  ( -255 <= Lo <= 255). [NOTE, min is -255, as eg 0 - 255 = -255].
    Hi:- Similar-ish to Hi in RT_YInRange. ( Lo   <= Hi <= 255). Lo and Hi Mode 6 only.
      Lo and Hi set the difference range that will be detected.
    
    ClipDiff modes except for modes 2 (MinMaxDifference), 5(Stdev), and 6 (InRange) return frames where 128 represents no difference
      [as in Subtract()].
    
    Mode 2  MinMaxDifference, returns frames relative to 0 (no difference), minmaxdifference differences greater than 255 will be 
      limited to 255.
    
    Mode 5  Stdev, returns frames relative to 0 (no stdev in difference), stdev differences greater than 255 will be limited to 255.
    EDIT: Above in blue, Dont think can be greater than 255
    
    Mode 6 (InRange) returns frames where 255 represents 100% population of sampled frames have difference in range lo -> hi.
      So eg, ClipDiffs(Last,Ref,delay=0,mode=6,lo= -255,hi= -1) set pixel/channels to 255 where all pixels/channels (from frames sampled so far)
      are less than the same pixel/channel reference frame. Where Lo=0 and Hi=0, would show 255 all pixels (from frames sampled so far) that are
      exactly the same as refernce frame pixels.
      There is no separate Lo/Hi for individual R,G,B channels, perhaps it would be better to move this over to Planar Y only.

    The plugin is compiled with debugging enabled, get DebugView to see what it is doing that is taking so much time :(      
     
*/

Ref = ClipBlend(delay=0).Trim(Framecount-1,-1)          # Ref is Average of all frames (Requires ClipBlend plugin)
#return Ref
ClipDiffs(Last,Ref,delay=0,mode=4,lo=0,hi=0,FN="ClipDiffs.Data")

T=StackHorizontal(Ref,Last)

Return T

# Show Amplified difference compared with Clipblend(0) ONLY INTENDED for Mode=4,delay=0 : should be same at last frame.
D=Levels(128-32,1.0,128+32,0,255)
B=StackHorizontal(ORG,D)
StackVertical(T,B)
Plugin supplied with VS 2008 build files.'
__________________
I sometimes post sober.
StainlessS@MediaFire ::: AND/OR ::: StainlessS@SendSpace

"Some infinities are bigger than other infinities", but how many of them are infinitely bigger ???

Last edited by StainlessS; 5th October 2015 at 22:14. Reason: Updated
StainlessS is offline   Reply With Quote
Old 5th October 2015, 01:07   #10  |  Link
jmac698
Registered User
 
Join Date: Jan 2006
Posts: 1,867
Thanks so much! I don't have time to try it right now but I look forward to it. You're one of the most helpful people in this community
jmac698 is offline   Reply With Quote
Old 5th October 2015, 03:51   #11  |  Link
StainlessS
HeartlessS Usurer
 
StainlessS's Avatar
 
Join Date: Dec 2009
Location: Over the rainbow
Posts: 10,980
Quote:
Originally Posted by jmac698 View Post
You're one of the most helpful people in ...
What you mean, one of !

A few words of note.
The ClipDiff thing is a little weird to get your head around.
Thresh applies to frames (whereas in eg YPlanemin it applies to Y pixels), so eg Thresh=1.0% for mode=MIN, means ignore up to 1% of frames (as the differences are counted separately for each and every pixel/channel in the frame). I've implemented almost exactly as in RT_Stats but on differences rather than Y levels.
Also, Mode=STDEV, on first frame will always produce gray $808080 frame as there is only a single difference to the corresponding pixel in reference frame (stdev of differences rather than eg Y levels).

Its real slow, perhaps a good idea to have Data file on Sata-3 Solid state if you have one.

Anyways, have a play and point out what mistakes I've made

EDIT: Oh yeah, MinMaxDifference is the MinMaxDifference of Differences.
__________________
I sometimes post sober.
StainlessS@MediaFire ::: AND/OR ::: StainlessS@SendSpace

"Some infinities are bigger than other infinities", but how many of them are infinitely bigger ???

Last edited by StainlessS; 5th October 2015 at 03:59.
StainlessS is offline   Reply With Quote
Old 5th October 2015, 08:13   #12  |  Link
ajk
Registered User
 
Join Date: Jan 2006
Location: Finland
Posts: 134
@jmac698

How far do you need the processing to reach? My median filter can extract the temporal minimum, maximum, average or median up to 25 frames at the moment, but that is just an arbitrary value and can be raised if needed.
ajk is offline   Reply With Quote
Old 5th October 2015, 17:26   #13  |  Link
StainlessS
HeartlessS Usurer
 
StainlessS's Avatar
 
Join Date: Dec 2009
Location: Over the rainbow
Posts: 10,980
Jmac, can you clarify, is the reference frame a separate frame not existing in the temporal clip ? (that is what I implemented in ClipDiffs).
__________________
I sometimes post sober.
StainlessS@MediaFire ::: AND/OR ::: StainlessS@SendSpace

"Some infinities are bigger than other infinities", but how many of them are infinitely bigger ???
StainlessS is offline   Reply With Quote
Old 5th October 2015, 21:51   #14  |  Link
StainlessS
HeartlessS Usurer
 
StainlessS's Avatar
 
Join Date: Dec 2009
Location: Over the rainbow
Posts: 10,980
Posts 8 and 9 updated.

Fixed a couple of problems, wish life were not so entertaining.

Made Mode = 3, MinMaxDifference relative 0, not 128. Same for mode=5, Stdev.

Check if I got em right JMac.

EDIT: This bit
Code:
    BYTE Diff[513];
    for (i=0; i<=512; i++) Diff[i] = max(0,min(255,i-129));
From Layer.cpp, Subtract Filter.
Use as in

Code:
    dstp[x] = Diff[z + 128 + 129];
To limit range, offset by 128, where z could be -255 -> 255.
__________________
I sometimes post sober.
StainlessS@MediaFire ::: AND/OR ::: StainlessS@SendSpace

"Some infinities are bigger than other infinities", but how many of them are infinitely bigger ???

Last edited by StainlessS; 5th October 2015 at 21:58.
StainlessS is offline   Reply With Quote
Reply

Thread Tools Search this Thread
Search this Thread:

Advanced Search
Display Modes

Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off

Forum Jump


All times are GMT +1. The time now is 17:32.


Powered by vBulletin® Version 3.8.11
Copyright ©2000 - 2024, vBulletin Solutions Inc.