void ZOffset3D( TH3* h, float z0) {
//     This boosts the height for each bin of a histogram by a given constant value of z
   int nbinx = h->GetNbinsX();
   int nbiny = h->GetNbinsY();
   int nbinz = h->GetNbinsZ();
   for (int i=1; i<=nbinx; i++) {
      for (int j=1; j<=nbiny; j++) {
          for (int k=1; k<=nbinz; k++) {
            int h0=(h->GetBinContent(i,j,k));
            if (h0>0){h->SetBinContent(i,j,k,z0+h0);}
          }
      }
   }
    h->SetMinimum(0);
}


void NewPalette(int top, int z0,int MaxCol, int MinCol){
    Int_t palette[top]; int q0; int co;
    
    for (int qq=0;qq<top+1;qq++){
        q0=qq-z0;
        if (qq<z0){
            palette[qq]=MinCol;

        }
        if (qq>=z0){
    //         Set the rest of the values between MinCol and Max
            co= (MinCol+(MaxCol-MinCol)*q0*pow(top-z0,-1));            
            palette[qq]=co;
        }
    }
    gStyle->SetPalette(top,palette);

}



void Many3DHistograms(int Nc, int N){
// Nc is the number of "layers" I want, N is the ammount of bins on each axis per layer
    
    TCanvas * c16 = new TCanvas("c16","c16",600,400);
    gStyle->SetCanvasPreferGL(true);
    gStyle->SetOptStat(kFALSE);

//     Bins for the histograms
    Float_t Xh[N+1],Yh[N+1],Zh[2*Nc];
    int cc2=0;
    int Z[Nc]; //The distances between histograms
    for (int cc=0;cc<Nc;cc++){
        Z[cc]=cc*10;
        Zh[cc2]=cc*10; Zh[cc2+1]=cc*10+1;
        cc2=cc2+2;
    }
    
    for (int ij=0;ij<N+1;ij++){
        Xh[ij]=ij; Yh[ij]=ij; 
    }
    

//     To make the straight line
//     Coordinates of the first and last points
    int x1=N*(gRandom->Rndm()); int y1=N*(gRandom->Rndm()); 
    int xL=N*(gRandom->Rndm()); int yL=N*(gRandom->Rndm()); 
    int z1=Z[0]-5; int zL=Z[Nc-1]+5;

    int elem=1+Nc;
    
//     Ri (i=x,y,z) contains all lines
    float Rx[2*elem],Ry[2*elem],Rz[2*elem];
    Rx[0]=x1; Rx[2*elem-1]=xL; 
    Ry[0]=y1; Ry[2*elem-1]=yL; 
    Rz[0]=z1; Rz[2*elem-1]=zL; 
    

    float sloped, slopeu;
// sloped is the slope for the lines below each layer, and slopeu is for the one on top of each layer
//     It'll be the same for all the objects ending with "-d" (for down) or "-u" (for up).

// I'll fill all the Ri lists    
    int ee2=1;
    for (int ee=0;ee<elem-1;ee++){
//         
        Rz[ee2]=10*ee; sloped=(Rz[ee2]-z1)*pow(zL-z1,-1);
        Ry[ee2]=y1+sloped*(yL-y1);  
        Rx[ee2]=x1+sloped*(xL-x1);
                
        Rz[ee2+1]=10*ee+1; slopeu=(Rz[ee2+1]-z1)*pow(zL-z1,-1);
        Ry[ee2+1]=y1+slopeu*(yL-y1);  
        Rx[ee2+1]=x1+slopeu*(xL-x1);
        

        ee2=ee2+2;        
    }

    int ind=0;

    int Nin=2000;//Nin events

    
//     I'll fill the histograms and lines
    for (int c=0;c<Nc;c++){
        int j=0; int k=0;
        int zc=Z[c]; 
        int xc; int yc; 

        TH3I * hc = new TH3I("ch","Layer", N,Xh, N,Yh, 2*Nc-1,Zh);//Represents the particles detected on each layer

        
       for (int i=0; i<Nin; i++){
            xc=N*(gRandom->Rndm());
            yc=N*(gRandom->Rndm()); 
            hc->Fill(xc,yc,zc);
        }   
        
        
//     To make the straight lines, two per layer. 
        float Rxd[2]={Rx[ind],Rx[ind+1]};
        float Rxu[2]={Rx[ind+2],Rx[ind+3]};
        float Ryd[2]={Ry[ind],Ry[ind+1]};
        float Ryu[2]={Ry[ind+2],Ry[ind+3]};
        float Rzd[2]={Rz[ind],Rz[ind+1]};
        float Rzu[2]={Rz[ind+2],Rz[ind+3]};
//         cout<<"Down: "<<Rzd[0]<<"   "<<Rzd[1]<<endl;
//         cout<<"Up: "<<Rzu[0]<<"   "<<Rzu[1]<<endl;
// 
        ind=ind+2;
        
        TPolyLine3D *gd = new TPolyLine3D(2, Rxd, Ryd, Rzd);
        TPolyLine3D *gu = new TPolyLine3D(2, Rxu, Ryu, Rzu);
        gd->SetLineColor(kRed); gd->SetLineWidth(7);
        gu->SetLineColor(kRed); gu->SetLineWidth(7);

        
        
//      "Offset in volume" for the histogram
        ZOffset3D(hc, 100);         

//      Drawing the histograms and lines.
        if (c==0){
            int T=hc->GetMaximum();
            NewPalette(T, 100,19,12);  //greys, works great
            gd->Draw(); hc->Draw("box2z fb same"); gu->Draw("same"); }
        else{
            gd->Draw("same"); hc->Draw("box2 fb same"); gu->Draw("same"); }
        c16->Modified();
    }
   
}


int N=10; int Nc=3;
Many3DHistograms(Nc,N);