Context Navigation

source: Sophya/trunk/SophyaExt/FitsIOServer/fitsfile.cc@ 1143

Visit:

Last change on this file since 1143 was 1143, checked in by ansari, 25 years ago
normalisation des DVList
File size: 32.3 KB

Rev	Line
[949]	1	#include "machdefs.h"
	2	#include <stdlib.h>
[839]	3	#include "fitsfile.h"
	4	#include "pexceptions.h"
	5	#include "strutil.h"
[903]	6	#include "anydataobj.h"
	7	#include "fitsspherehealpix.h"
[1136]	8
	9
	10	void FitsIOHandler::Read(char flnm[],int hdunum)
[839]	11	{
[1136]	12	// FitsInFile ifts(flnm, hdunum);
	13	FitsInFile ifts(flnm);
	14	Read(ifts, hdunum);
[839]	15	}
[1136]	16	void FitsIOHandler::Read(FitsInFile& is, int hdunum)
	17	{
	18	if (hdunum == 0) is.moveToFollowingHeader();
	19	else is.ReadFInit(hdunum);
	20	ReadFromFits(is);
	21	}
	22
	23
	24	void FitsIOHandler::Write(char flnm[], bool OldFile)
	25
	26	{
	27
	28	FitsOutFile of(flnm, OldFile);
	29	Write(of);
	30	}
	31
	32	void FitsIOHandler::Write(FitsOutFile& os)
	33	{
	34	WriteToFits(os);
	35	}
	36
	37
	38
[839]	39
	40	FitsFile::~FitsFile()
	41	{
	42	int status = 0;
[903]	43	if( fptr_ != NULL)
	44	{
	45	fits_close_file(fptr_,&status);
	46	delete fptr_;
	47	}
[839]	48	if( status ) printerror( status );
	49	}
[903]	50
[1136]	51
	52	void FitsFile::printerror(int &status)
	53	//*****************************************************/
	54	//* Print out cfitsio error messages and exit program */
	55	//*****************************************************/
[1045]	56	{
[1136]	57	if( status )
	58	{
	59	fits_report_error(stderr,status);
	60	throw IOExc("FitsFile:: error FITSIO status");
	61	}
	62	return;
	63	}
[903]	64
[1136]	65	void FitsFile::printerror(int& status, char* texte)
	66	//*****************************************************/
	67	//* Print out cfitsio error messages and exit program */
	68	//*****************************************************/
	69	{
	70	// print out cfitsio error messages and exit program
	71	// print error report
	72	fits_report_error(stderr, status);
	73	cout << " erreur:: " << texte << endl;
	74	throw IOExc("FitsFile:: error FITSIO status");
	75	}
	76
	77	void FitsFile::ResetStatus(int& status)
	78	{
	79	fits_status_ = status;
	80	status = 0;
	81	}
	82
	83	string FitsFile::getErrStatus(int status)
	84	{
	85	char text[31];
	86	fits_get_errstatus(status, text);
	87	return string(text);
	88	}
	89
	90	FitsInFile::FitsInFile()
	91	{
	92	InitNull();
	93	}
	94	//FitsInFile::FitsInFile(char flnm[], int hdunum)
	95	FitsInFile::FitsInFile(char flnm[])
	96	{
	97	InitNull();
	98	int status = 0;
	99	fits_open_file(&fptr_,flnm,READONLY,&status);
	100	if( status ) printerror( status );
	101	// ReadFInit(flnm, hdunum);
	102	}
	103
	104
	105	void FitsInFile::InitNull()
	106	{
[1045]	107	bitpix_ = 0;
	108	naxis_ = 0;
	109	nbData_ = 0;
	110	nrows_ = 0;
	111	nbcols_ = 0;
	112	naxisn_.clear();
	113	repeat_.clear();
	114	noms_.clear();
	115	taille_des_chaines_.clear();
	116	dvl_.Clear();
	117	}
	118
	119
[1136]	120	int FitsInFile::NbBlocks(char flnm[])
[903]	121	{
	122	int status = 0;
	123	int nbhdu = 0;
	124	fitsfile* fileptr;
	125	fits_open_file(&fileptr,flnm,READONLY,&status);
	126	if( status ) printerror( status, "NbBlocks: erreur ouverture fichier" );
	127	fits_get_num_hdus(fileptr, &nbhdu, &status);
	128	fits_close_file(fileptr,&status);
	129	return nbhdu;
	130	}
	131
[1136]	132	void FitsInFile::getBlockType(char flnm[], int hdunum, string& typeOfExtension, int& naxis, vector<int>& naxisn, string& dataType, DVList& dvl )
[903]	133	{
	134	int status = 0;
	135	fitsfile* fileptr;
	136	fits_open_file(&fileptr,flnm,READONLY,&status);
	137	if( status ) printerror( status, "getBlockType: erreur ouverture fichier" );
	138	// move to the specified HDU number
	139	int hdutype = 0;
	140	fits_movabs_hdu(fileptr,hdunum,&hdutype,&status);
	141	if( status ) printerror( status,"getBlockType: erreur movabs");
	142	if(hdutype == IMAGE_HDU)
	143	{
	144	typeOfExtension = "IMAGE";
	145	int bitpix;
	146	GetImageParameters (fileptr, bitpix, naxis, naxisn);
	147	if(bitpix == DOUBLE_IMG) dataType = "double";
	148	else
	149	if(bitpix == FLOAT_IMG) dataType = "float";
	150	else
	151	if(bitpix == LONG_IMG \|\| bitpix == SHORT_IMG ) dataType = "int";
	152	else
	153	{
	154	cout << " bitpix= " << bitpix << endl;
	155	throw PException(" FitsFile::getBlockType : unsupprted FITS data type");
	156	}
	157
	158	}
	159	else
	160	if(hdutype == ASCII_TBL \|\| hdutype == BINARY_TBL)
	161	{
	162	int nrows = 0;
	163	vector<string> noms;
	164	vector<char> types;
	165	vector<int> taille_des_chaines;
[971]	166	GetBinTabParameters(fileptr, naxis, nrows, naxisn, noms, types, taille_des_chaines);
	167	int k;
	168	for (k=0; k< naxisn.size(); k++) naxisn[k] *= nrows;
[903]	169	if(hdutype == ASCII_TBL)
	170	{
	171	typeOfExtension = "ASCII_TBL";
	172	dataType = "ASCII";
	173	}
	174	else
	175	{
	176	typeOfExtension = "BINARY_TBL";
	177	if(types[0] == 'D') dataType = "double";
	178	else
	179	if(types[0] == 'E') dataType = "float";
	180	else
	181	if(types[0] == 'I' ) dataType = "int";
	182	else
	183	if(types[0] == 'S' ) dataType = "char*";
	184	else
	185	{
	186	cout << " types[0]= " << types[0] << endl;
	187	throw PException(" FitsFile::getBlockType : unsupprted FITS data type");
	188	}
	189	}
	190	}
	191	else
	192	{
	193	cout << " hdutype= " << hdutype << endl;
	194	throw IOExc("FitsFile::getBlockType: this HDU type is unknown");
	195	}
	196
[971]	197	KeywordsIntoDVList(fileptr, dvl, hdunum);
[903]	198	fits_close_file(fileptr,&status);
	199	}
	200
