Спасибо всем большое. Будем думать

Господа! Предлагаю обратиться к высказанной уже хорошей идее - сплайн-аппроксимации.
Вот подпрограммка, которая аппроксимирует сплайном некоторую функцию в окрестности ее максимума. Т.е. мы сняли измерения с некоторымм грубым интервалом дискретизации, а "истинный" максимум находится где-то в стороне
от максимального измерения.

Код

void appr_max(double *sig, double *r, float tay, int maximum) {     float *g = new float[4];     float *g0 = new float[4];     float *g1 = new float[4];     float *g2 = new float[4];     float *g3 = new float[4]; int typ = 1; float t; float spl; // аппроксимация от sig[maximum-1] до sig[maximum]    for(int i=0; i<4; i++)        {       g0[i] = sig[(maximum-1)-1+i];       g[i] = g0[i];       if(typ == 0)           continue;       g1[i] = g0[i] - (sig[(maximum-1)-2+i] -                2*sig[(maximum-1)-1+i] + sig[(maximum-1)+i])/6;       g[i] = g1[i];       if(typ == 1)           continue;       float d = sig[(maximum-1)-3+i] - 4*sig[(maximum-1)-2+i] +                   6*sig[(maximum-1)-1+i] - 4*sig[(maximum-1)+i] +                                                    sig[(maximum-1)+1+i];       g2[i] = g1[i] + d/36;       g[i] = g2[i];       if(typ == 2)           continue;       g3[i] = g1[i] - 0.162*d;       g[i] = g3[i];    }     for(int nn = 1; nn < 1/tay; nn++) {        t = nn*tay;    spl = (g[3]*pow(t,3) + g[2]*(pow((1+t),3)-4*pow(t,3)) +    g[1]*(pow((2-t),3)-4*pow((1-t),3)) + g[0]*pow((1-t),3))/6;    r[nn] = spl; } // аппроксимация от sig[maximum] до sig[maximum+1]    for(int i=0; i<4; i++)        {       g0[i] = sig[maximum-1+i];       g[i] = g0[i];       if(typ == 0)           continue;       g1[i] = g0[i] - (sig[maximum-2+i] -                2*sig[maximum-1+i] + sig[maximum+i])/6;       g[i] = g1[i];       if(typ == 1)           continue;       float d = sig[maximum-3+i] - 4*sig[maximum-2+i] +                   6*sig[maximum-1+i] - 4*sig[maximum+i] +                                                    sig[maximum+1+i];       g2[i] = g1[i] + d/36;       g[i] = g2[i];       if(typ == 2)           continue;       g3[i] = g1[i] - 0.162*d;       g[i] = g3[i];    }     for(int mm = 1; mm < 1/tay; mm++) {        t = mm*tay;    spl = (g[3]*pow(t,3) + g[2]*(pow((1+t),3)-4*pow(t,3)) +    g[1]*(pow((2-t),3)-4*pow((1-t),3)) + g[0]*pow((1-t),3))/6;    r[(int)(1/tay+mm)] = spl; }     delete[] g; g = 0; delete[] g0; g0 = 0;   delete[] g1; g1 = 0;   delete[] g2; g2 = 0;   delete[] g3; g3 = 0; }

А вот как вызываем. Пусть в массиве max_array имеются отсчеты функции, максимум которой надо уточнить. В массив r будут записаны отсчёты сплайна - это и есть гладкая кривая, соединяющая наши точки-соседи максимума.

Код

// ищем номер максимального элемента max_array      int maximum = 0;      for(int i = 1;i < array_dim;i++)      {        if(max_array[i] > max_array[i-1])            maximum = i;      }     float tay = 0.01;  // это с какой дискретностью будем аппроксимировать                             // нашу функцию, хоть до пиксела (шутка)           // в этот массив пишутся значения  сплайна - это и есть     // отсчеты гладкой кривой     double *r = new double[2/tay+1];         // вызов функции     appr_max(max_array, r, tay, maximum);     // предположим, что хочется вывести график     for(int i = 0;i <= 2/tay;i++)        Сhart1->Series[0]->AddXY(i, r[i], "",clTeeColor);

Теперь несколько замечаний.
1. В задаче Семена можно пользоваться этим бредом, если предположить, что т.н. мой максимум - это его острый угол нарисованной кривой. Тут главное задать узлы интерполяции.
2. Как нет трудно заметить, в приведенной выше функции есть игра с некоторой переменной typ. Это тип сплайна. Тут есть сплайн минимального шаблона, сглаживающий сплайн. Ну вот хоть убей, не помню кто есть кто (давно дело было). По-момему, typ=1 настроен на СМШ. Он в принципе сглаживает прилично. Если Семен поиграется с этим "бредом", то может быть и поможет вспомнить.
3. Если увидели ошибки, сильно не ругайте. Не могу кнопку Edit нажать. Это админу относится - может постинги мои очень длинные?
4. Вот здесь немного для тех, кто не знает, о чем речь. А на языке "родных осин" можно провести такую аналогию. Это примерно, как рисование гладких кривых по нескольким точкам с помощью чертежного лекала.
5. Принят единичный интервал между измерениями, поэтому tau < 1.
6. Остальные замечания допишите сами, не могу больше, устал
7. Или спрашивайте

Вот вспомнил!
Кажись, есть в MATLAB к тулбоксу Spline Toolbox демка такая - тыкаешь мышкой несколько точек, а прога их соединяет кривой, вычисляемой на основе кубического сплайна. Довольно занятно!

В добавок к вышесказанному хочу рассказать, как я делал подобные "сплайны".
Пусть есть набор точек графика {xi,yi}.
Предполагаем погрешности точек равными, т.е. вес каждой точки(вклад в вид результирующего сплайна) одинаков.

Пусть мы приближаем график параболой, т.е. Yspl = ax^2+bx+c. a,b и с - на неизвестны. Выбираем поле перебора a,b и с, т.е. a будет меняется в диапазоне {amin..amax}, то же самое b и c. Теперь перебираем все значения a,b и с. Пусть для некоторых a,b и с имеем приближающий график Yspl, тогда оцениваем точность приближения, например по критерию X-квадрат(хи-квадрат):

X-квадрат(ну или почти он) = (Yspl[xi]-yi[xi])^2 (сумма по i).

Значение X-квадрат минимально в том наборе параметров a,b и c, где приближение максимально достоверно, те задача сводится к нахождению минимума Хи-квадрат на поле перебора a,b,c. Оценить достоверность найденного сплайна можно по таблицам Стюдента с заданным числом степеней свободы(тут их 3- a,b,c) и числом точек {yi(xi)}. Но это даже необязательно.
Естественно, что шаг перебора и тд можно менять в процессе приближения дабы ускорить процесс нахождения минимума Хи-квадрата(например, если приближаем многочленом 8-ой степени), можно сначала искать минимум по одному параметру и тд.

2 Chingachguk:
  Тебя не затруднит кинуть мне на мыло исходник твоей проги?

Конечно, нет

Только вот писал я такие вещи давно(криво тогда писал, непонятно). Например, есть даже что-то типа программы по приближению экспериментального распределения Гаусса(те нахождения сплайна для него) с графиками и интерфесом, но, повторюсь, очень все криво, к тому же на дикой смеси tp+asm. Скажи, и я все это вышлю !!!

С другой стороны, я могу по-быстрому состряпать программку на tp, находящую сплайн для графика параболы(то, что я выше написал). Изменить его для любого другого многочлена вам не составит труда.

Так что выбирай !

На АСМе я мало что пойму... Так, что нет другого выхода, кроме как попросить тебя сделать эту прогу... И автору темы, думаю, тоже будет интересно... Если тебя не затруднит, конечно...

{$G+}
{$N+}
const
 MaxPower = 10; {y(x)=a0+a1*x+a2*x^2+...aMaxPower*x^MaxPower}
type
 coff = array[0..MaxPower] of double;
const
 InFileName = 'in.dat'; OutFileName = 'out.dat';
 LatticeSize = 5;
 ImpLattices = 100;
 MaxPoints = 10;
 MaxErrorPercent = 25; {25%}
 StartX = 1; EndX = 10;
 Atest = 3; Btest = -5; Ctest = -8;
 ThisPower = 2; {current power - parabola}
 BreakValueOfXq = 1.0; {Look at the Student Table !!!}
var
 f:text;
 Data:array[1..MaxPoints,1..2] of double;
 i,NumberDataPoints,CurrLattice:longint;
 x,y,Xq,Xqmin,newamin,newamax:double;
 amin,amax,awork,afoundmin: coff;
 LatticeIndexes,NewIndexes:array[0..MaxPower] of longint;
{Get arg^power}
function get_pow(arg,power:double):double;
 begin
   if abs(arg) < 0.0001 then get_pow:=0
     else get_pow:=exp(power*ln(abs(arg)));
 end;
{Get value of LocCoff[i]*x^i, i=0..power}
function polinom(LocCoff:coff; power:longint; locx:double):double;
 var
   result:double;
   loci:longint;
 begin
   result:=LocCoff[0];
   for loci:=1 to power do
     begin
       result:=result+LocCoff[loci]*get_pow(locx,loci);
     end;
   polinom:=result;
 end;
Label NextApprox,NextLattice;
begin
{Create test data file}
 assign(f,InFileName);
 rewrite(f);
 randomize;
 for i:=0 to MaxPoints-1 do
   begin
     x:=StartX+((EndX-StartX+1)*i)/MaxPoints;
     y:=Atest*get_pow(x,2)+Btest*x+Ctest;
     y:=y+y*random(MaxErrorPercent)/100;
     writeln(f,x:4:0,' ',y:6:3);
   end;
 close(f);
{Read input data file: format: x y}
 NumberDataPoints:=1;
 assign(f,InFileName);
 reset(f);
 while not EOF(f) do
   begin
     readln(f,Data[NumberDataPoints,1],Data[NumberDataPoints,2]);
     inc(NumberDataPoints);
     if NumberDataPoints>MaxPoints then break;
   end;
 dec(NumberDataPoints);
 close(f);
{approximate Data by parabola}
{set initial coefficient's values}
 for i:=0 to ThisPower do
   begin
     amin[i]:=-MaxLongInt; {there ain't MaxDouble in tp7 ;)}
     amax[i]:=+MaxLongInt;
   end;
 CurrLattice:=1;
NextLattice:
{find Xq minimum at the current Lattice}
 Xqmin:=0;
 for i:=1 to NumberDataPoints
   do Xqmin:=Xqmin+get_pow(Data[i,2]-polinom(amin,ThisPower,Data[i,1]),2);
 for i:=0 to ThisPower do LatticeIndexes[i]:=0;
 repeat
   for i:=0 to ThisPower do
     awork[i]:=amin[i]+LatticeIndexes[i]*(amax[i]-amin[i])/(LatticeSize-1);
   Xq:=0;
   for i:=1 to NumberDataPoints do
     Xq:=Xq+get_pow(Data[i,2]-polinom(awork,ThisPower,Data[i,1]),2);
   if Xq<=Xqmin then
     begin
       Xqmin:=Xq;
       afoundmin:=awork;
       NewIndexes:=LatticeIndexes;
     end;
   i:=0;
NextApprox:
   inc(LatticeIndexes[i]);
   if LatticeIndexes[i]=LatticeSize then
     begin
       LatticeIndexes[i]:=0;
       inc(i);
       if i>ThisPower then break else goto NextApprox;
     end;
 until FALSE;
 for i:=0 to ThisPower do
   begin
     if (NewIndexes[i]=0) or (NewIndexes[i]=LatticeSize-1) then
       begin
         {point lie on the rage of lattice: no change size of lattice}
         newamin:=afoundmin[i]-(amax[i]-amin[i]);
         newamax:=afoundmin[i]+(amax[i]-amin[i]);
       end
     else
       begin
         {point lie inside of lattice: minimize size of lattice}
         newamin:=afoundmin[i]-((amax[i]-amin[i])/(LatticeSize-1));
         newamax:=afoundmin[i]+((amax[i]-amin[i])/(LatticeSize-1));
       end;
     amin[i]:=newamin;
     amax[i]:=newamax;
   end;
 inc(CurrLattice);
 if (Xqmin/NumberDataPoints)>(BreakValueOfXq*ThisPower) then
   begin
     if CurrLattice<=ImpLattices then goto NextLattice;
   end;
 write('y(x)=');
 for i:=0 to ThisPower do write('(',afoundmin[i]:4:3,')x^',i:1,'+');
 writeln;
 writeln('Xqmin=',Xqmin/NumberDataPoints);
{Create out data file: x ydata ysplain}
 assign(f,OutFileName);
 rewrite(f);
 randomize;
 for i:=1 to NumberDataPoints do
   begin
     y:=polinom(afoundmin,ThisPower,Data[i,1]);
     writeln(f,Data[i,1]:4:0,' ',Data[i,2]:6:3,' ',y:6:3);
   end;
 close(f);
end.

ЗЫ: Сам понимаешь, код приближает многочленом ЛЮБОЙ степени - какую задашь(параметр ThesePower)

Огромное спасибо! И не думал, что так длинно получится... (твой стиль программирования смахивает на АСМ )
А что это:

Код

y:=y+y*random(MaxErrorPercent)/100;

Погрешности?

Да, погрешности ! Ведь в первом постинге чел писал про неровный рисунок мышью ...

Вообще-то такая задача решается АНАЛИТИЧЕСКИ:

d
---     Xq  = 0
d(ai), i=0..Степени многочлена,

но как я уже писал в теме "Нужна ли программисту математика" ...
К тому же не для всех функций существует решение в аналитическом виде

Вот интересно, господинъ Семен сюда заглядывал? Так здорово поговорили...

Да, довольно интересная тема, предлагаю ее в FAQ по паскалю закинуть. Че скажешь, Vit?

Да, конечно закинем