//Copyright (c) 2011 ashelly.myopenid.com under <http://w...content-available-to-author-only...e.org/licenses/mit-license>
 
#include <stdlib.h>
#define inline
 
typedef int Item;
typedef struct Mediator_t
{
   Item* data;  //circular queue of values
   int*  pos;   //index into `heap` for each value
   int*  heap;  //max/median/min heap holding indexes into `data`.
   int   N;     //allocated size.
   int   idx;   //position in circular queue
   int   minCt; //count of items in min heap
   int   maxCt; //count of items in max heap
} Mediator;
 
/*--- Helper Functions ---*/
 
//returns 1 if heap[i] < heap[j]
inline int mmless(Mediator* m, int i, int j)
{
   return (m->data[m->heap[i]] < m->data[m->heap[j]]);
}
 
//swaps items i&j in heap, maintains indexes
int mmexchange(Mediator* m, int i, int j)
{
   int t = m->heap[i];
   m->heap[i]=m->heap[j];
   m->heap[j]=t;
   m->pos[m->heap[i]]=i;
   m->pos[m->heap[j]]=j;
   return 1;
}
 
//swaps items i&j if i<j;  returns true if swapped
inline int mmCmpExch(Mediator* m, int i, int j)
{
   return (mmless(m,i,j) && mmexchange(m,i,j));
}
 
//maintains minheap property for all items below i.
void minSortDown(Mediator* m, int i)
{
   for (i*=2; i <= m->minCt; i*=2)
   {  if (i < m->minCt && mmless(m, i+1, i)) { ++i; }
      if (!mmCmpExch(m,i,i/2)) { break; }
   }
}
 
//maintains maxheap property for all items below i. (negative indexes)
void maxSortDown(Mediator* m, int i)
{
   for (i*=2; i >= -m->maxCt; i*=2)
   {  if (i > -m->maxCt && mmless(m, i, i-1)) { --i; }
      if (!mmCmpExch(m,i/2,i)) { break; }
   }
}
 
//maintains minheap property for all items above i, including median
//returns true if median changed
inline int minSortUp(Mediator* m, int i)
{
   while (i>0 && mmCmpExch(m,i,i/2)) i/=2;
   return (i==0);
}
 
//maintains maxheap property for all items above i, including median
//returns true if median changed
inline int maxSortUp(Mediator* m, int i)
{
   while (i<0 && mmCmpExch(m,i/2,i))  i/=2;
   return (i==0);
}
 
/*--- Public Interface ---*/
 
 
//creates new Mediator: to calculate `nItems` running median. 
//mallocs single block of memory, caller must free.
Mediator* MediatorNew(int nItems)
{
   int size = sizeof(Mediator)+nItems*(sizeof(Item)+sizeof(int)*2);
   Mediator* m=  malloc(size);
   m->data= (Item*)(m+1);
   m->pos = (int*) (m->data+nItems);
   m->heap = m->pos+nItems + (nItems/2); //points to middle of storage.
   m->N=nItems;
   m->minCt = m->maxCt = m->idx = 0;
   while (nItems--)  //set up initial heap fill pattern: median,max,min,max,...
   {  m->pos[nItems]= ((nItems+1)/2) * ((nItems&1)?-1:1);
      m->heap[m->pos[nItems]]=nItems;
   }
   return m;
}
 
 
//Inserts item, maintains median in O(lg nItems)
void MediatorInsert(Mediator* m, Item v)
{
   int p = m->pos[m->idx];
   Item old = m->data[m->idx];
   m->data[m->idx]=v;
   m->idx = (m->idx+1) % m->N;
   if (p>0)         //new item is in minHeap
   {  if (m->minCt < (m->N-1)/2)  { m->minCt++; }
      else if (v>old) { minSortDown(m,p); return; }
      if (minSortUp(m,p) && mmCmpExch(m,0,-1)) { maxSortDown(m,-1); }
   }
   else if (p<0)   //new item is in maxheap
   {  if (m->maxCt < m->N/2) { m->maxCt++; }
      else if (v<old) { maxSortDown(m,p); return; }
      if (maxSortUp(m,p) && m->minCt && mmCmpExch(m,1,0)) { minSortDown(m,1); }
   }
   else //new item is at median
   {  if (m->maxCt && maxSortUp(m,-1)) { maxSortDown(m,-1); }
      if (m->minCt && minSortUp(m, 1)) { minSortDown(m, 1); }
   }
}
 
