Axel Naumann | 2 Jul 11:02 2012
Picon
Picon

Re: Dimension of THnSparse

Hi Satya,

Thanks for reporting! This is now fixed in the trunk and v5-34-patches. 
Note that this does not affect the integrity of the THnSparse, but only 
its interaction with TBrowser.

Cheers, Axel.

On 06/30/2012 12:18 PM, Satyajit Jena wrote:
> Dear Root Experts,
>
> I wonder about the maximum dimension of a THnSparse. When I create an
> object of THnSparse with more then 16 dimensions, it is having an problem
> while opening in TBrowser and it shows me only 16 axises.
>
> However nor it complains, neither gives any error during creation and
> filling.
>
> Could you please have a look. I am attaching a demo macro. Following are
> the
> errors with ndim = 20; but it works well when ndim < 17; generated by
> attached
> macro.
>
> With regards,
> satya
> -----------------------------------------------
> root [0] TBrowser v
> root [1] (class TFile*)0x103c09090
> Error in <TObjArray::operator[]>: index 16 out of bounds (size: 16,
(Continue reading)

N. Kulkarni | 2 Jul 17:11 2012

Re: Dimension of THnSparse

Please remove me from this mailing list. I don't know why am I getting these emails at all.
Regards,
au2442


On 6/30/12 6:18 AM, Satyajit Jena wrote:
Dear Root Experts,

I wonder about the maximum dimension of a THnSparse. When I create an
object of THnSparse with more then 16 dimensions, it is having an problem
while opening in TBrowser and it shows me only 16 axises.

However nor it complains, neither gives any error during creation and filling.

Could you please have a look. I am attaching a demo macro. Following are the 
errors with ndim = 20; but it works well when ndim < 17; generated by attached
macro.

With regards,
satya
-----------------------------------------------
root [0] TBrowser v
root [1] (class TFile*)0x103c09090
Error in <TObjArray::operator[]>: index 16 out of bounds (size: 16, this: 0x1042516b0)
Error in <TObjArray::operator[]>: index 17 out of bounds (size: 16, this: 0x1042516b0)
Error in <TObjArray::operator[]>: index 18 out of bounds (size: 16, this: 0x1042516b0)
Error in <TObjArray::operator[]>: index 19 out of bounds (size: 16, this: 0x1042516b0)


Paul Prichard | 2 Jul 17:13 2012
Picon

Re: Dimension of THnSparse

unsubscribe






Esch, Shorty | 2 Jul 19:02 2012

Tiff issue

Hi, has anybody a way to read a TIFF file grayscale into a 2d histogram in root?
Thanks, 
Shorty

Olivier Couet | 3 Jul 10:03 2012
Picon
Picon

Re: Tiff issue

You will find the examples about images here:

On Jul 2, 2012, at 7:02 PM, Esch, Shorty wrote:

Hi, has anybody a way to read a TIFF file grayscale into a 2d histogram in root?
Thanks,
Shorty






Valeri Onuchin | 3 Jul 13:29 2012
Picon
Picon

RE: Tiff issue


Ooops failed to send because of attachment size limit.

________________________________________
From: Valeri Onuchin
Sent: 03 July 2012 13:24
To: Esch, Shorty; roottalk (Mailing list for ROOT users.)
Subject: RE: Tiff issue

Hi Shorty,
the tiff format is not supported yet.
If you are working under windows
the attached codes demonstrate
- how to read a 16bit TIFF fie, convert it to "raw" format
- how to read the produced "raw" file in ROOT and
  fill 2D hist with scaled red component.

HTH.

Regards. Valeriy

________________________________________
From: Esch, Shorty [ernst <at> lanl.gov]
Sent: 02 July 2012 19:02
To: roottalk (Mailing list for ROOT users.)
Subject: Tiff issue

Hi, has anybody a way to read a TIFF file grayscale into a 2d histogram in root?
Thanks,
Shorty
Attachment (tif2raw.exe): application/octet-stream, 66 KiB
//$Id: doseHist.C 1092 2012-04-09 07:46:57Z onuchin $
// Author: Valeriy Onuchin   28.07.2011

/*************************************************************************
 *                                                                       *
 * Copyright (C) 2012,   Valeriy Onuchin                                 *
 * All rights reserved.                                                  *
 *                                                                       *
 *************************************************************************/

#include <stdio.h>

#include "TSystem.h"
#include "TH2.h"
#include "TString.h"
#include "TROOT.h"
#include "TCanvas.h"
#include "TStyle.h"
#include "TSpectrum2.h"
#include "TAxis.h"
#include "TFile.h"


//#include "LandscapeRed.h"
#include "PortraitRed.h"
#include "icru49.h"
double *dose_vs_LR;
//double *dose_vs_PR;

TH2D *dose2d = 0;
TH1D *dose1d = 0;
double *dose;
TCanvas *colz;
TCanvas *contz;
TCanvas *surf;
TCanvas *dose1dCanvas;
TCanvas *peakCanvas;
TH2 *dose2rebin = 0;
TH2 *dose2surf = 0;
TSpectrum2 *spectrum;
Float_t *xpeaks;
Float_t *ypeaks;

TFile *rootFile;

//______________________________________________________________________________
void doseHist(const char *file, const char *calibr = "PortraitRed.h")
{
   // read RAW data file and create hists from this

   int min, arrsz;
   Bool_t land = false;

   int ret = 0;//gROOT->LoadMacro(calibr);

   if (ret) {
      fprintf(stderr, "Failed to open calibr file: %s\n", calibr);
   }

   if (land) {
      min = int(dose_vs_LR[0]);
      arrsz = sizeof(dose_vs_LR)/sizeof(double) - 1;
   } else {
      min = int(dose_vs_PR[0]);
      arrsz = sizeof(dose_vs_PR)/sizeof(double) - 1;
   }

   TString fname = gSystem->UnixPathName(file);

   Bool_t tifile = kFALSE;
   TString tmp_file = "";
   TString rootfileName = gSystem->BaseName(fname);

   // three possible extensions replace with .root
   rootfileName.ReplaceAll(".tiff", ".root");
   rootfileName.ReplaceAll(".tif", ".root");
   rootfileName.ReplaceAll(".raw", ".root");

   if (fname.EndsWith(".tif", TString::kIgnoreCase) || 
       fname.EndsWith(".tiff", TString::kIgnoreCase)) {
      tifile = kTRUE;
      char *tmpname = gSystem->ConcatFileName(gSystem->TempDirectory(), gSystem->BaseName(fname));
      TString tmp_file = tmpname;
      tmp_file.ReplaceAll(".tif", ".raw");
      delete tmpname;
      
      TString tif2raw =gSystem->ExpandPathName("$(ROOTSYS)/macros/");
      tif2raw += "tif2raw.exe -f " + fname + " -o " + tmp_file;

      //printf("Convert tif file to raw: %s \n", tif2raw.Data());
      gSystem->Exec(tif2raw);
      fname = tmp_file;
   }

   //printf("%d\n", min);

   FILE *fp = fopen(fname.Data(), "rb");
 
   if (!fp) {
      fprintf(stderr, "Failed to open file: %s", fname.Data());
      return;
   }

   UInt_t w, h;

   fread(&w, sizeof(w), 1, fp);
   fread(&h, sizeof(h), 1, fp);

   //printf("%d %d\n", w, h);

   long sz = w*h*3;
   unsigned short *rgb = new unsigned short[sz];
   unsigned short r, g, b;

   fread(rgb, sizeof(short), sz, fp);
   fclose(fp);

   // delete temporary file
   if (tifile) {
      gSystem->Unlink(tmp_file.Data());
   }
   fname = file;

   TString title = "Dose distribution of ";
   title += file;

   TString name = "colz_";
   name += fname;

   // open ROOT file to save all hists 
   rootFile = new TFile(rootfileName, "RECREATE");

   float dpi = 25.4/300;
   dose2d = new TH2D(name, title, w, 0, w*dpi, h, 0, h*dpi);
   dose2d->SetStats(0); //delete stats box

   name = "1D_";
   name += fname;

   dose1d = new TH1D(name, title, 10000, 0, 50);

   //printf("Start filling hitogram %s\n", dose2d->GetName());

   Double_t value;
   long j = 0;
   long idx;

   for (int y = 0; y < h; y++) {
      for (int x = 0; x < w; x++) {
         r = rgb[j];
//
         idx = (r - min) > arrsz ? arrsz : r - min;

         if (idx < 1) {
            printf("x=%d y=%d r=%d\n", x, y, r);
         }
         idx = idx < 1 ? 1 : idx; 

         value = land ? dose_vs_LR[idx] : dose_vs_PR[idx];

         g = rgb[++j];
         b = rgb[++j];
         dose2d->SetBinContent(x, h - y - 1, value);
         dose1d->Fill(value);

         j++;
      }
   }

   delete [] rgb;

   gStyle->SetCanvasPreferGL();
   gStyle->SetOptStat(1001211);

   title = "2D dose ";
   title += fname;

   colz = new TCanvas(fname + " colz", title + " colz in Gy", w*2-28, h*2-18);
   colz->SetGridx();
   colz->SetGridy();
   colz->SetTickx();
   colz->SetTicky();
   //colz->SetCrosshair();
   colz->ToggleEventStatus();
   dose2d->SetTitle("Dose distribution of " + fname + " in Gy");
   dose2d->Draw("colz");
   dose2d->GetXaxis()->SetTitle("in mm");
   dose2d->GetYaxis()->SetTitle("in mm");
   colz->SetEditable(kFALSE);

   int rebin = (int)300./25.4;
   dose2rebin = dose2d->Rebin2D(10, 10, "rebinXY");
   dose2rebin->SetStats(0); //delete stats box
   dose2rebin->SetName("iso_" + fname);
   dose2rebin->SetTitle("Isodose field of " + fname + " in cGy");
   contz = new TCanvas(title + " contz", title + " contz", w*2-28, h*2-18);
   contz->SetGridx();
   contz->SetGridy();
   contz->SetTickx();
   contz->SetTicky();
   //contz->SetCrosshair();
   contz->ToggleEventStatus();
   dose2rebin->Draw("cont1z");
   dose2rebin->GetXaxis()->SetTitle("in mm");
   dose2rebin->GetYaxis()->SetTitle("in mm");
   //contz->SetEditable(kFALSE);

   dose2surf = (TH2*)dose2rebin->Clone("dose2surf");
   surf = new TCanvas(title + " surf", title + " surf in cGy");
   dose2surf->Draw("glsurf2");
   dose2surf->GetXaxis()->SetTitle("");
   dose2surf->GetYaxis()->SetTitle("");
   dose2surf->GetZaxis()->SetTitle("");
   surf->SetGridx();
   surf->SetGridy();
   surf->SetTickx();
   surf->SetTicky();
   surf->ToggleEventStatus();

   dose1dCanvas = new TCanvas(title + " 1D", "1D histogram of dose in Gy");
   dose1dCanvas->SetGridx();
   dose1dCanvas->SetGridy();
   dose1dCanvas->SetTickx();
   dose1dCanvas->SetTicky();
   dose1dCanvas->SetLogy();
   //dose1dCanvas->SetCrosshair();
   dose1dCanvas->ToggleEventStatus();
   dose1d->SetFillColor(2);
   //dose1d->SetName("Dose_distribution_of_" + fname);
   dose1d->SetTitle("Dose distribution of " + fname + " in Gy");
   dose1d->Draw();
   dose1d->GetXaxis()->SetTitle("in Gy");
   //dose1dCanvas->SetEditable(kFALSE);

   gStyle->SetCanvasPreferGL(0);

   spectrum = new TSpectrum2(10);
   TCanvas *tmp = new TCanvas();
   Int_t nfound = spectrum->Search(dose2rebin, 4, "cont1z noMarkov");
   delete tmp;

   xpeaks = spectrum->GetPositionX();
   ypeaks = spectrum->GetPositionY();

   if (nfound >= 1) {
      //peakCanvas =  new TCanvas(name + " slices", name + " silces", 1024, 800);
      //peakCanvas->SetCrosshair();
      //peakCanvas->ToggleEventStatus();

      //peakCanvas->Divide(2, nfound);
      //peakCanvas->cd();

      TH1D *xproj;
      TH1D *yproj;
      float idp = 300./25.4;

      for (int i = 0; i < nfound; i++) {
         //peakCanvas->cd(i+1);

         printf("x = %f y = %f\n", xpeaks[i], ypeaks[i]);

         TString prname = TString::Format("Y[%d] slice at X = ", i);

         peakCanvas = new TCanvas(prname);
         gPad->SetGridx();
         gPad->SetGridy();
         gPad->SetTickx();
         gPad->SetTicky();
         peakCanvas->ToggleEventStatus();

         int xbin = int(xpeaks[i]*idp);
         printf("i = %d xbin = %d\n" , i, xbin);

         yproj = dose2d->ProjectionY(TString::Format("X = %f mm", xpeaks[i]), xbin, xbin, "d");
         yproj->SetFillColor(9);
         yproj->SetTitle(prname + TString::Format(" %f mm", xpeaks[i]));

         //gPad->SetTitle(prname + TString::Format(" %f mm", xpeaks[i]));
         //yproj->Draw();
         yproj->GetXaxis()->SetTitle("in mm");
         yproj->GetYaxis()->SetTitle("in Gy");
         //peakCanvas->SetEditable(kFALSE);

         //peakCanvas->cd(i+2);

         prname = TString::Format("X[%d] slice at Y = ", i);

         peakCanvas = new TCanvas(prname);
         gPad->SetGridx();
         gPad->SetGridy();
         gPad->SetTickx();
         gPad->SetTicky();
         peakCanvas->ToggleEventStatus();

         int ybin = int(ypeaks[i]*idp);
         printf("i = %d ybin = %d\n" , i, ybin);

         xproj = dose2d->ProjectionX(TString::Format("Y = %f mm", ypeaks[i]), ybin, ybin, "d");
         xproj->SetFillColor(4);
         xproj->SetTitle(prname + TString::Format(" %f mm", ypeaks[i]));

         //gPad->SetTitle(prname + TString::Format(" %f mm", ypeaks[i]));
         //xproj->Draw();
         xproj->GetXaxis()->SetTitle("in mm");
         xproj->GetYaxis()->SetTitle("in Gy");
         //peakCanvas->SetEditable(kFALSE);
      }
   }

   peakCanvas->Update();
   rootFile->Write();
   
   printf("Stop\n");


   printf("\n\n==========================================================================\n");

   double maxE = dose2d->GetMaximum();
   double halfE = maxE*0.5;
   TAxis *ax = dose1d->GetXaxis();
   int b2 = ax->FindBin(maxE);
   int b1 = ax->FindBin(halfE);
   printf("Full integral = %E.    Half integral = %E\n", dose2d->Integral(), dose1d->Integral(b1, b2));
   int energy = 100; //100 MeV
   //float dE1 = 7.2;  // 7.2 MeV/cm at 100 MeV
   //float dE2 = 9.57; // 9.57 MeV/cm

   double iddp = (300./2.54)*(300./2.54); // (pixels per cantimeter)**2
   double intensity = dose2d->Integral()/iddp/dEdx_ICRU49[energy]/10.*6.24e9; // dEdx_ICRU49 array -
stopping power of protons in water in Mev/mm
   printf("Total Intensity = %E at proton Energy = %d (MeV)\n", intensity, energy);
}

//______________________________________________________________________________
void Intensity(int energy = 100)
{
   // calculates total proton beam intensity based on total dose field integral

   printf("\n\n==========================================================================\n");
   double maxE = dose2d->GetMaximum();
   double halfE = maxE*0.5;
   TAxis *ax = dose1d->GetXaxis();
   int b2 = ax->FindBin(maxE);
   int b1 = ax->FindBin(halfE);
   printf("Full integral = %E.    Half integral = %E\n", dose2d->Integral(),  dose1d->Integral(b1, b2));

   double iddp = (300./2.54)*(300./2.54); // (pixels per cantimeter)**2  
   double intensity = dose2d->Integral()/iddp/dEdx_ICRU49[energy]/10.*6.24e9; // dEdx_ICRU49 array -
stopping power of protons in water in Mev/mm
   printf("Total Intensity = %E at proton Energy = %d (MeV)\n", intensity, energy);
}
// $Id: tif2raw.cxx 945 2012-02-29 14:24:30Z onuchin $
// Author: Valeriy Onuchin   28.07.2011


/*************************************************************************
 *                                                                       *
 * Copyright (C) 2011,   Valeriy Onuchin                                 *
 * All rights reserved.                                                  *
 *                                                                       *
 *************************************************************************/

#include <windows.h>
#include <cctype>
#include <stdio.h>


#include "vcclr.h"          // for gcroot
using namespace System;
using namespace System::IO;
using namespace System::Drawing;
using namespace System::Windows::Media::Imaging;


#pragma warning(disable:4996) // CRT unsafe

#define LOCAL_DEBUG 1 // debugging


static char gUsageHelp[] = "\
Usage:   %s -f infile\n\
Index:\n\
     \"infile\"   is a TIFF image file (must have \".tif\" or \".tiff\" extension)\n\
                the output file will have the same name\n\
                as \"infile\" but with \".raw\" extension.\n\n\
Example:   tif2raw.exe -f scan111.tif\n\n";


static bool  gTickIsOn = false;

//______________________________________________________________________________
static void tickThreadStub(void *)
{
   //

   while (1) {
      ::Sleep(1000);

#ifdef LOCAL_DEBUG
      if (gTickIsOn) {
         fprintf(stderr, ".");
      }
#endif

   }
}


////////////////////////////////////////////////////////////////////////////////
//______________________________________________________________________________
static char *in2out(char *infile)
{
   // change file extension

   char *out = 0;

   out = strstr(infile, ".tif");

   if (out) {
      out[1] = 'r';
      out[2] = 'a';
      out[3] = 'w';
   } else {
      out = strstr(infile, ".tiff");
      
      if (out) {
         out[1] = 'r';
         out[2] = 'a';
         out[3] = 'w';
         out[4] = '\x0';
      } else {
         fprintf(stderr, "Input file must have extension \".dat\"");
         return 0;
      }
   }

   return infile;
}

//______________________________________________________________________________
static char *basename(char *pathname)
{
   static char *result;
   char *filename = NULL;
   DWORD len;

   free(result);
   len = 0;
   len = ::GetFullPathNameA(pathname, len, result, &filename);

   result = (char*)malloc(len);

   if (result) {
      if (::GetFullPathNameA(pathname, len, result, &filename) <= len) {
         return filename;
      }
   }
   return filename;
}



////////////////////////////////////////////////////////////////////////////////
//______________________________________________________________________________
int main(int argc, char *argv[])
{
   //

   char *optarg;
   char rawfile[512];
   char tiffile[512];

   memset(rawfile, 0, sizeof(rawfile));
   memset(tiffile, 0, sizeof(tiffile));

   bool hasinput = false;

   if ((argc > 1) && (!_stricmp(argv[1], "-h") || !_stricmp(argv[1], "/h") ||
                      !_stricmp(argv[1], "--help") ||
                      !_stricmp(argv[1], "-?") || !_stricmp(argv[1], "/?")) || (argc == 1)) {
      fprintf(stderr, gUsageHelp, basename(argv[0]));
      exit(0);
   }

   // check  flag
   if (!_stricmp(argv[1], "-f")) {
      optarg = argv[2];

      while (isspace(*optarg)) {
         optarg++;
      }

      strncpy(tiffile, optarg, sizeof(tiffile) - 1);
      tiffile[sizeof(tiffile) - 1] = '\x0';
      hasinput = true;
   }

   if (!hasinput) {
      fprintf(stderr, gUsageHelp, basename(argv[0]));
      exit(0);
   }

   HANDLE tickThread = ::CreateThread(NULL, NULL,
                                      (LPTHREAD_START_ROUTINE)tickThreadStub,
                                       0, NULL, NULL);

#ifdef LOCAL_DEBUG
   fprintf(stderr, "Reading TIFF file: %s\n", tiffile);
#endif

   gTickIsOn = true;

   FileStream^ stream;

   try {
      stream = gcnew FileStream(gcnew String(tiffile),
                                FileMode::Open,
                                FileAccess::Read,
                                FileShare::Read);

   } catch (Exception^ e)  {
      Console::WriteLine("Failed to read file. Reason: {0}", e->Message);
      throw;
   }

   Windows::Media::Imaging::TiffBitmapDecoder^ decoder;

   try {
      decoder = gcnew Windows::Media::Imaging::TiffBitmapDecoder(stream, 
                                                                 BitmapCreateOptions::PreservePixelFormat,
                                                                 BitmapCacheOption::Default);
   } catch (Exception^ e)  {
      Console::WriteLine("Failed to decode TIFF format. Reason: {0}", e->Message);
      throw;
   }

   Windows::Media::Imaging::BitmapFrame^ bitmap = decoder->Frames[0];

   int stride = bitmap->PixelWidth*6;
   int picSize = bitmap->PixelWidth*bitmap->PixelHeight*3; // w*h*(rgb  == 3 bytes)

   array<short>^ pixels = gcnew array<short>(picSize);

   if (bitmap->Format == Windows::Media::PixelFormats::Rgb48) {
      bitmap->CopyPixels(pixels, stride, 0);
   } else {
      Console::WriteLine("Wrong image format: {0}", bitmap->Format);
   }

   stream->Close();

   gTickIsOn = false;

   if (tiffile[0]) {
      char *out = in2out(tiffile);

      if (!out) {
#ifdef LOCAL_DEBUG
         fprintf(stderr, "Wrong output or input file name ... exiting\n");
#endif
         exit(-1);
      }
      memcpy(rawfile, out, sizeof(rawfile));
   }

#ifdef LOCAL_DEBUG
   fprintf(stderr, "Writing RAW data to file: %s, size = %d bytes\n\t\tWidth = %d. Height = %d\n", 
           rawfile, (picSize + 2)*2, bitmap->PixelWidth, bitmap->PixelHeight);
#endif

   FileStream^ fs  = gcnew FileStream(gcnew String(rawfile), FileMode::Create);
   BinaryWriter^ writer = gcnew BinaryWriter(fs);

   writer->Write(bitmap->PixelWidth);
   writer->Write(bitmap->PixelHeight);

   for (int i = 0; i < picSize; i++) {
      writer->Write(pixels[i]);
   }

   fs->Close();

#ifdef LOCAL_DEBUG
   fprintf(stderr, "Done\n\n");
#endif

   exit(0);
}

Timur Pocheptsov | 3 Jul 13:50 2012
Picon
Picon

RE: Tiff issue

I'd say it MUST fail since you are attaching 'exe' file.
________________________________________
From: Valeri Onuchin
Sent: 03 July 2012 13:29
To: ernst <at> lanl.gov; roottalk (Mailing list for ROOT users.)
Subject: RE: Tiff issue

Ooops failed to send because of attachment size limit.

________________________________________
From: Valeri Onuchin
Sent: 03 July 2012 13:24
To: Esch, Shorty; roottalk (Mailing list for ROOT users.)
Subject: RE: Tiff issue

Hi Shorty,
the tiff format is not supported yet.
If you are working under windows
the attached codes demonstrate
- how to read a 16bit TIFF fie, convert it to "raw" format
- how to read the produced "raw" file in ROOT and
  fill 2D hist with scaled red component.

HTH.

Regards. Valeriy

________________________________________
From: Esch, Shorty [ernst <at> lanl.gov]
Sent: 02 July 2012 19:02
To: roottalk (Mailing list for ROOT users.)
Subject: Tiff issue

Hi, has anybody a way to read a TIFF file grayscale into a 2d histogram in root?
Thanks,
Shorty

Valeri Onuchin | 3 Jul 13:53 2012
Picon
Picon

RE: Tiff issue


one more try.

________________________________________
From: Valeri Onuchin
Sent: 03 July 2012 13:29
To: Esch, Shorty; roottalk (Mailing list for ROOT users.)
Subject: RE: Tiff issue

Ooops failed to send because of attachment size limit.

________________________________________
From: Valeri Onuchin
Sent: 03 July 2012 13:24
To: Esch, Shorty; roottalk (Mailing list for ROOT users.)
Subject: RE: Tiff issue

Hi Shorty,
the tiff format is not supported yet.
If you are working under windows
the attached codes demonstrate
- how to read a 16bit TIFF fie, convert it to "raw" format
- how to read the produced "raw" file in ROOT and
  fill 2D hist with scaled red component.

HTH.

Regards. Valeriy

________________________________________
From: Esch, Shorty [ernst <at> lanl.gov]
Sent: 02 July 2012 19:02
To: roottalk (Mailing list for ROOT users.)
Subject: Tiff issue

Hi, has anybody a way to read a TIFF file grayscale into a 2d histogram in root?
Thanks,
Shorty
// $Id: tif2raw.cxx 945 2012-02-29 14:24:30Z onuchin $
// Author: Valeriy Onuchin   28.07.2011


/*************************************************************************
 *                                                                       *
 * Copyright (C) 2011,   Valeriy Onuchin                                 *
 * All rights reserved.                                                  *
 *                                                                       *
 *************************************************************************/

#include <windows.h>
#include <cctype>
#include <stdio.h>


#include "vcclr.h"          // for gcroot
using namespace System;
using namespace System::IO;
using namespace System::Drawing;
using namespace System::Windows::Media::Imaging;


#pragma warning(disable:4996) // CRT unsafe

#define LOCAL_DEBUG 1 // debugging


static char gUsageHelp[] = "\
Usage:   %s -f infile\n\
Index:\n\
     \"infile\"   is a TIFF image file (must have \".tif\" or \".tiff\" extension)\n\
                the output file will have the same name\n\
                as \"infile\" but with \".raw\" extension.\n\n\
Example:   tif2raw.exe -f scan111.tif\n\n";


static bool  gTickIsOn = false;

//______________________________________________________________________________
static void tickThreadStub(void *)
{
   //

   while (1) {
      ::Sleep(1000);

#ifdef LOCAL_DEBUG
      if (gTickIsOn) {
         fprintf(stderr, ".");
      }
#endif

   }
}


////////////////////////////////////////////////////////////////////////////////
//______________________________________________________________________________
static char *in2out(char *infile)
{
   // change file extension

   char *out = 0;

   out = strstr(infile, ".tif");

   if (out) {
      out[1] = 'r';
      out[2] = 'a';
      out[3] = 'w';
   } else {
      out = strstr(infile, ".tiff");
      
      if (out) {
         out[1] = 'r';
         out[2] = 'a';
         out[3] = 'w';
         out[4] = '\x0';
      } else {
         fprintf(stderr, "Input file must have extension \".dat\"");
         return 0;
      }
   }

   return infile;
}

//______________________________________________________________________________
static char *basename(char *pathname)
{
   static char *result;
   char *filename = NULL;
   DWORD len;

   free(result);
   len = 0;
   len = ::GetFullPathNameA(pathname, len, result, &filename);

   result = (char*)malloc(len);

   if (result) {
      if (::GetFullPathNameA(pathname, len, result, &filename) <= len) {
         return filename;
      }
   }
   return filename;
}



////////////////////////////////////////////////////////////////////////////////
//______________________________________________________________________________
int main(int argc, char *argv[])
{
   //

   char *optarg;
   char rawfile[512];
   char tiffile[512];

   memset(rawfile, 0, sizeof(rawfile));
   memset(tiffile, 0, sizeof(tiffile));

   bool hasinput = false;

   if ((argc > 1) && (!_stricmp(argv[1], "-h") || !_stricmp(argv[1], "/h") ||
                      !_stricmp(argv[1], "--help") ||
                      !_stricmp(argv[1], "-?") || !_stricmp(argv[1], "/?")) || (argc == 1)) {
      fprintf(stderr, gUsageHelp, basename(argv[0]));
      exit(0);
   }

   // check  flag
   if (!_stricmp(argv[1], "-f")) {
      optarg = argv[2];

      while (isspace(*optarg)) {
         optarg++;
      }

      strncpy(tiffile, optarg, sizeof(tiffile) - 1);
      tiffile[sizeof(tiffile) - 1] = '\x0';
      hasinput = true;
   }

   if (!hasinput) {
      fprintf(stderr, gUsageHelp, basename(argv[0]));
      exit(0);
   }

   HANDLE tickThread = ::CreateThread(NULL, NULL,
                                      (LPTHREAD_START_ROUTINE)tickThreadStub,
                                       0, NULL, NULL);

#ifdef LOCAL_DEBUG
   fprintf(stderr, "Reading TIFF file: %s\n", tiffile);
#endif

   gTickIsOn = true;

   FileStream^ stream;

   try {
      stream = gcnew FileStream(gcnew String(tiffile),
                                FileMode::Open,
                                FileAccess::Read,
                                FileShare::Read);

   } catch (Exception^ e)  {
      Console::WriteLine("Failed to read file. Reason: {0}", e->Message);
      throw;
   }

   Windows::Media::Imaging::TiffBitmapDecoder^ decoder;

   try {
      decoder = gcnew Windows::Media::Imaging::TiffBitmapDecoder(stream, 
                                                                 BitmapCreateOptions::PreservePixelFormat,
                                                                 BitmapCacheOption::Default);
   } catch (Exception^ e)  {
      Console::WriteLine("Failed to decode TIFF format. Reason: {0}", e->Message);
      throw;
   }

   Windows::Media::Imaging::BitmapFrame^ bitmap = decoder->Frames[0];

   int stride = bitmap->PixelWidth*6;
   int picSize = bitmap->PixelWidth*bitmap->PixelHeight*3; // w*h*(rgb  == 3 bytes)

   array<short>^ pixels = gcnew array<short>(picSize);

   if (bitmap->Format == Windows::Media::PixelFormats::Rgb48) {
      bitmap->CopyPixels(pixels, stride, 0);
   } else {
      Console::WriteLine("Wrong image format: {0}", bitmap->Format);
   }

   stream->Close();

   gTickIsOn = false;

   if (tiffile[0]) {
      char *out = in2out(tiffile);

      if (!out) {
#ifdef LOCAL_DEBUG
         fprintf(stderr, "Wrong output or input file name ... exiting\n");
#endif
         exit(-1);
      }
      memcpy(rawfile, out, sizeof(rawfile));
   }

#ifdef LOCAL_DEBUG
   fprintf(stderr, "Writing RAW data to file: %s, size = %d bytes\n\t\tWidth = %d. Height = %d\n", 
           rawfile, (picSize + 2)*2, bitmap->PixelWidth, bitmap->PixelHeight);
#endif

   FileStream^ fs  = gcnew FileStream(gcnew String(rawfile), FileMode::Create);
   BinaryWriter^ writer = gcnew BinaryWriter(fs);

   writer->Write(bitmap->PixelWidth);
   writer->Write(bitmap->PixelHeight);

   for (int i = 0; i < picSize; i++) {
      writer->Write(pixels[i]);
   }

   fs->Close();

#ifdef LOCAL_DEBUG
   fprintf(stderr, "Done\n\n");
#endif

   exit(0);
}

//$Id: doseHist.C 1092 2012-04-09 07:46:57Z onuchin $
// Author: Valeriy Onuchin   28.07.2011

/*************************************************************************
 *                                                                       *
 * Copyright (C) 2012,   Valeriy Onuchin                                 *
 * All rights reserved.                                                  *
 *                                                                       *
 *************************************************************************/

#include <stdio.h>

#include "TSystem.h"
#include "TH2.h"
#include "TString.h"
#include "TROOT.h"
#include "TCanvas.h"
#include "TStyle.h"
#include "TSpectrum2.h"
#include "TAxis.h"
#include "TFile.h"


//#include "LandscapeRed.h"
#include "PortraitRed.h"
#include "icru49.h"
double *dose_vs_LR;
//double *dose_vs_PR;

TH2D *dose2d = 0;
TH1D *dose1d = 0;
double *dose;
TCanvas *colz;
TCanvas *contz;
TCanvas *surf;
TCanvas *dose1dCanvas;
TCanvas *peakCanvas;
TH2 *dose2rebin = 0;
TH2 *dose2surf = 0;
TSpectrum2 *spectrum;
Float_t *xpeaks;
Float_t *ypeaks;

TFile *rootFile;

//______________________________________________________________________________
void doseHist(const char *file, const char *calibr = "PortraitRed.h")
{
   // read RAW data file and create hists from this

   int min, arrsz;
   Bool_t land = false;

   int ret = 0;//gROOT->LoadMacro(calibr);

   if (ret) {
      fprintf(stderr, "Failed to open calibr file: %s\n", calibr);
   }

   if (land) {
      min = int(dose_vs_LR[0]);
      arrsz = sizeof(dose_vs_LR)/sizeof(double) - 1;
   } else {
      min = int(dose_vs_PR[0]);
      arrsz = sizeof(dose_vs_PR)/sizeof(double) - 1;
   }

   TString fname = gSystem->UnixPathName(file);

   Bool_t tifile = kFALSE;
   TString tmp_file = "";
   TString rootfileName = gSystem->BaseName(fname);

   // three possible extensions replace with .root
   rootfileName.ReplaceAll(".tiff", ".root");
   rootfileName.ReplaceAll(".tif", ".root");
   rootfileName.ReplaceAll(".raw", ".root");

   if (fname.EndsWith(".tif", TString::kIgnoreCase) || 
       fname.EndsWith(".tiff", TString::kIgnoreCase)) {
      tifile = kTRUE;
      char *tmpname = gSystem->ConcatFileName(gSystem->TempDirectory(), gSystem->BaseName(fname));
      TString tmp_file = tmpname;
      tmp_file.ReplaceAll(".tif", ".raw");
      delete tmpname;
      
      TString tif2raw =gSystem->ExpandPathName("$(ROOTSYS)/macros/");
      tif2raw += "tif2raw.exe -f " + fname + " -o " + tmp_file;

      //printf("Convert tif file to raw: %s \n", tif2raw.Data());
      gSystem->Exec(tif2raw);
      fname = tmp_file;
   }

   //printf("%d\n", min);

   FILE *fp = fopen(fname.Data(), "rb");
 
   if (!fp) {
      fprintf(stderr, "Failed to open file: %s", fname.Data());
      return;
   }

   UInt_t w, h;

   fread(&w, sizeof(w), 1, fp);
   fread(&h, sizeof(h), 1, fp);

   //printf("%d %d\n", w, h);

   long sz = w*h*3;
   unsigned short *rgb = new unsigned short[sz];
   unsigned short r, g, b;

   fread(rgb, sizeof(short), sz, fp);
   fclose(fp);

   // delete temporary file
   if (tifile) {
      gSystem->Unlink(tmp_file.Data());
   }
   fname = file;

   TString title = "Dose distribution of ";
   title += file;

   TString name = "colz_";
   name += fname;

   // open ROOT file to save all hists 
   rootFile = new TFile(rootfileName, "RECREATE");

   float dpi = 25.4/300;
   dose2d = new TH2D(name, title, w, 0, w*dpi, h, 0, h*dpi);
   dose2d->SetStats(0); //delete stats box

   name = "1D_";
   name += fname;

   dose1d = new TH1D(name, title, 10000, 0, 50);

   //printf("Start filling hitogram %s\n", dose2d->GetName());

   Double_t value;
   long j = 0;
   long idx;

   for (int y = 0; y < h; y++) {
      for (int x = 0; x < w; x++) {
         r = rgb[j];
//
         idx = (r - min) > arrsz ? arrsz : r - min;

         if (idx < 1) {
            printf("x=%d y=%d r=%d\n", x, y, r);
         }
         idx = idx < 1 ? 1 : idx; 

         value = land ? dose_vs_LR[idx] : dose_vs_PR[idx];

         g = rgb[++j];
         b = rgb[++j];
         dose2d->SetBinContent(x, h - y - 1, value);
         dose1d->Fill(value);

         j++;
      }
   }

   delete [] rgb;

   gStyle->SetCanvasPreferGL();
   gStyle->SetOptStat(1001211);

   title = "2D dose ";
   title += fname;

   colz = new TCanvas(fname + " colz", title + " colz in Gy", w*2-28, h*2-18);
   colz->SetGridx();
   colz->SetGridy();
   colz->SetTickx();
   colz->SetTicky();
   //colz->SetCrosshair();
   colz->ToggleEventStatus();
   dose2d->SetTitle("Dose distribution of " + fname + " in Gy");
   dose2d->Draw("colz");
   dose2d->GetXaxis()->SetTitle("in mm");
   dose2d->GetYaxis()->SetTitle("in mm");
   colz->SetEditable(kFALSE);

   int rebin = (int)300./25.4;
   dose2rebin = dose2d->Rebin2D(10, 10, "rebinXY");
   dose2rebin->SetStats(0); //delete stats box
   dose2rebin->SetName("iso_" + fname);
   dose2rebin->SetTitle("Isodose field of " + fname + " in cGy");
   contz = new TCanvas(title + " contz", title + " contz", w*2-28, h*2-18);
   contz->SetGridx();
   contz->SetGridy();
   contz->SetTickx();
   contz->SetTicky();
   //contz->SetCrosshair();
   contz->ToggleEventStatus();
   dose2rebin->Draw("cont1z");
   dose2rebin->GetXaxis()->SetTitle("in mm");
   dose2rebin->GetYaxis()->SetTitle("in mm");
   //contz->SetEditable(kFALSE);

   dose2surf = (TH2*)dose2rebin->Clone("dose2surf");
   surf = new TCanvas(title + " surf", title + " surf in cGy");
   dose2surf->Draw("glsurf2");
   dose2surf->GetXaxis()->SetTitle("");
   dose2surf->GetYaxis()->SetTitle("");
   dose2surf->GetZaxis()->SetTitle("");
   surf->SetGridx();
   surf->SetGridy();
   surf->SetTickx();
   surf->SetTicky();
   surf->ToggleEventStatus();

   dose1dCanvas = new TCanvas(title + " 1D", "1D histogram of dose in Gy");
   dose1dCanvas->SetGridx();
   dose1dCanvas->SetGridy();
   dose1dCanvas->SetTickx();
   dose1dCanvas->SetTicky();
   dose1dCanvas->SetLogy();
   //dose1dCanvas->SetCrosshair();
   dose1dCanvas->ToggleEventStatus();
   dose1d->SetFillColor(2);
   //dose1d->SetName("Dose_distribution_of_" + fname);
   dose1d->SetTitle("Dose distribution of " + fname + " in Gy");
   dose1d->Draw();
   dose1d->GetXaxis()->SetTitle("in Gy");
   //dose1dCanvas->SetEditable(kFALSE);

   gStyle->SetCanvasPreferGL(0);

   spectrum = new TSpectrum2(10);
   TCanvas *tmp = new TCanvas();
   Int_t nfound = spectrum->Search(dose2rebin, 4, "cont1z noMarkov");
   delete tmp;

   xpeaks = spectrum->GetPositionX();
   ypeaks = spectrum->GetPositionY();

   if (nfound >= 1) {
      //peakCanvas =  new TCanvas(name + " slices", name + " silces", 1024, 800);
      //peakCanvas->SetCrosshair();
      //peakCanvas->ToggleEventStatus();

      //peakCanvas->Divide(2, nfound);
      //peakCanvas->cd();

      TH1D *xproj;
      TH1D *yproj;
      float idp = 300./25.4;

      for (int i = 0; i < nfound; i++) {
         //peakCanvas->cd(i+1);

         printf("x = %f y = %f\n", xpeaks[i], ypeaks[i]);

         TString prname = TString::Format("Y[%d] slice at X = ", i);

         peakCanvas = new TCanvas(prname);
         gPad->SetGridx();
         gPad->SetGridy();
         gPad->SetTickx();
         gPad->SetTicky();
         peakCanvas->ToggleEventStatus();

         int xbin = int(xpeaks[i]*idp);
         printf("i = %d xbin = %d\n" , i, xbin);

         yproj = dose2d->ProjectionY(TString::Format("X = %f mm", xpeaks[i]), xbin, xbin, "d");
         yproj->SetFillColor(9);
         yproj->SetTitle(prname + TString::Format(" %f mm", xpeaks[i]));

         //gPad->SetTitle(prname + TString::Format(" %f mm", xpeaks[i]));
         //yproj->Draw();
         yproj->GetXaxis()->SetTitle("in mm");
         yproj->GetYaxis()->SetTitle("in Gy");
         //peakCanvas->SetEditable(kFALSE);

         //peakCanvas->cd(i+2);

         prname = TString::Format("X[%d] slice at Y = ", i);

         peakCanvas = new TCanvas(prname);
         gPad->SetGridx();
         gPad->SetGridy();
         gPad->SetTickx();
         gPad->SetTicky();
         peakCanvas->ToggleEventStatus();

         int ybin = int(ypeaks[i]*idp);
         printf("i = %d ybin = %d\n" , i, ybin);

         xproj = dose2d->ProjectionX(TString::Format("Y = %f mm", ypeaks[i]), ybin, ybin, "d");
         xproj->SetFillColor(4);
         xproj->SetTitle(prname + TString::Format(" %f mm", ypeaks[i]));

         //gPad->SetTitle(prname + TString::Format(" %f mm", ypeaks[i]));
         //xproj->Draw();
         xproj->GetXaxis()->SetTitle("in mm");
         xproj->GetYaxis()->SetTitle("in Gy");
         //peakCanvas->SetEditable(kFALSE);
      }
   }

   peakCanvas->Update();
   rootFile->Write();
   
   printf("Stop\n");


   printf("\n\n==========================================================================\n");

   double maxE = dose2d->GetMaximum();
   double halfE = maxE*0.5;
   TAxis *ax = dose1d->GetXaxis();
   int b2 = ax->FindBin(maxE);
   int b1 = ax->FindBin(halfE);
   printf("Full integral = %E.    Half integral = %E\n", dose2d->Integral(), dose1d->Integral(b1, b2));
   int energy = 100; //100 MeV
   //float dE1 = 7.2;  // 7.2 MeV/cm at 100 MeV
   //float dE2 = 9.57; // 9.57 MeV/cm

   double iddp = (300./2.54)*(300./2.54); // (pixels per cantimeter)**2
   double intensity = dose2d->Integral()/iddp/dEdx_ICRU49[energy]/10.*6.24e9; // dEdx_ICRU49 array -
stopping power of protons in water in Mev/mm
   printf("Total Intensity = %E at proton Energy = %d (MeV)\n", intensity, energy);
}

//______________________________________________________________________________
void Intensity(int energy = 100)
{
   // calculates total proton beam intensity based on total dose field integral

   printf("\n\n==========================================================================\n");
   double maxE = dose2d->GetMaximum();
   double halfE = maxE*0.5;
   TAxis *ax = dose1d->GetXaxis();
   int b2 = ax->FindBin(maxE);
   int b1 = ax->FindBin(halfE);
   printf("Full integral = %E.    Half integral = %E\n", dose2d->Integral(),  dose1d->Integral(b1, b2));

   double iddp = (300./2.54)*(300./2.54); // (pixels per cantimeter)**2  
   double intensity = dose2d->Integral()/iddp/dEdx_ICRU49[energy]/10.*6.24e9; // dEdx_ICRU49 array -
stopping power of protons in water in Mev/mm
   printf("Total Intensity = %E at proton Energy = %d (MeV)\n", intensity, energy);
}
Timur Pocheptsov | 3 Jul 14:04 2012
Picon
Picon

RE: Tiff issue

Well, your first try was already successful :))) - that's why I was able to see exe file attached :)
________________________________________
From: Valeri Onuchin
Sent: 03 July 2012 13:53
To: ernst <at> lanl.gov; roottalk (Mailing list for ROOT users.)
Subject: RE: Tiff issue