[1136]	201	void FitsInFile::GetImageParameters (fitsfile* fileptr,int& bitpix,int& naxis,vector<int>& naxisn)
[839]	202	{
[1136]	203	int hdunum=0;
	204	cout << " Reading a FITS image in HDU : " << fits_get_hdu_num(fileptr,&hdunum) << endl;
	205	int status= 0;
	206
	207	// bits per pixels
	208	fits_read_key(fileptr,TINT,"BITPIX",&bitpix,NULL,&status);
	209	if( status ) printerror( status );
	210
	211	// number of dimensions in the FITS array
	212	naxis= 0;
	213	fits_read_key(fileptr,TINT,"NAXIS",&naxis,NULL,&status);
	214	if( status ) printerror( status );
	215	// read the NAXISn keywords to get image size
	216	long* naxes = new long[naxis] ;
	217	int nfound;
	218	fits_read_keys_lng(fileptr,"NAXIS",1,naxis,naxes,&nfound,&status);
	219	if( status ) printerror( status );
	220	if (nfound != naxis )
	221	cout << " WARNING : " << nfound << " axes found, expected naxis= " << naxis << endl;
	222	int k;
	223	for (k=0; k<naxis; k++)
	224	{
	225	naxisn.push_back( (int)naxes[k] );
	226	}
	227	delete [] naxes;
[971]	228	}
[839]	229
[971]	230
[1136]	231	void FitsInFile::ReadFInit(int hdunum)
[1045]	232	{
[1143]	233	InitNull();
[1136]	234	// int status = 0;
	235
	236	// fits_open_file(&fptr_,flnm,READONLY,&status);
	237	// if( status ) printerror( status );
	238
[1045]	239	int status = 0;
[1136]	240
[1045]	241	if (hdunum <= 1)
	242	{
	243	hdunum_ = 1;
	244	// presence of image ?
	245	int naxis= 0;
	246	fits_read_key(fptr_,TINT,"NAXIS",&naxis,NULL,&status);
	247	if( status ) printerror( status );
	248	if (naxis > 0 ) // there is an image
	249	{
	250	hdutype_ = IMAGE_HDU;
	251	GetImageParameters (fptr_, bitpix_, naxis_, naxisn_);
	252	nbData_ = 1;
	253	int k;
	254	for (k=0; k<naxis_; k++) if (naxisn_[k] > 0) nbData_ *= naxisn_[k];
	255	KeywordsIntoDVList(fptr_, dvl_,hdunum_);
	256	}
	257	else
	258	{
	259	throw PException(" first header : no image, probably error in hdunum");
	260	}
	261	//
	262	}
	263	else
	264	{
	265	hdunum_ = hdunum-1;
	266	int hdutype;
	267	fits_movabs_hdu(fptr_,hdunum_,&hdutype,&status);
	268	if( status ) printerror( status,":FitsFile::ReadF : erreur movabs");
	269	moveToFollowingHeader();
	270	}
	271	}
	272
[1143]	273	DVList FitsInFile::DVListFromPrimaryHeader() const
	274	{
	275	int status;
	276	DVList dvl;
	277	KeywordsIntoDVList(fptr_, dvl, 1);
	278	int hdutype = 0;
	279	if (hdunum_ > 0) fits_movabs_hdu(fptr_,hdunum_,&hdutype,&status);
	280	return dvl;
	281	}
[1136]	282
	283
	284	void FitsInFile::moveToFollowingHeader()
[971]	285	{
	286	int status = 0;
[839]	287	int hdutype;
[1143]	288	InitNull();
[971]	289	fits_movrel_hdu(fptr_, 1,&hdutype,&status);
	290	if( status ) printerror( status," lecture du header suivant" );
	291	hdunum_++;
[839]	292	hdutype_= hdutype;
	293	if(hdutype_ == IMAGE_HDU)
	294	{
[903]	295	GetImageParameters (fptr_, bitpix_, naxis_, naxisn_);
	296	nbData_ = 1;
[971]	297	int k;
	298	for (k=0; k<naxis_; k++) if (naxisn_[k] > 0) nbData_ *= naxisn_[k];
	299	KeywordsIntoDVList(fptr_, dvl_,hdunum_);
[839]	300	}
	301	if(hdutype_ == ASCII_TBL \|\| hdutype_ == BINARY_TBL)
	302	{
[903]	303	GetBinTabParameters(fptr_,nbcols_, nrows_,repeat_, noms_, types_, taille_des_chaines_);
[971]	304	KeywordsIntoDVList(fptr_, dvl_, hdunum_);
[839]	305	}
[971]	306	}
[839]	307
[1136]	308
	309	int FitsInFile::NbColsFromFits() const
[839]	310	{
[1136]	311	if(hdutype_ == BINARY_TBL) return nbcols_;
[971]	312	else
[1136]	313	if(hdutype_ == ASCII_TBL \|\| hdutype_ == IMAGE_HDU) return 1;
[839]	314	else
	315	{
[1136]	316	cout << " hdutype= " << hdutype_ << endl;
	317	throw PException("FitsFile::NbColsFromFits, this HDU is unknown");
[839]	318	}
	319	}
	320