//returns median item (or average of 2 when item count is even)
Item MediatorMedian(Mediator* m)
{
   Item v= m->data[m->heap[0]];
   if (m->minCt<m->maxCt) { v=(v+m->data[m->heap[-1]])/2; }
   return v;
}
 
 
 
 
/*--- Test Code ---*/
#include <stdio.h>
 
void PrintMaxHeap(Mediator* m)
{
   int i;
   if(m->maxCt) printf("Max: %3d",m->data[m->heap[-1]]);
   for (i=2;i<=m->maxCt;++i)
   {  printf("|%3d ",m->data[m->heap[-i]]);
      if(++i<=m->maxCt) printf("%3d",m->data[m->heap[-i]]);
   }
   printf("\n");
}
void PrintMinHeap(Mediator* m)
{
   int i;
   if(m->minCt) printf("Min: %3d",m->data[m->heap[1]]);
   for (i=2;i<=m->minCt;++i)
   {  printf("|%3d ",m->data[m->heap[i]]);
      if(++i<=m->minCt) printf("%3d",m->data[m->heap[i]]);
   }
   printf("\n");
}
void ShowTree(Mediator* m)
{
   PrintMaxHeap(m);
   printf("Mid: %3d\n",m->data[m->heap[0]]);
   PrintMinHeap(m);
   printf("\n");
}
 
int main(int argc, char* argv[])
{
   int i,v;
   Mediator* m = MediatorNew(15);
 
   for (i=0;i<30;i++)
   {
      v = rand()&127;
      printf("Inserting %3d \n",v);
      MediatorInsert(m,v);
      v=MediatorMedian(m);
      printf("Median = %3d.\n\n",v);
      ShowTree(m);
   }
}

