// @(#)root/graf:$Id$
// Author: Anna Kreshuk 18/11/2005

/*************************************************************************
 * Copyright (C) 1995-2005, Rene Brun and Fons Rademakers.               *
 * All rights reserved.                                                  *
 *                                                                       *
 * For the licensing terms see $ROOTSYS/LICENSE.                         *
 * For the list of contributors see $ROOTSYS/README/CREDITS.             *
 *************************************************************************/

#include "TGraphQQ.h"
#include "TAxis.h"
#include "TF1.h"
#include "TMath.h"
#include "TVirtualPad.h"
#include "TLine.h"

ClassImp(TGraphQQ)

//______________________________________________________________________________
//
// This class allows to draw quantile-quantile plots
//
// Plots can be drawn for 2 datasets or for a dataset and a theoretical
// distribution function
//
// 2 datasets:
//   Quantile-quantile plots are used to determine whether 2 samples come from
//   the same distribution.
//   A qq-plot draws the quantiles of one dataset against the quantile of the
//   the other. The quantiles of the dataset with fewer entries are on Y axis,
//   with more entries - on X axis.
//   A straight line, going through 0.25 and 0.75 quantiles is also plotted
//   for reference. It represents a robust linear fit, not sensitive to the
//   extremes of the datasets.
//   If the datasets come from the same distribution, points of the plot should
//   fall approximately on the 45 degrees line. If they have the same
//   distribution function, but location or scale different parameters,
//   they should still fall on the straight line, but not the 45 degrees one.
//   The greater their departure from the straight line, the more evidence there
//   is, that the datasets come from different distributions.
//   The advantage of qq-plot is that it not only shows that the underlying
//   distributions are different, but, unlike the analytical methods, it also
//   gives information on the nature of this difference: heavier tails,
//   different location/scale, different shape, etc.
//
//   Some examples of qqplots of 2 datasets:
//Begin_Html
/*
<img src="gif/qqplots.gif">
*/
//End_Html
//
// 1 dataset:
//   Quantile-quantile plots are used to determine if the dataset comes from the
//   specified theoretical distribution, such as normal.
//   A qq-plot draws quantiles of the dataset against quantiles of the specified
//   theoretical distribution.
//   (NOTE, that density, not CDF should be specified)
//   A straight line, going through 0.25 and 0.75 quantiles can also be plotted
//   for reference. It represents a robust linear fit, not sensitive to the
//   extremes of the dataset.
//   As in the 2 datasets case, departures from straight line indicate departures
//   from the specified distribution.
//
//   " The correlation coefficient associated with the linear fit to the data
//     in the probability plot (qq plot in our case) is a measure of the
//     goodness of the fit.
//     Estimates of the location and scale parameters  of the distribution
//     are given by the intercept and slope. Probability plots can be generated
//     for several competing distributions to see which provides the best fit,
//     and the probability plot generating the highest correlation coefficient
//     is the best choice since it generates the straightest probability plot."
//   From "Engineering statistic handbook",
//   http://www.itl.nist.gov/div898/handbook/eda/section3/probplot.htm
//
//   Example of a qq-plot of a dataset from N(3, 2) distribution and
//           TMath::Gaus(0, 1) theoretical function. Fitting parameters
//           are estimates of the distribution mean and sigma.
//
//Begin_Html
/*
<img src="gif/qqnormal.gif">
*/
//End_Html//
//
//
// References:
// http://www.itl.nist.gov/div898/handbook/eda/section3/qqplot.htm
// http://www.itl.nist.gov/div898/handbook/eda/section3/probplot.htm
//



//______________________________________________________________________________
TGraphQQ::TGraphQQ()
{
   //default constructor

   fF   = 0;
   fY0  = 0;
   fNy0 = 0;
   fXq1 = 0.;
   fXq2 = 0.;
   fYq1 = 0.;
   fYq2 = 0.;

}


//______________________________________________________________________________
TGraphQQ::TGraphQQ(Int_t n, Double_t *x)
   : TGraph(n)
{
   //Creates a quantile-quantile plot of dataset x.
   //Theoretical distribution function can be defined later by SetFunction method

   fNy0 = 0;
   fXq1 = 0.;
   fXq2 = 0.;
   fYq1 = 0.;
   fYq2 = 0.;

   Int_t *index = new Int_t[n];
   TMath::Sort(n, x, index, kFALSE);
   for (Int_t i=0; i<fNpoints; i++)
      fY[i] = x[index[i]];
   fF=0;
   fY0=0;
   delete [] index;
}