[1136]	321	int FitsInFile::NentriesFromFits(int nocol) const
[839]	322	{
[1136]	323	if(hdutype_ == BINARY_TBL ) return nrows_*repeat_[nocol];
	324	else
	325	if(hdutype_ == ASCII_TBL) return nrows_;
	326	else
	327	if(hdutype_ == IMAGE_HDU) return nbData_;
	328	else
[839]	329	{
[1136]	330	cout << "hdutype= " << hdutype_ << endl;
	331	throw PException("FitsFile::NentriesFromFits, this HDU is unknown");
[839]	332	}
	333	}
	334
[1136]	335	char FitsInFile::ColTypeFromFits(int nocol) const
[839]	336	{
[1136]	337	if(hdutype_ != ASCII_TBL && hdutype_ != BINARY_TBL)
[839]	338	{
[1136]	339	throw IOExc("FitsFile::TypeFromFits, this HDU is not an ASCII table nor a binary table");
[839]	340	}
[1136]	341	return types_[nocol];
[839]	342	}
[1136]	343	string FitsInFile::ColNameFromFits(int nocol) const
[839]	344	{
[1136]	345	if(hdutype_ != ASCII_TBL && hdutype_ != BINARY_TBL)
[1045]	346	{
[1136]	347	throw IOExc("FitsFile::TypeFromFits, this HDU is not an ASCII table nor a binary table");
[1045]	348	}
[1136]	349	return noms_[nocol];
[839]	350	}
	351
[1136]	352	int FitsInFile::ColStringLengthFromFits(int nocol) const
[839]	353	{
[1136]	354	if(hdutype_ != ASCII_TBL && hdutype_ != BINARY_TBL)
	355	{
	356	throw IOExc("FitsFile::TypeFromFits, this HDU is not an ASCII table nor a binary table");
	357	}
	358	int index=-1;
	359	int k;
	360	for (k=0; k<=nocol; k++)
	361	{
	362	if (types_[k] == 'S') index++;
	363	}
	364	return taille_des_chaines_[index];
[839]	365	}
[1136]	366	void FitsInFile::GetBinTabLine(int NoLine, double* ddata, float* fdata, int* idata, char ** cdata)
[839]	367	{
	368	int status= 0;
[1136]	369	int anull;
	370	double dnull= 0.;
	371	float fnull= 0.;
	372	int inull= 0;
	373	char* cnull= "";
	374	int dcount = 0.;
	375	int fcount = 0.;
	376	int icount = 0;
	377	int ccount =0;
	378	int ncol;
	379	long nels=1;
	380	for (ncol=0; ncol<nbcols_; ncol++)
[861]	381	{
[1136]	382	switch (types_[ncol])
	383	{
	384	case 'D' :
	385	fits_read_col(fptr_,TDOUBLE,ncol+1,NoLine+1,1,1,&dnull,&ddata[dcount++],&anull,&status);
	386	break;
	387	case 'E' :
	388	fits_read_col(fptr_,TFLOAT,ncol+1,NoLine+1,1,1,&fnull,&fdata[fcount++],&anull,&status);
	389	break;
	390	case 'I' :
	391	fits_read_col(fptr_,TINT,ncol+1,NoLine+1,1,1,&inull,&idata[icount++],
	392	&anull,&status);
	393	break;
	394	case 'S' :
	395	fits_read_col(fptr_,TSTRING,ncol+1,NoLine+1,1,1,cnull,&cdata[ccount++],&anull,&status);
	396	break;
	397	}
	398	if (status)
	399	{
	400	ResetStatus(status);
	401	break;
	402	}
[861]	403	}
[903]	404	}
[839]	405
[1136]	406	void FitsInFile::GetBinTabLine(int NoLine, float* fdata)
[903]	407	{
[1136]	408	int status= 0;
	409	int anull;
	410	float fnull= 0.;
	411	long nels=1;
	412	int ncol;
	413	for (ncol=0; ncol<nbcols_; ncol++)
[861]	414	{
[1136]	415	fits_read_col(fptr_,TFLOAT,ncol+1,NoLine+1,1,1,&fnull,&fdata[ncol],&anull,&status);
	416	if (status)
[903]	417	{
[1136]	418	ResetStatus(status);
	419	break;
[903]	420	}
[1136]	421	}
	422	}
[839]	423
[903]	424
[1136]	425	void FitsInFile::GetBinTabFCol(double* valeurs,int nentries, int NoCol) const
[839]	426	{
	427	int status= 0;
	428	int DTYPE;
	429	long repeat,width;
	430	fits_get_coltype(fptr_, NoCol+1,&DTYPE,&repeat,&width,&status);
	431	if( DTYPE != TDOUBLE)
	432	{
[1045]	433	if (DTYPE == TFLOAT) cout << " WARNING: reading double from float : conversion will be made by fitsio library" << endl;
	434	else
	435	throw IOExc("FitsFile::GetBinTabFCol, tentative de lecture non double");
[839]	436	}
	437	long nels=nentries;
[971]	438	int anull;
[839]	439	// no checking for undefined pixels
[971]	440	double dnull= 0.;
	441	// fits_read_key(fptr_,TDOUBLE,"BAD_DATA",&dnull,NULL,&status);
	442	// if (status != 0)
	443	// {
	444	// dnull = -1.6375e30; // default value
	445	// status = 0;
	446	// }
	447	if (nentries != nrows_*repeat)
	448	{
	449	cout << " found " << nentries << " pixels, expected: " << nrows_*repeat << endl;
	450	throw PException(" FitsFile:::GetBinTabFCol ");
	451	}
[839]	452	fits_read_col(fptr_,TDOUBLE,NoCol+1,1,1,nels,&dnull,valeurs,
	453	&anull,&status);
	454	if( status ) printerror( status,"erreur lecture de colonne" );
[971]	455
[839]	456	}
	457
