source: trunk/Cosy/caos/Ring.cc@ 14751

#include "Ring.h"

#include <iostream>

#include <math.h>

#include "Led.h"
#include "Leds.h"

#include "Rings.h"
#include "MString.h"

ClassImp(Ring);

using namespace std;

Ring::Ring(Double_t x, Double_t y) :
    fX(x), fY(y), fR(0), fPhi(0), fDx(-1), fDy(-1), fDr(-1), fDphi(-1)
{
}

bool Ring::CalcCenter(const Leds &leds, Int_t i, Int_t j, Int_t k)
{
    if (leds.At(i)==NULL)
    {
        cout << "Ring::CalcCenter: Led i=" << i << " is NULL." << endl;
        return kFALSE;
    }
    if (leds.At(j)==NULL)
    {
        cout << "Ring::CalcCenter: Led j=" << j << " is NULL." << endl;
        return kFALSE;
    }
    if (leds.At(k)==NULL)
    {
        cout << "Ring::CalcCenter: Led k=" << k << " is NULL." << endl;
        return kFALSE;
    }

    Double_t h1 = leds(i).GetY()- leds(j).GetY();

    if (h1==0)
    {
        Swap(j, k);
        h1 = leds(i).GetY()- leds(j).GetY();
        if (h1==0)
        {
            cout << "Ring::CalcCenter: h1==0" <<endl;
            return kFALSE;
        }
    }

    Double_t h2 = leds(j).GetY() - leds(k).GetY();

    if (h2==0)
    {
        Swap(i, j);
        h2 = leds(j).GetY() - leds(k).GetY();
        if (h2==0)
        {
            cout << "Ring::CalcCenter: h2==0" << endl;
            return kFALSE;
        }
    }

    const Double_t w1 = leds(i).GetX() - leds(j).GetX();
    const Double_t w2 = leds(j).GetX() - leds(k).GetX();

    const Double_t m1 = -w1/h1;
    const Double_t m2 = -w2/h2;

    if (m2-m1==0)
    {
        cout << "Ring::CalcCenter: All three points in a row! (m2-m1==0)" << endl;
        return kFALSE;
    }

    fX = ((m2*(leds(j).GetX() + leds(k).GetX()) + leds(i).GetY() - leds(k).GetY()       -m1*(leds(i).GetX() + leds(j).GetX()))/(m2-m1)/2);
    fY = ((m2*(leds(i).GetY() + leds(j).GetY()) +m1*m2*(leds(k).GetX() - leds(i).GetX())-m1*(leds(j).GetY() + leds(k).GetY()))/(m2-m1)/2);

    fR = hypot(fX-leds(i).GetX(), fY-leds(i).GetY());

    fMag = (leds(i).GetMag() + leds(j).GetMag() + leds(k).GetMag())/3;

    return kTRUE;
}

void Ring::InterpolCenters(const Rings &rings)
{
    const int n=rings.GetEntries();

    fX = 0;
    fY = 0;
    fR = 0;

    fDx=0;
    fDy=0;
    fDr=0;

    fMag=0;

    if (n<1)
        return;

    for (int i=0; i<n; i++)
    {
        const Ring &ring = rings(i);

        fX   += ring.GetX();
        fY   += ring.GetY();
        fR   += ring.GetR();
        fMag += ring.GetMag();
    }

    fX   /= n;
    fY   /= n;
    fR   /= n;
    fMag /= n;

    if (n<2)
        return;

    //
    // deviation of x- and y coordinate and radius
    //
    for (int i=0; i<n; i++)
    {
        const Ring &ring = rings(i);

        fDx += sqr(ring.GetX()-fX);
        fDy += sqr(ring.GetY()-fY);
        fDr += sqr(ring.GetR()-fR);
    }

    fDx=sqrt(fDx)/n;
    fDy=sqrt(fDy)/n;
    fDr=sqrt(fDr)/n;
}

void Ring::Print(Option_t *o) const
{
    cout << "Ring: ";
    cout << "x="   << MString::Format("%5.1f", fX) << "+-" << MString::Format("%.1f", fDx) << ", ";
    cout << "y="   << MString::Format("%5.1f", fY) << "+-" << MString::Format("%.1f", fDy) << ", ";
    cout << "r="   << MString::Format("%5.1f", fR) << "+-" << MString::Format("%.1f", fDr) << ", ";
    cout << "phi=" << fPhi              << "+-" << fDphi << endl;
}