Changeset 1805 for trunk/MagicSoft/Cosy/main/MBending.cc
- Timestamp:
- 03/02/03 20:46:12 (22 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/MagicSoft/Cosy/main/MBending.cc
r1699 r1805 2 2 3 3 #include <fstream.h> 4 #include <iomanip.h> 4 5 5 6 #include <TMinuit.h> 6 7 8 #include "timer.h" 7 9 8 10 ClassImp(MBending); … … 10 12 void MBending::Reset() 11 13 { 12 fIa = 0; 13 fIe = 0; 14 fCa = 0; 15 fAn = 0; 16 fAw = 0; 17 fNrx = 0; 18 fNry = 0; 19 fCrx = 0; 20 fCry = 0; 21 fNpae = 0; 22 fEces = 0; 23 fAces = 0; 24 fEcec = 0; 25 fAcec = 0; 14 Clear(); 26 15 } 27 16 … … 104 93 val *= kDeg2Rad; 105 94 106 if (str=="IA") fIa = val; 107 if (str=="IE") fIe = val; 108 if (str=="CA") fCa = val; 109 if (str=="AN") fAn = val; 110 if (str=="AW") fAw = val; 111 if (str=="NRX") fNrx = val; 112 if (str=="NRY") fNry = val; 113 if (str=="CRX") fCrx = val; 114 if (str=="CRY") fCry = val; 115 if (str=="NPAE") fNpae = val; 116 if (str=="ECES") fEces = val; 117 if (str=="ACES") fAces = val; 118 if (str=="ECEC") fEcec = val; 119 if (str=="ACEC") fAcec = val; 95 if (str=="IA") fIa = val; 96 if (str=="IE") fIe = val; 97 if (str=="CA") fCa = val; 98 if (str=="AN") fAn = val; 99 if (str=="AW") fAw = val; 100 if (str=="NRX") fNrx = val; 101 if (str=="NRY") fNry = val; 102 if (str=="CRX") fCrx = val; 103 if (str=="CRY") fCry = val; 104 if (str=="NPAE") fNpae = val; 105 if (str=="ECES") fEces = val; 106 if (str=="ACES") fAces = val; 107 if (str=="ECEC") fEcec = val; 108 if (str=="ACEC") fAcec = val; 109 if (str=="MAGIC1") fMagic1 = val; 110 if (str=="MAGIC2") fMagic2 = val; 120 111 121 112 fin >> val; … … 125 116 } 126 117 127 void MBending::Save(const char *name) {} 118 void MBending::Save(const char *name) 119 { 120 /* 121 ! MMT 1987 July 8 122 ! T 36 7.3622 41.448 -0.0481 123 ! IA -37.5465 20.80602 124 ! IE -13.9180 1.25217 125 ! NPAE +7.0751 26.44763 126 ! CA -6.9149 32.05358 127 ! AN +0.5053 1.40956 128 ! AW -2.2016 1.37480 129 ! END 130 */ 131 132 ofstream fout(name); 133 if (!fout) 134 { 135 cout << "Error: Cannot open file '" << name << "'" << endl; 136 return; 137 } 138 139 Timer t; 140 t.Now(); 141 142 fout << "MAGIC1 " << t.GetTimeStr() << endl; 143 fout << "S 00 000000 000000 0000000" << endl; 144 fout << setprecision(8); 145 fout << " IA " << kRad2Deg*fIa << " -1" << endl; 146 fout << " IE " << kRad2Deg*fIe << " -1" << endl; 147 fout << " CA " << kRad2Deg*fNpae << " -1" << endl; 148 fout << " NPAE " << kRad2Deg*fCa << " -1" << endl; 149 fout << " AN " << kRad2Deg*fAn << " -1" << endl; 150 fout << " AW " << kRad2Deg*fAw << " -1" << endl; 151 fout << " NRX " << kRad2Deg*fNrx << " -1" << endl; 152 fout << " NRY " << kRad2Deg*fNry << " -1" << endl; 153 fout << " CRX " << kRad2Deg*fCrx << " -1" << endl; 154 fout << " CRY " << kRad2Deg*fCry << " -1" << endl; 155 fout << " ECES " << kRad2Deg*fEces << " -1" << endl; 156 fout << " ACES " << kRad2Deg*fAces << " -1" << endl; 157 fout << " ECEC " << kRad2Deg*fEcec << " -1" << endl; 158 fout << " ACEC " << kRad2Deg*fAcec << " -1" << endl; 159 fout << " MAGIC1 " << kRad2Deg*fMagic1 << " -1" << endl; 160 fout << " MAGIC2 " << kRad2Deg*fMagic2 << " -1" << endl; 161 fout << "END" << endl; 162 } 128 163 129 164 AltAz MBending::Correct(const AltAz &aa) const … … 133 168 AltAz p = aa; 134 169 135 const AltAz CES(-fEces* kDeg2Rad*sin(p.Alt()), -fAces*kDeg2Rad*sin(p.Az()));136 const AltAz CEC(-fEcec* kDeg2Rad*cos(p.Alt()), -fAcec*kDeg2Rad*cos(p.Az()));170 const AltAz CES(-fEces*sin(p.Alt()), -fAces*sin(p.Az())); 171 const AltAz CEC(-fEcec*cos(p.Alt()), -fAcec*cos(p.Az())); 137 172 p += CES; 138 173 p += CEC; 139 174 140 const AltAz CRX(-fCrx* kDeg2Rad*sin(p.Az()-p.Alt()), fCrx*kDeg2Rad*cos(p.Az()-p.Alt())/cos(p.Alt()));141 const AltAz CRY(-fCry* kDeg2Rad*cos(p.Az()-p.Alt()), -fCry*kDeg2Rad*sin(p.Az()-p.Alt())/cos(p.Alt()));175 const AltAz CRX(-fCrx*sin(p.Az()-p.Alt()), fCrx*cos(p.Az()-p.Alt())/cos(p.Alt())); 176 const AltAz CRY(-fCry*cos(p.Az()-p.Alt()), -fCry*sin(p.Az()-p.Alt())/cos(p.Alt())); 142 177 p += CRX; 143 178 p += CRY; … … 148 183 p += NRY; 149 184 150 const AltAz AW( fAw*kDeg2Rad*sin(p.Az()), -fAw*kDeg2Rad*cos(p.Az())*tan(p.Alt())); 151 const AltAz AN(-fAn*kDeg2Rad*cos(p.Az()), -fAn*kDeg2Rad*sin(p.Az())*tan(p.Alt())); 185 const AltAz MAGIC(-fMagic1*sin(p.Az()-fMagic2), 0); 186 p += MAGIC; 187 188 const AltAz AW( fAw*sin(p.Az()), -fAw*kDeg2Rad*cos(p.Az())*tan(p.Alt())); 189 const AltAz AN(-fAn*cos(p.Az()), -fAn*sin(p.Az())*tan(p.Alt())); 152 190 p += AW; 153 191 p += AN; 154 192 155 const AltAz CA(0, -fCa*kDeg2Rad/cos(p.Alt())); 193 //zd1 += 11.22; 194 //zd1 += 1.15*sin((az0+0.05)*3.1415/180); 195 196 const AltAz CA(0, -fCa/cos(p.Alt())); 156 197 p += CA; 157 198 158 const AltAz NPAE(0, -fNpae* kDeg2Rad*tan(p.Alt()));199 const AltAz NPAE(0, -fNpae*tan(p.Alt())); 159 200 p += NPAE; 160 201 161 const AltAz I(fIe /**kDeg2Rad*/, fIa/**kDeg2Rad*/);202 const AltAz I(fIe, fIa); 162 203 p += I; 163 204 … … 171 212 AltAz p = aa; 172 213 173 const AltAz CES(-fEces*kDeg2Rad*sin(p.Alt()), -fAces*kDeg2Rad*sin(p.Az())); 174 const AltAz CEC(-fEcec*kDeg2Rad*cos(p.Alt()), -fAcec*kDeg2Rad*cos(p.Az())); 175 p -= CES; 176 p -= CEC; 177 178 const AltAz CRX(-fCrx*kDeg2Rad*sin(p.Az()-p.Alt()), fCrx*kDeg2Rad*cos(p.Az()-p.Alt())/cos(p.Alt())); 179 const AltAz CRY(-fCry*kDeg2Rad*cos(p.Az()-p.Alt()), -fCry*kDeg2Rad*sin(p.Az()-p.Alt())/cos(p.Alt())); 180 p -= CRX; 181 p -= CRY; 182 183 const AltAz NRX(fNrx*kDeg2Rad*sin(p.Alt()), -fNrx*kDeg2Rad); 184 const AltAz NRY(fNry*kDeg2Rad*cos(p.Alt()), -fNry*kDeg2Rad*tan(p.Alt())); 185 p -= NRX; 186 p -= NRY; 187 188 const AltAz AW( fAw*kDeg2Rad*sin(p.Az()), -fAw*kDeg2Rad*cos(p.Az())*tan(p.Alt())); 189 const AltAz AN(-fAn*kDeg2Rad*cos(p.Az()), -fAn*kDeg2Rad*sin(p.Az())*tan(p.Alt())); 190 p -= AW; 191 p -= AN; 214 const AltAz I(fIe/**kDeg2Rad*/, fIa/**kDeg2Rad*/); 215 p -= I; 216 217 const AltAz NPAE(0, -fNpae*kDeg2Rad*tan(p.Alt())); 218 p -= NPAE; 192 219 193 220 const AltAz CA(0, -fCa*kDeg2Rad/cos(p.Alt())); 194 221 p -= CA; 195 222 196 const AltAz NPAE(0, -fNpae*kDeg2Rad*tan(p.Alt())); 197 p -= NPAE; 198 199 const AltAz I(fIe/**kDeg2Rad*/, fIa/**kDeg2Rad*/); 200 p -= I; 223 const AltAz AN(-fAn*kDeg2Rad*cos(p.Az()), -fAn*kDeg2Rad*sin(p.Az())*tan(p.Alt())); 224 const AltAz AW( fAw*kDeg2Rad*sin(p.Az()), 0/*-fAw*kDeg2Rad*cos(p.Az())*tan(p.Alt())*/); 225 p -= AN; 226 p -= AW; 227 228 const AltAz NRY(fNry*kDeg2Rad*cos(p.Alt()), -fNry*kDeg2Rad*tan(p.Alt())); 229 const AltAz NRX(fNrx*kDeg2Rad*sin(p.Alt()), -fNrx*kDeg2Rad); 230 p -= NRY; 231 p -= NRX; 232 233 const AltAz CRY(-fCry*kDeg2Rad*cos(p.Az()-p.Alt()), -fCry*kDeg2Rad*sin(p.Az()-p.Alt())/cos(p.Alt())); 234 const AltAz CRX(-fCrx*kDeg2Rad*sin(p.Az()-p.Alt()), fCrx*kDeg2Rad*cos(p.Az()-p.Alt())/cos(p.Alt())); 235 p -= CRY; 236 p -= CRX; 237 238 const AltAz CEC(-fEcec*kDeg2Rad*cos(p.Alt()), -fAcec*kDeg2Rad*cos(p.Az())); 239 const AltAz CES(-fEces*kDeg2Rad*sin(p.Alt()), -fAces*kDeg2Rad*sin(p.Az())); 240 p -= CEC; 241 p -= CES; 201 242 202 243 return p; … … 221 262 } 222 263 223 void MBending::SetParameters(const Double_t *par, Int_t n =14)264 void MBending::SetParameters(const Double_t *par, Int_t n) 224 265 { 225 266 Clear(); … … 227 268 switch (n) 228 269 { 270 case 16: 271 fMagic2=par[15]/kRad2Deg; // Magic User Defined 272 case 15: 273 fMagic1=par[14]/kRad2Deg; // Magic User Defined 229 274 case 14: 230 fAcec =par[13]; // Azimuth Centering Error (cos)275 fAcec =par[13]/kRad2Deg; // Azimuth Centering Error (cos) 231 276 case 13: 232 fEcec =par[12]; // Elevation Centering Error (cos)277 fEcec =par[12]/kRad2Deg; // Elevation Centering Error (cos) 233 278 case 12: 234 fAces =par[11]; // Azimuth Centering Error (sin)279 fAces =par[11]/kRad2Deg; // Azimuth Centering Error (sin) 235 280 case 11: 236 fEces =par[10]; // Elevation Centering Error (sin)281 fEces =par[10]/kRad2Deg; // Elevation Centering Error (sin) 237 282 case 10: 238 fCry =par[9] ;// Alt/Az Coude Displacement (E-W)283 fCry =par[9]/kRad2Deg; // Alt/Az Coude Displacement (E-W) 239 284 case 9: 240 fCrx =par[8] ;// Alt/Az Coude Displacement (N-S)285 fCrx =par[8]/kRad2Deg; // Alt/Az Coude Displacement (N-S) 241 286 case 8: 242 fNry =par[7] ;// Nasmyth rotator displacement, vertical287 fNry =par[7]/kRad2Deg; // Nasmyth rotator displacement, vertical 243 288 case 7: 244 fNrx =par[6] ;// Nasmyth rotator displacement, horizontan289 fNrx =par[6]/kRad2Deg; // Nasmyth rotator displacement, horizontan 245 290 case 6: 246 fAw =par[5] ;// Azimuth Axis Misalignment (E-W)291 fAw =par[5]/kRad2Deg; // Azimuth Axis Misalignment (E-W) 247 292 case 5: 248 fAn =par[4] ;// Azimuth Axis Misalignment (N-S)293 fAn =par[4]/kRad2Deg; // Azimuth Axis Misalignment (N-S) 249 294 case 4: 250 fCa =par[3] ;// Left-Right Collimation Error295 fCa =par[3]/kRad2Deg; // Left-Right Collimation Error 251 296 case 3: 252 fNpae =par[2] ;// Az-El Nonperpendicularity297 fNpae =par[2]/kRad2Deg; // Az-El Nonperpendicularity 253 298 case 2: 254 fIe =par[1] ;// Index Error in Elevation299 fIe =par[1]/kRad2Deg; // Index Error in Elevation 255 300 case 1: 256 fIa =par[0] ;// Index Error in Azimuth257 } 258 } 259 260 void MBending::GetParameters(Double_t *par, Int_t n =14) const301 fIa =par[0]/kRad2Deg; // Index Error in Azimuth 302 } 303 } 304 305 void MBending::GetParameters(Double_t *par, Int_t n) const 261 306 { 262 307 switch (n) 263 308 { 309 case 16: 310 par[15]=fMagic2*kRad2Deg; // 311 case 15: 312 par[14]=fMagic1*kRad2Deg; // 264 313 case 14: 265 par[13]=fAcec ; // Azimuth Centering Error (cos)314 par[13]=fAcec*kRad2Deg; // Azimuth Centering Error (cos) 266 315 case 13: 267 par[12]=fEcec ; // Elevation Centering Error (cos)316 par[12]=fEcec*kRad2Deg; // Elevation Centering Error (cos) 268 317 case 12: 269 par[11]=fAces 318 par[11]=fAces*kRad2Deg; // Azimuth Centering Error (sin) 270 319 case 11: 271 par[10]=fEces 320 par[10]=fEces*kRad2Deg; // Elevation Centering Error (sin) 272 321 case 10: 273 par[9]=fCry 322 par[9]=fCry*kRad2Deg; // Alt/Az Coude Displacement (E-W) 274 323 case 9: 275 par[8]=fCrx 324 par[8]=fCrx*kRad2Deg; // Alt/Az Coude Displacement (N-S) 276 325 case 8: 277 par[7]=fNry 326 par[7]=fNry*kRad2Deg; // Nasmyth rotator displacement, vertical 278 327 case 7: 279 par[6]=fNrx 328 par[6]=fNrx*kRad2Deg; // Nasmyth rotator displacement, horizontan 280 329 case 6: 281 par[5]=fAw 330 par[5]=fAw*kRad2Deg; // Azimuth Axis Misalignment (E-W) 282 331 case 5: 283 par[4]=fAn 332 par[4]=fAn*kRad2Deg; // Azimuth Axis Misalignment (N-S) 284 333 case 4: 285 par[3]=fCa 334 par[3]=fCa*kRad2Deg; // Left-Right Collimation Error 286 335 case 3: 287 par[2]=fNpae 336 par[2]=fNpae*kRad2Deg; // Az-El Nonperpendicularity 288 337 case 2: 289 par[1]=fIe 338 par[1]=fIe*kRad2Deg; // Index Error in Elevation 290 339 case 1: 291 par[0]=fIa 340 par[0]=fIa*kRad2Deg; // Index Error in Azimuth 292 341 } 293 342 } … … 302 351 switch (n) 303 352 { 353 case 16: 354 m.mnparm(15,"MAGIC2", fMagic2*kRad2Deg, 1, -360, 360, ierflg); 355 // cout << "Init 3 CA: " << fCa << endl; 356 case 15: 357 m.mnparm(14,"MAGIC1", fMagic1*kRad2Deg, 1, -360, 360, ierflg); 358 // cout << "Init 3 CA: " << fCa << endl; 304 359 case 14: 360 m.mnparm(13,"ACEC", fAcec*kRad2Deg, 1, -360, 360, ierflg); 361 // cout << "Init 3 CA: " << fCa << endl; 305 362 case 13: 363 m.mnparm(12,"ECEC", fEcec*kRad2Deg, 1, -360, 360, ierflg); 364 // cout << "Init 3 CA: " << fCa << endl; 306 365 case 12: 366 m.mnparm(11,"ACES", fAcec*kRad2Deg, 1, -360, 360, ierflg); 367 // cout << "Init 3 CA: " << fCa << endl; 307 368 case 11: 369 m.mnparm(10,"ECES", fEcec*kRad2Deg, 1, -360, 360, ierflg); 370 // cout << "Init 3 CA: " << fCa << endl; 308 371 case 10: 372 m.mnparm(9, "CRY", fCry*kRad2Deg, 1, -360, 360, ierflg); 373 // cout << "Init 3 CA: " << fCa << endl; 309 374 case 9: 375 m.mnparm(8, "CRX", fCrx*kRad2Deg, 1, -360, 360, ierflg); 376 // cout << "Init 3 CA: " << fCa << endl; 310 377 case 8: 378 m.mnparm(7, "NRY", fNry*kRad2Deg, 1, -360, 360, ierflg); 379 // cout << "Init 3 CA: " << fCa << endl; 311 380 case 7: 381 m.mnparm(6, "NRX", fNrx*kRad2Deg, 1, -360, 360, ierflg); 382 // cout << "Init 3 CA: " << fCa << endl; 383 case 6: 384 m.mnparm(5, "AW", fAw*kRad2Deg, 1, -360, 360, ierflg); 385 // cout << "Init 3 CA: " << fCa << endl; 312 386 case 5: 387 m.mnparm(4, "AN", fAn*kRad2Deg, 1, -360, 360, ierflg); 388 // cout << "Init 3 CA: " << fCa << endl; 313 389 case 4: 314 m.mnparm(3, "CA", fCa,1, -360, 360, ierflg);390 m.mnparm(3, "CA", fCa*kRad2Deg, 1, -360, 360, ierflg); 315 391 // cout << "Init 3 CA: " << fCa << endl; 316 392 case 3: 317 m.mnparm(2, "NPAE", fNpae,1, -360, 360, ierflg);393 m.mnparm(2, "NPAE", fNpae*kRad2Deg, 1, -360, 360, ierflg); 318 394 // cout << "Init 2 NPAE: " << fNpae << endl; 319 395 case 2: 320 m.mnparm(1, "IE", fIe,1, -360, 360, ierflg);396 m.mnparm(1, "IE", fIe*kRad2Deg, 1, -360, 360, ierflg); 321 397 // cout << "Init 1 IE: " << fIe << endl; 322 398 case 1: 323 m.mnparm(0, "IA", fIa,1, -360, 360, ierflg);399 m.mnparm(0, "IA", fIa*kRad2Deg, 1, -360, 360, ierflg); 324 400 // cout << "Init 0 IA: " << fIa << endl; 325 401 } … … 335 411 switch (n) 336 412 { 413 case 16: 414 m.GetParameter(15, fMagic2, err); 415 fMagic2 /= kRad2Deg; 416 case 15: 417 m.GetParameter(14, fMagic1, err); 418 fMagic1 /= kRad2Deg; 337 419 case 14: 420 m.GetParameter(13, fAcec, err); 421 fAcec /= kRad2Deg; 338 422 case 13: 423 m.GetParameter(12, fEcec, err); 424 fEcec /= kRad2Deg; 339 425 case 12: 426 m.GetParameter(11, fAces, err); 427 fAces /= kRad2Deg; 340 428 case 11: 429 m.GetParameter(10, fEces, err); 430 fEces /= kRad2Deg; 341 431 case 10: 432 m.GetParameter(9, fCry, err); 433 fCry /= kRad2Deg; 342 434 case 9: 435 m.GetParameter(8, fCrx, err); 436 fCrx /= kRad2Deg; 343 437 case 8: 438 m.GetParameter(7, fNry, err); 439 fNry /= kRad2Deg; 344 440 case 7: 441 m.GetParameter(6, fNrx, err); 442 fNrx /= kRad2Deg; 345 443 case 6: 444 m.GetParameter(5, fAw, err); 445 fAw /= kRad2Deg; 346 446 case 5: 447 m.GetParameter(4, fAn, err); 448 fAn /= kRad2Deg; 347 449 case 4: 348 450 m.GetParameter(3, fCa, err); 451 fCa /= kRad2Deg; 349 452 case 3: 350 453 m.GetParameter(2, fNpae, err); 454 fNpae /= kRad2Deg; 351 455 case 2: 352 456 m.GetParameter(1, fIe, err); 457 fIe /= kRad2Deg; 353 458 case 1: 354 459 m.GetParameter(0, fIa, err); 460 fIa /= kRad2Deg; 355 461 } 356 462 } … … 365 471 switch (n) 366 472 { 473 case 16: 474 m.GetParameter(15, par, err); 475 cout << " 15 MAGIC2: " << par << " +- " << err << endl; 476 case 15: 477 m.GetParameter(14, par, err); 478 cout << " 14 MAGIC1: " << par << " +- " << err << endl; 367 479 case 14: 480 m.GetParameter(13, par, err); 481 cout << " 13 ACEC: " << par << " +- " << err << endl; 368 482 case 13: 483 m.GetParameter(12, par, err); 484 cout << " 12 ECEC: " << par << " +- " << err << endl; 369 485 case 12: 486 m.GetParameter(11, par, err); 487 cout << " 11 ACES: " << par << " +- " << err << endl; 370 488 case 11: 489 m.GetParameter(10, par, err); 490 cout << " 10 ECES: " << par << " +- " << err << endl; 371 491 case 10: 492 m.GetParameter(9, par, err); 493 cout << " 9 CRY: " << par << " +- " << err << endl; 372 494 case 9: 495 m.GetParameter(8, par, err); 496 cout << " 8 CRX: " << par << " +- " << err << endl; 373 497 case 8: 498 m.GetParameter(7, par, err); 499 cout << " 7 NRY: " << par << " +- " << err << endl; 374 500 case 7: 501 m.GetParameter(6, par, err); 502 cout << " 6 NRX: " << par << " +- " << err << endl; 375 503 case 6: 504 m.GetParameter(5, par, err); 505 cout << " 5 AW: " << par << " +- " << err << endl; 376 506 case 5: 507 m.GetParameter(4, par, err); 508 cout << " 4 AN: " << par << " +- " << err << endl; 377 509 case 4: 378 510 m.GetParameter(3, par, err); 379 cout << " 3 CA: " << par << " +- " << err << endl;511 cout << " 3 CA: " << par << " +- " << err << endl; 380 512 case 3: 381 513 m.GetParameter(2, par, err); 382 cout << " 2 NPAE: " << par << " +- " << err << endl;514 cout << " 2 NPAE: " << par << " +- " << err << endl; 383 515 case 2: 384 516 m.GetParameter(1, par, err); 385 cout << " 1 IE: " << par << " +- " << err << endl;517 cout << " 1 IE: " << par << " +- " << err << endl; 386 518 case 1: 387 519 m.GetParameter(0, par, err); 388 cout << " 0 IA: " << par << " +- " << err << endl;389 } 390 } 520 cout << " 0 IA: " << par << " +- " << err << endl; 521 } 522 }
Note:
See TracChangeset
for help on using the changeset viewer.