[1136]	458	void FitsInFile::GetBinTabFCol(float* valeurs,int nentries, int NoCol) const
[839]	459	{
	460	int status= 0;
	461	int DTYPE;
	462	long repeat,width;
	463	fits_get_coltype(fptr_, NoCol+1,&DTYPE,&repeat,&width,&status);
	464	if( DTYPE != TFLOAT)
	465	{
[1045]	466	if (DTYPE == TDOUBLE) cout << " WARNING: reading float from double : conversion will be made by fitsio library" << endl;
	467	else
	468	throw IOExc("FitsFile::GetBinTabFCol, tentative de lecture non float");
[839]	469	}
	470	long nels=nentries;
[971]	471	int anull;
[839]	472	// no checking for undefined pixels
	473	float fnull= 0.;
[971]	474	// fits_read_key(fptr_,TFLOAT,"BAD_DATA",&fnull,NULL,&status);
	475	// if (status != 0)
	476	// {
	477	// fnull = -1.6375e30; // default value
	478	// status = 0;
	479	// }
	480	if (nentries != nrows_*repeat)
	481	{
	482	cout << " found " << nentries << " pixels, expected: " << nrows_*repeat << endl;
	483	throw PException(" FitsFile:::GetBinTabFCol ");
	484	}
[839]	485	fits_read_col(fptr_,TFLOAT,NoCol+1,1,1,nels,&fnull,valeurs,
	486	&anull,&status);
	487	if( status ) printerror( status,"erreur lecture de colonne" );
	488	}
[1136]	489
	490	void FitsInFile::GetBinTabFCol(int* valeurs,int nentries, int NoCol) const
[839]	491	{
	492	int status= 0;
	493	int DTYPE;
	494	long repeat,width;
	495	fits_get_coltype(fptr_, NoCol+1,&DTYPE,&repeat,&width,&status);
	496	if( DTYPE != TLONG && DTYPE != TINT && DTYPE != TSHORT )
	497	{
	498	throw IOExc("FitsFile::GetBinTabFCol, tentative de lecture non entier");
	499	}
	500	long nels=nentries;
	501	// no checking for undefined pixels
	502	int anull;
	503	int inull= 0;
[971]	504	// fits_read_key(fptr_,TINT,"BAD_DATA",&inull,NULL,&status);
	505	// if (status != 0)
	506	// {
	507	// inull = -999999; // default value
	508	// status = 0;
	509	// }
	510	if (nentries != nrows_*repeat)
	511	{
	512	cout << " found " << nentries << " pixels, expected: " << nrows_*repeat << endl;
	513	throw PException(" FitsFile:::GetBinTabFCol ");
	514	}
[839]	515	fits_read_col(fptr_,TINT,NoCol+1,1,1,nels,&inull,valeurs,
	516	&anull,&status);
	517	if( status ) printerror( status,"erreur lecture de colonne" );
	518	}
[1136]	519
	520	void FitsInFile::GetBinTabFCol(char** valeurs, int nentries, int NoCol) const
[839]	521	{
	522	int status= 0;
	523	int DTYPE;
	524	long repeat,width;
	525	fits_get_coltype(fptr_, NoCol+1,&DTYPE,&repeat,&width,&status);
	526	if( DTYPE != TSTRING )
	527	{
	528	throw IOExc("FitsFile::GetBinTabFCol, tentative de lecture non float");
	529	}
	530	long nels=nentries;
	531	// no checking for undefined pixels
	532	int anull;
[971]	533	char* cnull= "";
	534	if (nentries != nrows_*repeat/width)
	535	{
	536	cout << " found " << nentries << " pixels, expected: " << nrows_*repeat/width << endl;
	537	throw PException(" FitsFile:::GetBinTabFCol ");
	538	}
[839]	539	long frow=1;
	540	long felem=1;
	541	fits_read_col(fptr_,TSTRING,NoCol+1,frow,felem,nels,cnull,valeurs,
	542	&anull,&status);
	543	if( status ) printerror( status,"erreur lecture de colonne" );
	544	}
[1045]	545
[1136]	546	void FitsInFile::GetSingleColumn(double* map, int nentries) const
	547	{
	548	int status = 0;
	549	if(hdutype_ == IMAGE_HDU)
[1045]	550	{
[1136]	551
	552	if(bitpix_ != DOUBLE_IMG)
[1047]	553	{
[1136]	554	cout << " The data type on fits file is not double...";
	555	cout << " Conversion to double achieved by cfitsio lib" << endl;
[1047]	556	}
[1136]	557
	558	// no checking for undefined pixels
	559	int anull;
	560	double dnull= 0.;
	561
	562	long nels= nentries;
	563	fits_read_img(fptr_,TDOUBLE,1,nels,&dnull,map,&anull,&status);
	564	if( status ) printerror( status );
[1045]	565	}
[1136]	566	else
	567	if(hdutype_ == ASCII_TBL \|\| hdutype_ == BINARY_TBL)
	568	{
	569	GetBinTabFCol(map,nentries, 0);
	570	}
	571	else
	572	{
	573	cout << " hdutype= " << hdutype_ << endl;
	574	throw IOExc("FitsFile::GetSingleColumn this HDU is unknown");
	575	}
[1045]	576	}
	577
[1136]	578	void FitsInFile::GetSingleColumn(float* map, int nentries) const
[1045]	579	{
[1136]	580	int status = 0;
	581	if(hdutype_ == IMAGE_HDU)
[1045]	582	{
[1136]	583	if(bitpix_ != FLOAT_IMG)
[1047]	584	{
[1136]	585	cout << " The data type on fits file is not float ";
	586	cout << " Conversion to float achieved by cfitsio lib" << endl;
[1047]	587	}
[1136]	588	// no checking for undefined pixels
	589	int anull;
	590	float fnull= 0.;
	591
	592	long nels= nentries;
	593	fits_read_img(fptr_,TFLOAT,1,nels,&fnull, map,&anull,&status);
	594	if( status ) printerror( status );
[1045]	595	}
[839]	596	else
[1136]	597	if(hdutype_ == ASCII_TBL \|\| hdutype_ == BINARY_TBL)
	598	{
	599	GetBinTabFCol(map,nentries, 0);
	600	}
[839]	601	else
	602	{
[1136]	603	cout << " hdutype= " << hdutype_ << endl;
	604	throw IOExc("FitsFile::GetSingleColumn this HDU is unknown");
[839]	605	}
	606	}
	607