//______________________________________________________________________________
TGraphQQ::TGraphQQ(Int_t n, Double_t *x, TF1 *f)
   : TGraph(n)
{
   //Creates a quantile-quantile plot of dataset x against function f

   fNy0 = 0;

   Int_t *index = new Int_t[n];
   TMath::Sort(n, x, index, kFALSE);
   for (Int_t i=0; i<fNpoints; i++)
      fY[i] = x[index[i]];
   delete [] index;
   fF = f;
   fY0=0;
   MakeFunctionQuantiles();
}


//______________________________________________________________________________
TGraphQQ::TGraphQQ(Int_t nx, Double_t *x, Int_t ny, Double_t *y)
{
   //Creates a quantile-quantile plot of dataset x against dataset y
   //Parameters nx and ny are respective array sizes

   fNy0 = 0;
   fXq1 = 0.;
   fXq2 = 0.;
   fYq1 = 0.;
   fYq2 = 0.;

   nx<=ny ? fNpoints=nx : fNpoints=ny;

   if (!CtorAllocate()) return;
   fF=0;
   Int_t *index = new Int_t[TMath::Max(nx, ny)];
   TMath::Sort(nx, x, index, kFALSE);
   if (nx <=ny){
      for (Int_t i=0; i<fNpoints; i++)
         fY[i] = x[index[i]];
      TMath::Sort(ny, y, index, kFALSE);
      if (nx==ny){
         for (Int_t i=0; i<fNpoints; i++)
            fX[i] = y[index[i]];
         fY0 = 0;
         Quartiles();
      } else {
         fNy0 = ny;
         fY0 = new Double_t[ny];
         for (Int_t i=0; i<ny; i++)
            fY0[i] = y[i];
         MakeQuantiles();
      }
   } else {
      fNy0 = nx;
      fY0 = new Double_t[nx];
      for (Int_t i=0; i<nx; i++)
         fY0[i] = x[index[i]];
      TMath::Sort(ny, y, index, kFALSE);
      for (Int_t i=0; i<ny; i++)
         fY[i] = y[index[i]];
      MakeQuantiles();
   }


   delete [] index;
}


//______________________________________________________________________________
TGraphQQ::~TGraphQQ()
{
   //Destroys a TGraphQQ

   if (fY0)
      delete [] fY0;
   if (fF)
      fF = 0;
}


//______________________________________________________________________________
void TGraphQQ::MakeFunctionQuantiles()
{
   //Computes quantiles of theoretical distribution function

   if (!fF) return;
   TString s = fF->GetTitle();
   Double_t pk;
   if (s.Contains("TMath::Gaus") || s.Contains("gaus")){
      //use plotting positions optimal for normal distribution
      for (Int_t k=1; k<=fNpoints; k++){
         pk = (k-0.375)/(fNpoints+0.25);
         fX[k-1]=TMath::NormQuantile(pk);
      }
   } else {
      Double_t *prob = new Double_t[fNpoints];
      if (fNpoints > 10){
         for (Int_t k=1; k<=fNpoints; k++)
            prob[k-1] = (k-0.5)/fNpoints;
      } else {
         for (Int_t k=1; k<=fNpoints; k++)
            prob[k-1] = (k-0.375)/(fNpoints+0.25);
      }
      //fF->GetQuantiles(fNpoints, prob, fX);
      fF->GetQuantiles(fNpoints, fX, prob);
      delete [] prob;
   }

   Quartiles();
}


//______________________________________________________________________________

void TGraphQQ::MakeQuantiles()
{
   //When sample sizes are not equal, computes quantiles of the bigger sample
   //by linear interpolation

   if (!fY0) return;

   Double_t pi, pfrac;
   Int_t pint;
   for (Int_t i=0; i<fNpoints-1; i++){
      pi = (fNy0-1)*Double_t(i)/Double_t(fNpoints-1);
      pint = TMath::FloorNint(pi);
      pfrac = pi - pint;
      fX[i] = (1-pfrac)*fY0[pint]+pfrac*fY0[pint+1];
   }
   fX[fNpoints-1]=fY0[fNy0-1];

   Quartiles();
}