one more try.

________________________________________
From: Valeri Onuchin
Sent: 03 July 2012 13:29
To: Esch, Shorty; roottalk (Mailing list for ROOT users.)
Subject: RE: Tiff issue

Ooops failed to send because of attachment size limit.

________________________________________
From: Valeri Onuchin
Sent: 03 July 2012 13:24
To: Esch, Shorty; roottalk (Mailing list for ROOT users.)
Subject: RE: Tiff issue

Hi Shorty,
the tiff format is not supported yet.
If you are working under windows
the attached codes demonstrate
- how to read a 16bit TIFF fie, convert it to "raw" format
- how to read the produced "raw" file in ROOT and
  fill 2D hist with scaled red component.

HTH.

Regards. Valeriy

________________________________________
From: Esch, Shorty [ernst <at> lanl.gov]
Sent: 02 July 2012 19:02
To: roottalk (Mailing list for ROOT users.)
Subject: Tiff issue

Hi, has anybody a way to read a TIFF file grayscale into a 2d histogram in root?
Thanks,
Shorty

Jakub Čermák | 3 Jul 14:42 2012
Picon

MLP in TMVA

Hello,

I’m trying to use the MLP neural network in TMVA package. It’s working on some data, but not with some other (usually some more complex data) - I get values like this:

<Neuron NSynapses="7">

          -1.#IND000000000000e+000 -1.#IND000000000000e+000 -1.#IND000000000000e+000 -1.#IND000000000000e+000 -1.#IND000000000000e+000 -1.#IND000000000000e+000 -1.#IND000000000000e+000

</Neuron>

 

What can be the problem? Why the values are-1.#IND000000000000e+000?

My code is here: http://pastebin.com/7SpqXEaR

 

Thanks for some advice, I’m becoming really desperate about this

 

 

S pozdravem / Best regards

 

Jakub Čermák

 


Gmane