Ly9Db3B5cmlnaHQgKGMpIDIwMTEgYXNoZWxseS5teW9wZW5pZC5jb20gdW5kZXIgPGh0dHA6Ly93Li4uY29udGVudC1hdmFpbGFibGUtdG8tYXV0aG9yLW9ubHkuLi5lLm9yZy9saWNlbnNlcy9taXQtbGljZW5zZT4KCiNpbmNsdWRlIDxzdGRsaWIuaD4KI2RlZmluZSBpbmxpbmUKCnR5cGVkZWYgaW50IEl0ZW07CnR5cGVkZWYgc3RydWN0IE1lZGlhdG9yX3QKewogICBJdGVtKiBkYXRhOyAgLy9jaXJjdWxhciBxdWV1ZSBvZiB2YWx1ZXMKICAgaW50KiAgcG9zOyAgIC8vaW5kZXggaW50byBgaGVhcGAgZm9yIGVhY2ggdmFsdWUKICAgaW50KiAgaGVhcDsgIC8vbWF4L21lZGlhbi9taW4gaGVhcCBob2xkaW5nIGluZGV4ZXMgaW50byBgZGF0YWAuCiAgIGludCAgIE47ICAgICAvL2FsbG9jYXRlZCBzaXplLgogICBpbnQgICBpZHg7ICAgLy9wb3NpdGlvbiBpbiBjaXJjdWxhciBxdWV1ZQogICBpbnQgICBtaW5DdDsgLy9jb3VudCBvZiBpdGVtcyBpbiBtaW4gaGVhcAogICBpbnQgICBtYXhDdDsgLy9jb3VudCBvZiBpdGVtcyBpbiBtYXggaGVhcAp9IE1lZGlhdG9yOwoKLyotLS0gSGVscGVyIEZ1bmN0aW9ucyAtLS0qLwoKLy9yZXR1cm5zIDEgaWYgaGVhcFtpXSA8IGhlYXBbal0KaW5saW5lIGludCBtbWxlc3MoTWVkaWF0b3IqIG0sIGludCBpLCBpbnQgaikKewogICByZXR1cm4gKG0tPmRhdGFbbS0+aGVhcFtpXV0gPCBtLT5kYXRhW20tPmhlYXBbal1dKTsKfQoKLy9zd2FwcyBpdGVtcyBpJmogaW4gaGVhcCwgbWFpbnRhaW5zIGluZGV4ZXMKaW50IG1tZXhjaGFuZ2UoTWVkaWF0b3IqIG0sIGludCBpLCBpbnQgaikKewogICBpbnQgdCA9IG0tPmhlYXBbaV07CiAgIG0tPmhlYXBbaV09bS0+aGVhcFtqXTsKICAgbS0+aGVhcFtqXT10OwogICBtLT5wb3NbbS0+aGVhcFtpXV09aTsKICAgbS0+cG9zW20tPmhlYXBbal1dPWo7CiAgIHJldHVybiAxOwp9CgovL3N3YXBzIGl0ZW1zIGkmaiBpZiBpPGo7ICByZXR1cm5zIHRydWUgaWYgc3dhcHBlZAppbmxpbmUgaW50IG1tQ21wRXhjaChNZWRpYXRvciogbSwgaW50IGksIGludCBqKQp7CiAgIHJldHVybiAobW1sZXNzKG0saSxqKSAmJiBtbWV4Y2hhbmdlKG0saSxqKSk7Cn0KCi8vbWFpbnRhaW5zIG1pbmhlYXAgcHJvcGVydHkgZm9yIGFsbCBpdGVtcyBiZWxvdyBpLgp2b2lkIG1pblNvcnREb3duKE1lZGlhdG9yKiBtLCBpbnQgaSkKewogICBmb3IgKGkqPTI7IGkgPD0gbS0+bWluQ3Q7IGkqPTIpCiAgIHsgIGlmIChpIDwgbS0+bWluQ3QgJiYgbW1sZXNzKG0sIGkrMSwgaSkpIHsgKytpOyB9CiAgICAgIGlmICghbW1DbXBFeGNoKG0saSxpLzIpKSB7IGJyZWFrOyB9CiAgIH0KfQoKLy9tYWludGFpbnMgbWF4aGVhcCBwcm9wZXJ0eSBmb3IgYWxsIGl0ZW1zIGJlbG93IGkuIChuZWdhdGl2ZSBpbmRleGVzKQp2b2lkIG1heFNvcnREb3duKE1lZGlhdG9yKiBtLCBpbnQgaSkKewogICBmb3IgKGkqPTI7IGkgPj0gLW0tPm1heEN0OyBpKj0yKQogICB7ICBpZiAoaSA+IC1tLT5tYXhDdCAmJiBtbWxlc3MobSwgaSwgaS0xKSkgeyAtLWk7IH0KICAgICAgaWYgKCFtbUNtcEV4Y2gobSxpLzIsaSkpIHsgYnJlYWs7IH0KICAgfQp9CgovL21haW50YWlucyBtaW5oZWFwIHByb3BlcnR5IGZvciBhbGwgaXRlbXMgYWJvdmUgaSwgaW5jbHVkaW5nIG1lZGlhbgovL3JldHVybnMgdHJ1ZSBpZiBtZWRpYW4gY2hhbmdlZAppbmxpbmUgaW50IG1pblNvcnRVcChNZWRpYXRvciogbSwgaW50IGkpCnsKICAgd2hpbGUgKGk+MCAmJiBtbUNtcEV4Y2gobSxpLGkvMikpIGkvPTI7CiAgIHJldHVybiAoaT09MCk7Cn0KCi8vbWFpbnRhaW5zIG1heGhlYXAgcHJvcGVydHkgZm9yIGFsbCBpdGVtcyBhYm92ZSBpLCBpbmNsdWRpbmcgbWVkaWFuCi8vcmV0dXJucyB0cnVlIGlmIG1lZGlhbiBjaGFuZ2VkCmlubGluZSBpbnQgbWF4U29ydFVwKE1lZGlhdG9yKiBtLCBpbnQgaSkKewogICB3aGlsZSAoaTwwICYmIG1tQ21wRXhjaChtLGkvMixpKSkgIGkvPTI7CiAgIHJldHVybiAoaT09MCk7Cn0KCi8qLS0tIFB1YmxpYyBJbnRlcmZhY2UgLS0tKi8KCgovL2NyZWF0ZXMgbmV3IE1lZGlhdG9yOiB0byBjYWxjdWxhdGUgYG5JdGVtc2AgcnVubmluZyBtZWRpYW4uIAovL21hbGxvY3Mgc2luZ2xlIGJsb2NrIG9mIG1lbW9yeSwgY2FsbGVyIG11c3QgZnJlZS4KTWVkaWF0b3IqIE1lZGlhdG9yTmV3KGludCBuSXRlbXMpCnsKICAgaW50IHNpemUgPSBzaXplb2YoTWVkaWF0b3IpK25JdGVtcyooc2l6ZW9mKEl0ZW0pK3NpemVvZihpbnQpKjIpOwogICBNZWRpYXRvciogbT0gIG1hbGxvYyhzaXplKTsKICAgbS0+ZGF0YT0gKEl0ZW0qKShtKzEpOwogICBtLT5wb3MgPSAoaW50KikgKG0tPmRhdGErbkl0ZW1zKTsKICAgbS0+aGVhcCA9IG0tPnBvcytuSXRlbXMgKyAobkl0ZW1zLzIpOyAvL3BvaW50cyB0byBtaWRkbGUgb2Ygc3RvcmFnZS4KICAgbS0+Tj1uSXRlbXM7CiAgIG0tPm1pbkN0ID0gbS0+bWF4Q3QgPSBtLT5pZHggPSAwOwogICB3aGlsZSAobkl0ZW1zLS0pICAvL3NldCB1cCBpbml0aWFsIGhlYXAgZmlsbCBwYXR0ZXJuOiBtZWRpYW4sbWF4LG1pbixtYXgsLi4uCiAgIHsgIG0tPnBvc1tuSXRlbXNdPSAoKG5JdGVtcysxKS8yKSAqICgobkl0ZW1zJjEpPy0xOjEpOwogICAgICBtLT5oZWFwW20tPnBvc1tuSXRlbXNdXT1uSXRlbXM7CiAgIH0KICAgcmV0dXJuIG07Cn0KCgovL0luc2VydHMgaXRlbSwgbWFpbnRhaW5zIG1lZGlhbiBpbiBPKGxnIG5JdGVtcykKdm9pZCBNZWRpYXRvckluc2VydChNZWRpYXRvciogbSwgSXRlbSB2KQp7CiAgIGludCBwID0gbS0+cG9zW20tPmlkeF07CiAgIEl0ZW0gb2xkID0gbS0+ZGF0YVttLT5pZHhdOwogICBtLT5kYXRhW20tPmlkeF09djsKICAgbS0+aWR4ID0gKG0tPmlkeCsxKSAlIG0tPk47CiAgIGlmIChwPjApICAgICAgICAgLy9uZXcgaXRlbSBpcyBpbiBtaW5IZWFwCiAgIHsgIGlmIChtLT5taW5DdCA8IChtLT5OLTEpLzIpICB7IG0tPm1pbkN0Kys7IH0KICAgICAgZWxzZSBpZiAodj5vbGQpIHsgbWluU29ydERvd24obSxwKTsgcmV0dXJuOyB9CiAgICAgIGlmIChtaW5Tb3J0VXAobSxwKSAmJiBtbUNtcEV4Y2gobSwwLC0xKSkgeyBtYXhTb3J0RG93bihtLC0xKTsgfQogICB9CiAgIGVsc2UgaWYgKHA8MCkgICAvL25ldyBpdGVtIGlzIGluIG1heGhlYXAKICAgeyAgaWYgKG0tPm1heEN0IDwgbS0+Ti8yKSB7IG0tPm1heEN0Kys7IH0KICAgICAgZWxzZSBpZiAodjxvbGQpIHsgbWF4U29ydERvd24obSxwKTsgcmV0dXJuOyB9CiAgICAgIGlmIChtYXhTb3J0VXAobSxwKSAmJiBtLT5taW5DdCAmJiBtbUNtcEV4Y2gobSwxLDApKSB7IG1pblNvcnREb3duKG0sMSk7IH0KICAgfQogICBlbHNlIC8vbmV3IGl0ZW0gaXMgYXQgbWVkaWFuCiAgIHsgIGlmIChtLT5tYXhDdCAmJiBtYXhTb3J0VXAobSwtMSkpIHsgbWF4U29ydERvd24obSwtMSk7IH0KICAgICAgaWYgKG0tPm1pbkN0ICYmIG1pblNvcnRVcChtLCAxKSkgeyBtaW5Tb3J0RG93bihtLCAxKTsgfQogICB9Cn0KCi8vcmV0dXJucyBtZWRpYW4gaXRlbSAob3IgYXZlcmFnZSBvZiAyIHdoZW4gaXRlbSBjb3VudCBpcyBldmVuKQpJdGVtIE1lZGlhdG9yTWVkaWFuKE1lZGlhdG9yKiBtKQp7CiAgIEl0ZW0gdj0gbS0+ZGF0YVttLT5oZWFwWzBdXTsKICAgaWYgKG0tPm1pbkN0PG0tPm1heEN0KSB7IHY9KHYrbS0+ZGF0YVttLT5oZWFwWy0xXV0pLzI7IH0KICAgcmV0dXJuIHY7Cn0KCgoKCi8qLS0tIFRlc3QgQ29kZSAtLS0qLwojaW5jbHVkZSA8c3RkaW8uaD4KCnZvaWQgUHJpbnRNYXhIZWFwKE1lZGlhdG9yKiBtKQp7CiAgIGludCBpOwogICBpZihtLT5tYXhDdCkgcHJpbnRmKCJNYXg6ICUzZCIsbS0+ZGF0YVttLT5oZWFwWy0xXV0pOwogICBmb3IgKGk9MjtpPD1tLT5tYXhDdDsrK2kpCiAgIHsgIHByaW50ZigifCUzZCAiLG0tPmRhdGFbbS0+aGVhcFstaV1dKTsKICAgICAgaWYoKytpPD1tLT5tYXhDdCkgcHJpbnRmKCIlM2QiLG0tPmRhdGFbbS0+aGVhcFstaV1dKTsKICAgfQogICBwcmludGYoIlxuIik7Cn0Kdm9pZCBQcmludE1pbkhlYXAoTWVkaWF0b3IqIG0pCnsKICAgaW50IGk7CiAgIGlmKG0tPm1pbkN0KSBwcmludGYoIk1pbjogJTNkIixtLT5kYXRhW20tPmhlYXBbMV1dKTsKICAgZm9yIChpPTI7aTw9bS0+bWluQ3Q7KytpKQogICB7ICBwcmludGYoInwlM2QgIixtLT5kYXRhW20tPmhlYXBbaV1dKTsKICAgICAgaWYoKytpPD1tLT5taW5DdCkgcHJpbnRmKCIlM2QiLG0tPmRhdGFbbS0+aGVhcFtpXV0pOwogICB9CiAgIHByaW50ZigiXG4iKTsKfQp2b2lkIFNob3dUcmVlKE1lZGlhdG9yKiBtKQp7CiAgIFByaW50TWF4SGVhcChtKTsKICAgcHJpbnRmKCJNaWQ6ICUzZFxuIixtLT5kYXRhW20tPmhlYXBbMF1dKTsKICAgUHJpbnRNaW5IZWFwKG0pOwogICBwcmludGYoIlxuIik7Cn0KCmludCBtYWluKGludCBhcmdjLCBjaGFyKiBhcmd2W10pCnsKICAgaW50IGksdjsKICAgTWVkaWF0b3IqIG0gPSBNZWRpYXRvck5ldygxNSk7CgogICBmb3IgKGk9MDtpPDMwO2krKykKICAgewogICAgICB2ID0gcmFuZCgpJjEyNzsKICAgICAgcHJpbnRmKCJJbnNlcnRpbmcgJTNkIFxuIix2KTsKICAgICAgTWVkaWF0b3JJbnNlcnQobSx2KTsKICAgICAgdj1NZWRpYXRvck1lZGlhbihtKTsKICAgICAgcHJpbnRmKCJNZWRpYW4gPSAlM2QuXG5cbiIsdik7CiAgICAgIFNob3dUcmVlKG0pOwogICB9Cn0=