[1136]	608	void FitsInFile::GetSingleColumn( int* map, int nentries) const
[839]	609	{
[1136]	610	int status = 0;
	611	if(hdutype_ == IMAGE_HDU)
[839]	612	{
[1136]	613	if(bitpix_ != LONG_IMG)
	614	{
	615	cout << " The data type on fits file is not int ";
	616	cout << " Conversion to float achieved by cfitsio lib" << endl;
	617	}
	618	// no checking for undefined pixels
	619	int anull;
	620	float fnull= 0.;
	621
	622	long nels= nentries;
	623	fits_read_img(fptr_,TINT,1,nels,&fnull,map,&anull,&status);
	624	if( status ) printerror( status );
[839]	625	}
	626	else
[1136]	627	if(hdutype_ == ASCII_TBL \|\| hdutype_ == BINARY_TBL)
	628	{
	629	GetBinTabFCol(map,nentries, 0);
	630	}
[839]	631	else
[1136]	632	{
	633	cout << " hdutype= " << hdutype_ << endl;
	634	throw IOExc("FitsFile::GetSingleColumn this HDU is unknown");
	635	}
[839]	636	}
	637
[1136]	638	void FitsInFile::GetBinTabParameters(fitsfile* fileptr, int& nbcols, int& nrows,
[903]	639	vector<int>& repeat,
	640	vector<string>& noms,
	641	vector<char>& types,
	642	vector<int>& taille_des_chaines)
[839]	643	{
	644	int status= 0;
[903]	645	int hdunum=0;
	646	int hdutype=0;
	647	fits_get_hdu_num(fileptr,&hdunum);
	648	fits_get_hdu_type(fileptr, &hdutype, &status);
	649
	650	if(hdutype != ASCII_TBL && hdutype != BINARY_TBL)
[839]	651	{
[903]	652	throw IOExc("FitsFile::GetBinTabParameters this HDU is not an ASCII table nor a binary table");
[839]	653	}
[903]	654	if(hdutype == ASCII_TBL)
	655	cout << " Reading a FITS ascii table in HDU : " << hdunum << endl;
	656	if(hdutype == BINARY_TBL)
	657	cout << " Reading a FITS binary table in HDU : " << hdunum << endl;
[839]	658
	659	// get the number of columns
[903]	660	fits_get_num_cols(fileptr, &nbcols,&status);
[839]	661	if( status ) printerror( status );
	662
	663	// get the number of rows
	664	long naxis2= 0;
[903]	665	fits_get_num_rows(fileptr,&naxis2,&status);
[839]	666	if( status ) printerror( status );
[903]	667	nrows = (int)naxis2;
[839]	668
	669	// get the datatype, names and the repeat count
[903]	670	noms.clear();
	671	noms.reserve(nbcols);
	672	types.clear();
	673	types.reserve(nbcols);
	674	repeat.clear();
	675	repeat.reserve(nbcols);
	676	taille_des_chaines.clear();
[839]	677	char *ttype = new char[nbcols];
[923]	678	int ii;
	679	for (ii=0; ii < nbcols; ii++) ttype[ii]=new char[FLEN_VALUE];
[839]	680	int nfound;
[903]	681	fits_read_keys_str(fileptr, "TTYPE",1,nbcols,ttype,&nfound, &status);
[839]	682	if( status ) printerror( status,"erreur lecture des noms de colonne");
	683	int rept=0;
[923]	684	for(ii = 0; ii < nbcols; ii++)
[839]	685	{
	686	int DTYPE;
	687	long width;
[903]	688	long repete = 0;
	689	fits_get_coltype(fileptr,ii+1,&DTYPE,&repete,&width,&status);
[839]	690	if( status ) printerror( status,"erreur lecture type de colonne");
[903]	691	rept = repete;
	692	noms.push_back(string(ttype[ii]));
[839]	693	switch (DTYPE)
	694	{
	695	case TDOUBLE :
[903]	696	types.push_back('D');
[839]	697	break;
	698	case TFLOAT :
[903]	699	types.push_back('E');
[839]	700	break;
	701	case TLONG :
[903]	702	types.push_back('I');
[839]	703	break;
	704	case TINT :
[903]	705	types.push_back('I');
[839]	706	break;
	707	case TSHORT :
[903]	708	types.push_back('I');
[839]	709	break;
	710	case TSTRING :
[903]	711	types.push_back('S');
	712	taille_des_chaines.push_back(width);
[839]	713	rept/=width;
	714	break;
	715	default :
	716	cout << " field " << ii+1 << " DTYPE= " << DTYPE << endl;
	717	throw IOExc("FitsFile:: unknown type of field");
	718	}
[903]	719	repeat.push_back(rept);
[1136]	720	}
	721	for (ii=0; ii < nbcols; ii++) delete [] ttype[ii];
	722	delete [] ttype;
	723	}
	724
	725	void FitsInFile::KeywordsIntoDVList(fitsfile* fileptr, DVList& dvl, int hdunum)
	726	{
	727	int status = 0;
	728	int hdutype;
	729	fits_movabs_hdu(fileptr,hdunum,&hdutype,&status);
	730	if( status ) printerror( status,":KeywordsIntoDVList : erreur movabs");
	731	// get number of keywords
	732	int nkeys,keypos;
	733	fits_get_hdrpos(fileptr,&nkeys,&keypos,&status);
	734	if( status ) printerror( status );
	735
	736	// put keywords in a DVList object
	737	char keyname[LEN_KEYWORD]= "";
	738	char strval[FLEN_VALUE]= "";
	739	char dtype;
	740	char card[FLEN_CARD];
	741	char *comkey = "COMMENT";
[1143]	742	char comment[FLEN_COMMENT];
[1136]	743
	744	// shift with the number of mandatory keywords
[1143]	745	// int num= 8;
	746	int num= 0;
	747	// primary header
	748	if (hdunum == 1)
	749	{
	750	// read NAXIS
	751	int naxis=0;
	752	fits_read_key(fileptr,TINT,"NAXIS",&naxis,NULL,&status);
	753	// number of mandatory keywords
	754	num = naxis+3;
	755	}
	756	// extensions
	757	else
	758	{
	759	if (hdutype == IMAGE_HDU)
	760	{
	761	// read NAXIS
	762	int naxis=0;
	763	fits_read_key(fileptr,TINT,"NAXIS",&naxis,NULL,&status);
	764	// number of mandatory keywords
	765	num = naxis+5;
	766	}
	767	else
	768	if(hdutype == ASCII_TBL \|\| hdutype == BINARY_TBL)
	769	{
	770	// number of mandatory keywords
	771	num = 8;
	772	}
	773	}
[1136]	774	int j;
	775	for(j = num+1; j <= nkeys; j++)
	776	{
	777	fits_read_keyn(fileptr,j,card,strval,NULL,&status);
	778	if(status) printerror(status);
	779
	780	strncpy(keyname,card,LEN_KEYWORD-1);
	781
	782	if(strncmp(keyname,comkey,LEN_KEYWORD-1) != 0 && strlen(keyname) != 0
	783	&& strlen(strval) != 0)
	784	{
	785	fits_get_keytype(strval,&dtype,&status);
	786	if(status) printerror(status);
	787
	788	strip(keyname, 'B',' ');
	789	strip(strval, 'B',' ');
	790	strip(strval, 'B','\'');
	791
	792	switch( dtype )
	793	{
	794	case 'C':
[1143]	795	fits_read_key(fileptr,TSTRING,keyname,strval,comment,&status);
	796	dvl[keyname]= strval;
	797	dvl.SetComment(keyname, comment);
[1136]	798	break;
	799	case 'I':
	800	int ival;
[1143]	801	fits_read_key(fileptr,TINT,keyname,&ival,comment,&status);
[1136]	802	dvl[keyname]= (int_4) ival; // Portage mac DY
[1143]	803	dvl.SetComment(keyname, comment);
[1136]	804	break;
	805	case 'L':
	806	int ilog;
[1143]	807	fits_read_key(fileptr,TLOGICAL,keyname,&ilog,comment,&status);
[1136]	808	dvl[keyname]= (int_4) ilog;
[1143]	809	dvl.SetComment(keyname, comment);
[1136]	810	break;
	811	case 'F':
	812	double dval;
[1143]	813	fits_read_key(fileptr,TDOUBLE,keyname,&dval,comment,&status);
[1136]	814	dvl[keyname]= dval;
[1143]	815	dvl.SetComment(keyname, comment);
[1136]	816	break;
	817	}
	818
	819	}
[839]	820	}
[1136]	821	// dvl_.Print();
	822	}
	823
	824	FitsOutFile::FitsOutFile()
	825	{
	826	InitNull();
[903]	827	}
[839]	828
[1136]	829	FitsOutFile::FitsOutFile(char flnm[], bool OldFile)
	830	{
[839]	831
[1136]	832	InitNull();
	833	int status = 0;
[839]	834
[1136]	835	// create new FITS file
	836	if (!OldFile)
	837	{
	838	fits_create_file(&fptr_,flnm,&status);
	839	if( status ) printerror(status,"file already exists");
	840	}
	841	else
	842	{
	843	fits_open_file(&fptr_,flnm,READWRITE,&status);
	844	if( status ) printerror(status,"file does not exist");
	845	fits_get_num_hdus(fptr_, &hdunum_, &status);
	846	int hdutype;
	847	fits_movabs_hdu(fptr_,hdunum_,&hdutype,&status);
	848	if( status ) printerror( status,":FitsFile::WriteF : erreur movabs");
	849
	850	}
	851
	852	}
	853
	854
	855