//______________________________________________________________________________
void TGraphQQ::Quartiles()
{
   // compute quartiles
   // a quartile is a 25 per cent or 75 per cent quantile

   Double_t prob[]={0.25, 0.75};
   Double_t x[2];
   Double_t y[2];
   TMath::Quantiles(fNpoints, 2, fY, y, prob, kTRUE);
   if (fY0)
      TMath::Quantiles(fNy0, 2, fY0, x, prob, kTRUE);
   else if (fF) {
      TString s = fF->GetTitle();
      if (s.Contains("TMath::Gaus") || s.Contains("gaus")){
         x[0] = TMath::NormQuantile(0.25);
         x[1] = TMath::NormQuantile(0.75);
      } else
         fF->GetQuantiles(2, x, prob);
   }
   else
      TMath::Quantiles(fNpoints, 2, fX, x, prob, kTRUE);

   fXq1=x[0]; fXq2=x[1]; fYq1=y[0]; fYq2=y[1];
}


//______________________________________________________________________________
void TGraphQQ::SetFunction(TF1 *f)
{
   //Sets the theoretical distribution function (density!)
   //and computes its quantiles

   fF = f;
   MakeFunctionQuantiles();
}
 TGraphQQ.cxx:1
 TGraphQQ.cxx:2
 TGraphQQ.cxx:3
 TGraphQQ.cxx:4
 TGraphQQ.cxx:5
 TGraphQQ.cxx:6
 TGraphQQ.cxx:7
 TGraphQQ.cxx:8
 TGraphQQ.cxx:9
 TGraphQQ.cxx:10
 TGraphQQ.cxx:11
 TGraphQQ.cxx:12
 TGraphQQ.cxx:13
 TGraphQQ.cxx:14
 TGraphQQ.cxx:15
 TGraphQQ.cxx:16
 TGraphQQ.cxx:17
 TGraphQQ.cxx:18
 TGraphQQ.cxx:19
 TGraphQQ.cxx:20
 TGraphQQ.cxx:21
 TGraphQQ.cxx:22
 TGraphQQ.cxx:23
 TGraphQQ.cxx:24
 TGraphQQ.cxx:25
 TGraphQQ.cxx:26
 TGraphQQ.cxx:27
 TGraphQQ.cxx:28
 TGraphQQ.cxx:29
 TGraphQQ.cxx:30
 TGraphQQ.cxx:31
 TGraphQQ.cxx:32
 TGraphQQ.cxx:33
 TGraphQQ.cxx:34
 TGraphQQ.cxx:35
 TGraphQQ.cxx:36
 TGraphQQ.cxx:37
 TGraphQQ.cxx:38
 TGraphQQ.cxx:39
 TGraphQQ.cxx:40
 TGraphQQ.cxx:41
 TGraphQQ.cxx:42
 TGraphQQ.cxx:43
 TGraphQQ.cxx:44
 TGraphQQ.cxx:45
 TGraphQQ.cxx:46
 TGraphQQ.cxx:47
 TGraphQQ.cxx:48
 TGraphQQ.cxx:49
 TGraphQQ.cxx:50
 TGraphQQ.cxx:51
 TGraphQQ.cxx:52
 TGraphQQ.cxx:53
 TGraphQQ.cxx:54
 TGraphQQ.cxx:55
 TGraphQQ.cxx:56
 TGraphQQ.cxx:57
 TGraphQQ.cxx:58
 TGraphQQ.cxx:59
 TGraphQQ.cxx:60
 TGraphQQ.cxx:61
 TGraphQQ.cxx:62
 TGraphQQ.cxx:63
 TGraphQQ.cxx:64
 TGraphQQ.cxx:65
 TGraphQQ.cxx:66
 TGraphQQ.cxx:67
 TGraphQQ.cxx:68
 TGraphQQ.cxx:69
 TGraphQQ.cxx:70
 TGraphQQ.cxx:71
 TGraphQQ.cxx:72
 TGraphQQ.cxx:73
 TGraphQQ.cxx:74
 TGraphQQ.cxx:75
 TGraphQQ.cxx:76
 TGraphQQ.cxx:77
 TGraphQQ.cxx:78
 TGraphQQ.cxx:79
 TGraphQQ.cxx:80
 TGraphQQ.cxx:81
 TGraphQQ.cxx:82
 TGraphQQ.cxx:83
 TGraphQQ.cxx:84
 TGraphQQ.cxx:85
 TGraphQQ.cxx:86
 TGraphQQ.cxx:87
 TGraphQQ.cxx:88
 TGraphQQ.cxx:89
 TGraphQQ.cxx:90
 TGraphQQ.cxx:91
 TGraphQQ.cxx:92
 TGraphQQ.cxx:93
 TGraphQQ.cxx:94
 TGraphQQ.cxx:95
 TGraphQQ.cxx:96
 TGraphQQ.cxx:97
 TGraphQQ.cxx:98
 TGraphQQ.cxx:99
 TGraphQQ.cxx:100
 TGraphQQ.cxx:101
 TGraphQQ.cxx:102
 TGraphQQ.cxx:103
 TGraphQQ.cxx:104
 TGraphQQ.cxx:105
 TGraphQQ.cxx:106
 TGraphQQ.cxx:107
 TGraphQQ.cxx:108
 TGraphQQ.cxx:109
 TGraphQQ.cxx:110
 TGraphQQ.cxx:111
 TGraphQQ.cxx:112
 TGraphQQ.cxx:113
 TGraphQQ.cxx:114
 TGraphQQ.cxx:115
 TGraphQQ.cxx:116
 TGraphQQ.cxx:117
 TGraphQQ.cxx:118
 TGraphQQ.cxx:119
 TGraphQQ.cxx:120
 TGraphQQ.cxx:121
 TGraphQQ.cxx:122
 TGraphQQ.cxx:123
 TGraphQQ.cxx:124
 TGraphQQ.cxx:125
 TGraphQQ.cxx:126
 TGraphQQ.cxx:127
 TGraphQQ.cxx:128
 TGraphQQ.cxx:129
 TGraphQQ.cxx:130
 TGraphQQ.cxx:131
 TGraphQQ.cxx:132
 TGraphQQ.cxx:133
 TGraphQQ.cxx:134
 TGraphQQ.cxx:135
 TGraphQQ.cxx:136
 TGraphQQ.cxx:137
 TGraphQQ.cxx:138
 TGraphQQ.cxx:139
 TGraphQQ.cxx:140
 TGraphQQ.cxx:141
 TGraphQQ.cxx:142
 TGraphQQ.cxx:143
 TGraphQQ.cxx:144
 TGraphQQ.cxx:145
 TGraphQQ.cxx:146
 TGraphQQ.cxx:147
 TGraphQQ.cxx:148
 TGraphQQ.cxx:149
 TGraphQQ.cxx:150
 TGraphQQ.cxx:151
 TGraphQQ.cxx:152
 TGraphQQ.cxx:153
 TGraphQQ.cxx:154
 TGraphQQ.cxx:155
 TGraphQQ.cxx:156
 TGraphQQ.cxx:157
 TGraphQQ.cxx:158
 TGraphQQ.cxx:159
 TGraphQQ.cxx:160
 TGraphQQ.cxx:161
 TGraphQQ.cxx:162
 TGraphQQ.cxx:163
 TGraphQQ.cxx:164
 TGraphQQ.cxx:165
 TGraphQQ.cxx:166
 TGraphQQ.cxx:167
 TGraphQQ.cxx:168
 TGraphQQ.cxx:169
 TGraphQQ.cxx:170
 TGraphQQ.cxx:171
 TGraphQQ.cxx:172
 TGraphQQ.cxx:173
 TGraphQQ.cxx:174
 TGraphQQ.cxx:175
 TGraphQQ.cxx:176
 TGraphQQ.cxx:177
 TGraphQQ.cxx:178
 TGraphQQ.cxx:179
 TGraphQQ.cxx:180
 TGraphQQ.cxx:181
 TGraphQQ.cxx:182
 TGraphQQ.cxx:183
 TGraphQQ.cxx:184
 TGraphQQ.cxx:185
 TGraphQQ.cxx:186
 TGraphQQ.cxx:187
 TGraphQQ.cxx:188
 TGraphQQ.cxx:189
 TGraphQQ.cxx:190
 TGraphQQ.cxx:191
 TGraphQQ.cxx:192
 TGraphQQ.cxx:193
 TGraphQQ.cxx:194
 TGraphQQ.cxx:195
 TGraphQQ.cxx:196
 TGraphQQ.cxx:197
 TGraphQQ.cxx:198
 TGraphQQ.cxx:199
 TGraphQQ.cxx:200
 TGraphQQ.cxx:201
 TGraphQQ.cxx:202
 TGraphQQ.cxx:203
 TGraphQQ.cxx:204
 TGraphQQ.cxx:205
 TGraphQQ.cxx:206
 TGraphQQ.cxx:207
 TGraphQQ.cxx:208
 TGraphQQ.cxx:209
 TGraphQQ.cxx:210
 TGraphQQ.cxx:211
 TGraphQQ.cxx:212
 TGraphQQ.cxx:213
 TGraphQQ.cxx:214
 TGraphQQ.cxx:215
 TGraphQQ.cxx:216
 TGraphQQ.cxx:217
 TGraphQQ.cxx:218
 TGraphQQ.cxx:219
 TGraphQQ.cxx:220
 TGraphQQ.cxx:221
 TGraphQQ.cxx:222
 TGraphQQ.cxx:223
 TGraphQQ.cxx:224
 TGraphQQ.cxx:225
 TGraphQQ.cxx:226
 TGraphQQ.cxx:227
 TGraphQQ.cxx:228
 TGraphQQ.cxx:229
 TGraphQQ.cxx:230
 TGraphQQ.cxx:231
 TGraphQQ.cxx:232
 TGraphQQ.cxx:233
 TGraphQQ.cxx:234
 TGraphQQ.cxx:235
 TGraphQQ.cxx:236
 TGraphQQ.cxx:237
 TGraphQQ.cxx:238
 TGraphQQ.cxx:239
 TGraphQQ.cxx:240
 TGraphQQ.cxx:241
 TGraphQQ.cxx:242
 TGraphQQ.cxx:243
 TGraphQQ.cxx:244
 TGraphQQ.cxx:245
 TGraphQQ.cxx:246
 TGraphQQ.cxx:247
 TGraphQQ.cxx:248
 TGraphQQ.cxx:249
 TGraphQQ.cxx:250
 TGraphQQ.cxx:251
 TGraphQQ.cxx:252
 TGraphQQ.cxx:253
 TGraphQQ.cxx:254
 TGraphQQ.cxx:255
 TGraphQQ.cxx:256
 TGraphQQ.cxx:257
 TGraphQQ.cxx:258
 TGraphQQ.cxx:259
 TGraphQQ.cxx:260
 TGraphQQ.cxx:261
 TGraphQQ.cxx:262
 TGraphQQ.cxx:263
 TGraphQQ.cxx:264
 TGraphQQ.cxx:265
 TGraphQQ.cxx:266
 TGraphQQ.cxx:267
 TGraphQQ.cxx:268
 TGraphQQ.cxx:269
 TGraphQQ.cxx:270
 TGraphQQ.cxx:271
 TGraphQQ.cxx:272
 TGraphQQ.cxx:273
 TGraphQQ.cxx:274
 TGraphQQ.cxx:275
 TGraphQQ.cxx:276
 TGraphQQ.cxx:277
 TGraphQQ.cxx:278
 TGraphQQ.cxx:279
 TGraphQQ.cxx:280
 TGraphQQ.cxx:281
 TGraphQQ.cxx:282
 TGraphQQ.cxx:283
 TGraphQQ.cxx:284
 TGraphQQ.cxx:285
 TGraphQQ.cxx:286
 TGraphQQ.cxx:287
 TGraphQQ.cxx:288
 TGraphQQ.cxx:289
 TGraphQQ.cxx:290
 TGraphQQ.cxx:291
 TGraphQQ.cxx:292
 TGraphQQ.cxx:293
 TGraphQQ.cxx:294
 TGraphQQ.cxx:295
 TGraphQQ.cxx:296
 TGraphQQ.cxx:297
 TGraphQQ.cxx:298
 TGraphQQ.cxx:299
 TGraphQQ.cxx:300
 TGraphQQ.cxx:301
 TGraphQQ.cxx:302
 TGraphQQ.cxx:303
 TGraphQQ.cxx:304
 TGraphQQ.cxx:305