[1143]	856	void FitsOutFile::makeHeaderImageOnFits(char type, int nbdim, int* naxisn, DVList &dvl)
[1136]	857	{
	858	int status = 0;
	859	long naxis = nbdim;
	860	long* naxes = new long[nbdim];
	861	if (hdunum_ == 1)
	862	{
[1143]	863	if (imageOnPrimary_ == false)
	864	{
	865	fits_create_img(fptr_,FLOAT_IMG,0,naxes,&status);
	866	writeSignatureOnFits(); }
[1136]	867	else hdunum_--;
	868	}
	869	int k;
	870	for (k=0; k< nbdim; k++) naxes[k] = (long)naxisn[k];
	871	if (type == 'D')
	872	fits_create_img(fptr_,DOUBLE_IMG,naxis,naxes,&status);
	873	else
	874	if (type == 'E')
	875	fits_create_img(fptr_,FLOAT_IMG,naxis,naxes,&status);
	876	else
	877	if (type == 'I')
	878	fits_create_img(fptr_,LONG_IMG,naxis,naxes,&status);
	879	else
	880	{
	881	cout << " type of data: " << type << endl;
	882	throw PException("FitsFile:::makeHeaderImageOnFits:unprogrammed type of data ");
	883	}
	884
	885	hdunum_++;
[1143]	886
	887	// write supplementary keywords
	888	// dvl.Print();
	889	addKeywordsOfDVList(dvl);
	890
[1136]	891	delete [] naxes;
	892	if( status ) printerror( status, "erreur creation HDU IMAGE" );
	893
	894	}
	895	void FitsOutFile::putImageToFits(int nbData, double* map) const
	896	{
	897	int status = 0;
	898	long npix= nbData;
	899	fits_write_img(fptr_,TDOUBLE,1,npix,map,&status);
	900	if( status ) printerror( status, "erreur ecriture putImageToFits" );
[1143]	901	// writeSignatureOnFits();
[1136]	902	}
	903
	904	void FitsOutFile::putImageToFits(int nbData, float* map) const
	905	{
	906	int status = 0;
	907	long npix= nbData;
	908	fits_write_img(fptr_,TFLOAT,1,npix, map,&status);
	909	if( status ) printerror( status, "erreur ecriture putImageToFits" );
[1143]	910	// writeSignatureOnFits();
[1136]	911
	912	}
	913	void FitsOutFile::putImageToFits( int nbData, int* map) const
	914	{
	915	int status = 0;
	916
	917	long npix= nbData;
	918	fits_write_img(fptr_,TINT,1,npix,map,&status);
	919	if( status ) printerror( status, "erreur ecriture putImageToFits" );
[1143]	920	// writeSignatureOnFits();
[1136]	921	}
	922
	923
	924
	925	void FitsOutFile::makeHeaderBntblOnFits( char* fieldType, char** Noms, int nentries, int tfields, DVList &dvl, char* extname, vector<int> taille_des_chaines)
[839]	926	{
	927	int status = 0;
	928	long nrows;
	929	if (strlen(fieldType) != tfields)
	930	{
	931	cout << " nombre de champs :" << tfields << "nombre de types: " << strlen(fieldType) << endl;
[1136]	932	throw ParmError("FitsFile:: fields and types don't match");
[839]	933
	934	}
	935	char ** ttype= new char*[tfields];
	936	char ** tform= new char*[tfields];
	937	char largeur[FLEN_VALUE];
	938	int noColString=0;
[971]	939	int k;
	940	for (k=0; k<tfields;k++)
[839]	941	{
	942	char format[FLEN_VALUE];
	943
	944	if(nentries < 1024)
	945	{
	946	nrows= nentries;
	947	if (fieldType[k] == 'A')
	948	{
	949	sprintf(largeur,"%d",taille_des_chaines[noColString++]);
	950	strcpy(format,largeur);
	951	}
	952	else strcpy(format,"1");
	953	}
	954	else
	955	{
	956	nrows = nentries/1024;
	957	if(nentries%1024 != 0) nrows++;
	958	if (fieldType[k] == 'A')
	959	{
[1136]	960	char largaux[FLEN_VALUE];
	961	sprintf(largeur,"%d",taille_des_chaines[noColString]);
	962	sprintf(largaux,"%d",1024*taille_des_chaines[noColString]);
	963	noColString++;
	964	strcpy(format, largaux);
[839]	965	}
	966	else strcpy(format,"1024");
	967	}
	968	strncat(format,&fieldType[k],1);
	969	if (fieldType[k] == 'A')
	970	{
	971	strcat(format,largeur);
	972	}
	973	ttype[k]= new char[FLEN_VALUE];
	974	strcpy(ttype[k],Noms[k]);
	975	tform[k]= new char[FLEN_VALUE];
	976	strcpy(tform[k],format);
	977	}
	978	// value of the EXTNAME keyword
	979	char extn[FLEN_VALUE];
	980	strncpy(extn,extname,FLEN_VALUE);
	981
	982	// create a new empty binary table onto the FITS file
	983	// physical units if they exist, are defined in the DVList object
	984	// so the NULL pointer is given for the tunit parameters.
	985	nrows=0;
	986	fits_create_tbl(fptr_,BINARY_TBL,nrows,tfields,ttype,tform,
	987	NULL,extn,&status);
	988	if( status ) printerror( status );
[1026]	989	if ( hdunum_ == 0 ) hdunum_ = 2;
	990	else hdunum_++;
[971]	991	int ii;
	992	for(ii = 0; ii < tfields; ii++)
[839]	993	{
	994	delete [] ttype[ii];
	995	delete [] tform[ii];
	996	}
	997	delete [] ttype;
	998	delete [] tform;
	999	//
	1000	// write supplementary keywords
[1143]	1001	addKeywordsOfDVList(dvl);
[839]	1002
	1003	}
	1004
[1136]	1005	void FitsOutFile::putColToFits(int nocol, int nentries, double* donnees) const
[839]	1006	{
	1007	int status = 0;
[971]	1008	int hdutype;
[839]	1009	fits_movabs_hdu(fptr_,hdunum_,&hdutype,&status);
	1010	if( status ) printerror(status,"putColToFits: le movabs a foire");
	1011	fits_get_hdu_type(fptr_, &hdutype, &status);
[867]	1012	if(hdutype != ASCII_TBL && hdutype != BINARY_TBL)
	1013	{
	1014	cout << " hdunum= " << hdunum_ << " hdutype= " << hdutype << endl;
	1015	throw IOExc("FitsFile::putColToFits, this HDU is not an ASCII table nor a binary table");
	1016	}
[839]	1017	int code;
	1018	long repeat, width;
	1019	fits_get_coltype(fptr_, nocol+1, &code, &repeat,&width, &status);
	1020	if( code != TDOUBLE)
	1021	{
	1022	cout << " WARNING : types don't match (putColToFits) : on fits file= " << code << " to be written= DOUBLE " << endl;
	1023	}
	1024	fits_write_col(fptr_,TDOUBLE,nocol+1,1,1,nentries, donnees ,&status);
	1025	if( status ) printerror( status,"erreur ecriture du fichier fits" );
	1026	}
[1136]	1027	void FitsOutFile::putColToFits(int nocol, int nentries, float* donnees) const
[839]	1028	{
	1029	int status = 0;
	1030	int hdutype;
	1031	fits_movabs_hdu(fptr_,hdunum_,&hdutype,&status);
	1032	if( status ) printerror(status,"putColToFits: le movabs a foire");
	1033	fits_get_hdu_type(fptr_, &hdutype, &status);
	1034	if(hdutype != ASCII_TBL && hdutype != BINARY_TBL)
	1035	{
	1036	cout << " hdunum= " << hdunum_ << " hdutype= " << hdutype << endl;
	1037	throw IOExc("FitsFile::putColToFits, this HDU is not an ASCII table nor a binary table");
	1038	}
	1039	if(hdutype == ASCII_TBL && nocol>0)
	1040	{
	1041	throw IOExc("FitsFile::putColToFits, this HDU is an ASCII table, nocol>0 forbidden");
	1042	}
	1043	int code;
	1044	long repeat, width;
	1045	fits_get_coltype(fptr_, nocol+1, &code, &repeat,&width, &status);
	1046	if( code != TFLOAT)
	1047	{
	1048	cout << " WARNING : types don't match (putColToFits) : on fits file= " << code << " (FITS code), to be written= FLOAT " << endl;
	1049	}
	1050	fits_write_col(fptr_,TFLOAT,nocol+1,1,1,nentries, donnees ,&status);
	1051	if( status ) printerror( status,"erreur ecriture du fichier fits" );
	1052	}
[1136]	1053	void FitsOutFile::putColToFits(int nocol, int nentries, int* donnees) const
[839]	1054	{
	1055	int status = 0;
	1056	int hdutype;
	1057	fits_movabs_hdu(fptr_,hdunum_,&hdutype,&status);
	1058	if( status ) printerror(status,"putColToFits: le movabs a foire");
	1059	fits_get_hdu_type(fptr_, &hdutype, &status);
	1060	if(hdutype != ASCII_TBL && hdutype != BINARY_TBL)
	1061	{
	1062	cout << " hdunum= " << hdunum_ << " hdutype= " << hdutype << endl;
	1063	throw IOExc("FitsFile::putColToFits, this HDU is not an ASCII table nor a binary table");
	1064	}
	1065	if(hdutype == ASCII_TBL && nocol>0)
	1066	{
	1067	throw IOExc("FitsFile::putColToFits, this HDU is an ASCII table, nocol>0 forbidden");
	1068	}
	1069	int code;
	1070	long repeat, width;
	1071	fits_get_coltype(fptr_, nocol+1, &code, &repeat,&width, &status);
	1072	if( code != TLONG && code != TINT && code != TSHORT )
	1073	{
[971]	1074	cout << " WARNING : types don't match (putColToFits) : on fits file= " << code << " (FITS code), to be written= INT " << endl;
[839]	1075	}
	1076	fits_write_col(fptr_,TINT,nocol+1,1,1,nentries, donnees ,&status);
[971]	1077	if( status ) printerror( status," ecriture du fichier fits" );
[839]	1078	}
[1136]	1079	void FitsOutFile::putColToFits(int nocol, int nentries, char** donnees) const
[839]	1080	{
	1081	int status = 0;
	1082	int hdutype;
	1083	fits_movabs_hdu(fptr_,hdunum_,&hdutype,&status);
	1084	if( status ) printerror(status,"putColToFits: le movabs a foire");
	1085	fits_get_hdu_type(fptr_, &hdutype, &status);
	1086	if(hdutype != ASCII_TBL && hdutype != BINARY_TBL)
	1087	{
	1088	cout << " hdunum= " << hdunum_ << " hdutype= " << hdutype << endl;
	1089	throw IOExc("FitsFile::putColToFits, this HDU is not an ASCII table nor a binary table");
	1090	}
	1091	if(hdutype == ASCII_TBL && nocol>0)
	1092	{
	1093	throw IOExc("FitsFile::putColToFits, this HDU is an ASCII table, nocol>0 forbidden");
	1094	}
	1095	int code;
	1096	long repeat, width;
	1097	fits_get_coltype(fptr_, nocol+1, &code, &repeat,&width, &status);
	1098	if( code != TSTRING)
	1099	{
	1100	cout << " WARNING : types don't match (putColToFits) : on fits file= " << code << " (FITS code), to be written= char** " << endl;
	1101	}
	1102	fits_write_col(fptr_,TSTRING,nocol+1,1,1,nentries, donnees ,&status);
	1103	if( status ) printerror( status,"erreur ecriture du fichier fits" );
	1104	}
	1105
[1143]	1106	void FitsOutFile::DVListIntoPrimaryHeader(DVList& dvl) const
	1107	{
	1108	int status = 0;
	1109	int hdutype;
	1110	fits_movabs_hdu(fptr_,1,&hdutype,&status);
	1111	addKeywordsOfDVList(dvl);
	1112	fits_movabs_hdu(fptr_,hdunum_,&hdutype,&status);
	1113	}
[839]	1114
[1143]	1115
[1136]	1116	void FitsOutFile::writeSignatureOnFits() const
[839]	1117	{
	1118	int status = 0;
	1119	char keyname[LEN_KEYWORD];
	1120	char strval[FLEN_VALUE];
	1121	char comment[FLEN_COMMENT];
[971]	1122	strncpy(keyname, "CREATOR", LEN_KEYWORD);
	1123	keyname[LEN_KEYWORD-1] = '\0';
	1124	strcpy(strval, "SOPHYA");
	1125	strcpy(comment," SOPHYA Package - FITSIOServer ");
	1126	fits_write_key(fptr_, TSTRING, keyname, &strval, comment, &status);
	1127	if( status ) printerror( status );
[1143]	1128	fits_write_date(fptr_, &status);
[971]	1129	fits_write_comment(fptr_,"..............................................", &status);
	1130	fits_write_comment(fptr_, " SOPHYA package - FITSIOSever ", &status);
	1131	fits_write_comment(fptr_, " (C) LAL/IN2P3-CNRS Orsay, FRANCE 2000", &status);
	1132	fits_write_comment(fptr_, " (C) DAPNIA/CEA Saclay, FRANCE 2000", &status);
	1133	fits_write_comment(fptr_,"..............................................", &status);
[1045]	1134	if( status ) printerror( status, "erreur writeSignatureOnFits" );
[839]	1135	}
	1136
[903]	1137
[1143]	1138	void FitsOutFile::addKeywordsOfDVList(DVList& dvl) const
	1139	{
	1140	int status = 0;
	1141	fits_write_comment(fptr_,"---------- keywords from SOPHYA ---------", &status);
	1142	if (hdunum_ == 1) writeSignatureOnFits();
	1143	DVList::ValList::const_iterator it;
	1144	for(it = dvl.Begin(); it != dvl.End(); it++)
	1145	{
	1146	char keytype= (*it).second.elval.typ;
	1147	char keyname[10];
	1148	strncpy(keyname,(*it).first.substr(0,64).c_str(),10);
	1149	char comment[FLEN_COMMENT];
	1150	char strval[FLEN_VALUE]= "";
	1151	char *comkey = "COMMENT";
	1152	fits_read_keyword(fptr_, keyname, strval, NULL, &status);
	1153	if (status != 0 \|\| strncmp(keyname,comkey,LEN_KEYWORD-1) == 0 )
	1154	{
	1155	status = 0;
	1156	switch (keytype)
	1157	{
	1158	case 'I' :
	1159	{
	1160	int ival=(*it).second.elval.iv;
	1161	strncpy(comment,(*it).second.elcomm.c_str(),FLEN_COMMENT );
	1162	fits_write_key(fptr_,TINT,keyname,&ival,comment,&status);
	1163	break;
	1164	}
	1165	case 'D' :
	1166	{
	1167	double dval=(*it).second.elval.dv;
	1168	strncpy(comment,(*it).second.elcomm.c_str(),FLEN_COMMENT );
	1169	fits_write_key(fptr_,TDOUBLE,keyname,&dval,comment,&status);
	1170	break;
	1171	}
	1172	case 'S' :
	1173	{
	1174	char strval[128];
	1175	strncpy(strval,(*it).second.elval.strv->c_str(),127);
	1176	strncpy(comment,(*it).second.elcomm.c_str(),FLEN_COMMENT );
	1177	fits_write_key(fptr_,TSTRING,keyname,&strval,comment,&status);
	1178	break;
	1179	}
	1180	}
	1181	}
	1182	if( status ) printerror( status,"fitsfile: probleme ecriture mot-cle du dvlist" );
	1183	}
	1184	fits_write_comment(fptr_,"--------------------------------------", &status);
	1185	}
[903]	1186
	1187
[839]	1188

Note: See TracBrowser for help on using the repository browser.

Download in other formats: